Actual source code: state.c
2: /*
3: Provides utility routines for manulating any type of PETSc object.
4: */
5: #include <petscsys.h> /*I "petscsys.h" I*/
9: /*@C
10: PetscObjectStateQuery - Gets the state of any PetscObject,
11: regardless of the type.
13: Not Collective
15: Input Parameter:
16: . obj - any PETSc object, for example a Vec, Mat or KSP. This must be
17: cast with a (PetscObject), for example,
18: PetscObjectStateQuery((PetscObject)mat,&state);
20: Output Parameter:
21: . state - the object state
23: Notes: object state is an integer which gets increased every time
24: the object is changed. By saving and later querying the object state
25: one can determine whether information about the object is still current.
26: Currently, state is maintained for Vec and Mat objects.
28: Level: advanced
30: seealso: PetscObjectStateIncrease(), PetscObjectSetState()
32: Concepts: state
34: @*/
35: PetscErrorCode PetscObjectStateQuery(PetscObject obj,PetscInt *state)
36: {
40: *state = obj->state;
41: return(0);
42: }
46: /*@C
47: PetscObjectSetState - Sets the state of any PetscObject,
48: regardless of the type.
50: Not Collective
52: Input Parameter:
53: + obj - any PETSc object, for example a Vec, Mat or KSP. This must be
54: cast with a (PetscObject), for example,
55: PetscObjectSetState((PetscObject)mat,state);
56: - state - the object state
58: Notes: This function should be used with extreme caution. There is
59: essentially only one use for it: if the user calls Mat(Vec)GetRow(Array),
60: which increases the state, but does not alter the data, then this
61: routine can be used to reset the state.
63: Level: advanced
65: seealso: PetscObjectStateQuery(),PetscObjectStateIncrease()
67: Concepts: state
69: @*/
70: PetscErrorCode PetscObjectSetState(PetscObject obj,PetscInt state)
71: {
74: obj->state = state;
75: return(0);
76: }
78: PetscInt globalcurrentstate = 0;
79: PetscInt globalmaxstate = 10;
83: /*@C
84: PetscObjectComposedDataRegister - Get an available id for
85: composed data
87: Not Collective
89: Output parameter:
90: . id - an identifier under which data can be stored
92: Level: developer
94: seealso: PetscObjectComposedDataSetInt()
96: @*/
97: PetscErrorCode PetscObjectComposedDataRegister(PetscInt *id)
98: {
100: *id = globalcurrentstate++;
101: if (globalcurrentstate > globalmaxstate) globalmaxstate += 10;
102: return(0);
103: }
107: PetscErrorCode PetscObjectComposedDataIncreaseInt(PetscObject obj)
108: {
109: PetscInt *ar = obj->intcomposeddata,*new_ar;
110: PetscInt *ir = obj->intcomposedstate,*new_ir,n = obj->int_idmax,new_n,i;
114: new_n = globalmaxstate;
115: PetscMalloc(new_n*sizeof(PetscInt),&new_ar);
116: PetscMemzero(new_ar,new_n*sizeof(PetscInt));
117: PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
118: PetscMemzero(new_ir,new_n*sizeof(PetscInt));
119: if (n) {
120: for (i=0; i<n; i++) {
121: new_ar[i] = ar[i]; new_ir[i] = ir[i];
122: }
123: PetscFree(ar);
124: PetscFree(ir);
125: }
126: obj->int_idmax = new_n;
127: obj->intcomposeddata = new_ar; obj->intcomposedstate = new_ir;
128: return(0);
129: }
133: PetscErrorCode PetscObjectComposedDataIncreaseIntstar(PetscObject obj)
134: {
135: PetscInt **ar = obj->intstarcomposeddata,**new_ar;
136: PetscInt *ir = obj->intstarcomposedstate,*new_ir,n = obj->intstar_idmax,new_n,i;
140: new_n = globalmaxstate;
141: PetscMalloc(new_n*sizeof(PetscInt*),&new_ar);
142: PetscMemzero(new_ar,new_n*sizeof(PetscInt*));
143: PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
144: PetscMemzero(new_ir,new_n*sizeof(PetscInt));
145: if (n) {
146: for (i=0; i<n; i++) {
147: new_ar[i] = ar[i]; new_ir[i] = ir[i];
148: }
149: PetscFree(ar);
150: PetscFree(ir);
151: }
152: obj->intstar_idmax = new_n;
153: obj->intstarcomposeddata = new_ar; obj->intstarcomposedstate = new_ir;
154: return(0);
155: }
159: PetscErrorCode PetscObjectComposedDataIncreaseReal(PetscObject obj)
160: {
161: PetscReal *ar = obj->realcomposeddata,*new_ar;
162: PetscInt *ir = obj->realcomposedstate,*new_ir,n = obj->real_idmax,new_n,i;
166: new_n = globalmaxstate;
167: PetscMalloc(new_n*sizeof(PetscReal),&new_ar);
168: PetscMemzero(new_ar,new_n*sizeof(PetscReal));
169: PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
170: PetscMemzero(new_ir,new_n*sizeof(PetscInt));
171: if (n) {
172: for (i=0; i<n; i++) {
173: new_ar[i] = ar[i]; new_ir[i] = ir[i];
174: }
175: PetscFree(ar);
176: PetscFree(ir);
177: }
178: obj->real_idmax = new_n;
179: obj->realcomposeddata = new_ar; obj->realcomposedstate = new_ir;
180: return(0);
181: }
185: PetscErrorCode PetscObjectComposedDataIncreaseRealstar(PetscObject obj)
186: {
187: PetscReal **ar = obj->realstarcomposeddata,**new_ar;
188: PetscInt *ir = obj->realstarcomposedstate,*new_ir,n = obj->realstar_idmax,new_n,i;
192: new_n = globalmaxstate;
193: PetscMalloc(new_n*sizeof(PetscReal*),&new_ar);
194: PetscMemzero(new_ar,new_n*sizeof(PetscReal*));
195: PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
196: PetscMemzero(new_ir,new_n*sizeof(PetscInt));
197: if (n) {
198: for (i=0; i<n; i++) {
199: new_ar[i] = ar[i]; new_ir[i] = ir[i];
200: }
201: PetscFree(ar);
202: PetscFree(ir);
203: }
204: obj->realstar_idmax = new_n;
205: obj->realstarcomposeddata = new_ar; obj->realstarcomposedstate = new_ir;
206: return(0);
207: }
211: PetscErrorCode PetscObjectComposedDataIncreaseScalar(PetscObject obj)
212: {
213: PetscScalar *ar = obj->scalarcomposeddata,*new_ar;
214: PetscInt *ir = obj->scalarcomposedstate,*new_ir,n = obj->scalar_idmax,new_n,i;
218: new_n = globalmaxstate;
219: PetscMalloc(new_n*sizeof(PetscScalar),&new_ar);
220: PetscMemzero(new_ar,new_n*sizeof(PetscScalar));
221: PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
222: PetscMemzero(new_ir,new_n*sizeof(PetscInt));
223: if (n) {
224: for (i=0; i<n; i++) {
225: new_ar[i] = ar[i]; new_ir[i] = ir[i];
226: }
227: PetscFree(ar);
228: PetscFree(ir);
229: }
230: obj->scalar_idmax = new_n;
231: obj->scalarcomposeddata = new_ar; obj->scalarcomposedstate = new_ir;
232: return(0);
233: }
237: PetscErrorCode PetscObjectComposedDataIncreaseScalarstar(PetscObject obj)
238: {
239: PetscScalar **ar = obj->scalarstarcomposeddata,**new_ar;
240: PetscInt *ir = obj->scalarstarcomposedstate,*new_ir,n = obj->scalarstar_idmax,new_n,i;
244: new_n = globalmaxstate;
245: PetscMalloc(new_n*sizeof(PetscScalar*),&new_ar);
246: PetscMemzero(new_ar,new_n*sizeof(PetscScalar*));
247: PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
248: PetscMemzero(new_ir,new_n*sizeof(PetscInt));
249: if (n) {
250: for (i=0; i<n; i++) {
251: new_ar[i] = ar[i]; new_ir[i] = ir[i];
252: }
253: PetscFree(ar);
254: PetscFree(ir);
255: }
256: obj->scalarstar_idmax = new_n;
257: obj->scalarstarcomposeddata = new_ar; obj->scalarstarcomposedstate = new_ir;
258: return(0);
259: }