Actual source code: state.c

  1: /*
  2:      Provides utility routines for manulating any type of PETSc object.
  3: */
 4:  #include petsc.h

  8: /*@C
  9:    PetscObjectGetState - Gets the state of any PetscObject, 
 10:    regardless of the type.

 12:    Not Collective

 14:    Input Parameter:
 15: .  obj - any PETSc object, for example a Vec, Mat or KSP. This must be
 16:          cast with a (PetscObject), for example, 
 17:          PetscObjectGetState((PetscObject)mat,&state);

 19:    Output Parameter:
 20: .  state - the object state

 22:    Notes: object state is an integer which gets increased every time
 23:    the object is changed. By saving and later querying the object state
 24:    one can determine whether information about the object is still current.
 25:    Currently, state is maintained for Vec and Mat objects.

 27:    Level: advanced

 29:    seealso: PetscObjectIncreaseState, PetscObjectSetState

 31:    Concepts: state

 33: @*/
 34: PetscErrorCode PetscObjectGetState(PetscObject obj,PetscInt *state)
 35: {
 37:   if (!obj) SETERRQ(PETSC_ERR_ARG_CORRUPT,"Null object");
 38:   *state = obj->state;
 39:   return(0);
 40: }

 44: /*@C
 45:    PetscObjectSetState - Sets the state of any PetscObject, 
 46:    regardless of the type.

 48:    Not Collective

 50:    Input Parameter:
 51: +  obj - any PETSc object, for example a Vec, Mat or KSP. This must be
 52:          cast with a (PetscObject), for example, 
 53:          PetscObjectSetState((PetscObject)mat,state);
 54: -  state - the object state

 56:    Notes: This function should be used with extreme caution. There is 
 57:    essentially only one use for it: if the user calls Mat(Vec)GetRow(Array),
 58:    which increases the state, but does not alter the data, then this 
 59:    routine can be used to reset the state.

 61:    Level: advanced

 63:    seealso: PetscObjectGetState,PetscObjectIncreaseState

 65:    Concepts: state

 67: @*/
 68: PetscErrorCode PetscObjectSetState(PetscObject obj,PetscInt state)
 69: {
 71:   if (!obj) SETERRQ(PETSC_ERR_ARG_CORRUPT,"Null object");
 72:   obj->state = state;
 73:   return(0);
 74: }

 78: /*@C
 79:    PetscObjectIncreaseState - Increases the state of any PetscObject, 
 80:    regardless of the type.

 82:    Not Collective

 84:    Input Parameter:
 85: .  obj - any PETSc object, for example a Vec, Mat or KSP. This must be
 86:          cast with a (PetscObject), for example, 
 87:          PetscObjectIncreaseState((PetscObject)mat);

 89:    Notes: object state is an integer which gets increased every time
 90:    the object is changed. By saving and later querying the object state
 91:    one can determine whether information about the object is still current.
 92:    Currently, state is maintained for Vec and Mat objects.

 94:    This routine is mostly for internal use by PETSc; a developer need only
 95:    call it after explicit access to an object's internals. Routines such
 96:    as VecSet or MatScale already call this routine. It is also called, as a 
 97:    precaution, in VecRestoreArray, MatRestoreRow, MatRestoreArray.

 99:    Level: developer

101:    seealso: PetscObjectGetState

103:    Concepts: state

105: @*/
106: PetscErrorCode PetscObjectIncreaseState(PetscObject obj)
107: {
109:   if (!obj) SETERRQ(PETSC_ERR_ARG_CORRUPT,"Null object");
110:   obj->state++;
111:   return(0);
112: }

114: PetscInt globalcurrentstate = 0, globalmaxstate = 10;
115: PetscErrorCode PetscRegisterComposedData(PetscInt *id)
116: {
118:   *id = globalcurrentstate++;
119:   if (globalcurrentstate > globalmaxstate) globalmaxstate += 10;
120:   return(0);
121: }

123: PetscErrorCode PetscObjectIncreaseIntComposedData(PetscObject obj)
124: {
125:   PetscInt       *ar = obj->intcomposeddata,*new_ar;
126:   PetscInt       *ir = obj->intcomposedstate,*new_ir,n = obj->int_idmax,new_n,i;

130:   new_n = globalmaxstate;
131:   PetscMalloc(new_n*sizeof(PetscInt),&new_ar);
132:   PetscMemzero(new_ar,new_n*sizeof(PetscInt));
133:   PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
134:   PetscMemzero(new_ir,new_n*sizeof(PetscInt));
135:   if (n) {
136:     for (i=0; i<n; i++) {
137:       new_ar[i] = ar[i]; new_ir[i] = ir[i];
138:     }
139:     PetscFree(ar);
140:     PetscFree(ir);
141:   }
142:   obj->int_idmax = new_n;
143:   obj->intcomposeddata = new_ar; obj->intcomposedstate = new_ir;
144:   return(0);
145: }
146: PetscErrorCode PetscObjectIncreaseIntstarComposedData(PetscObject obj)
147: {
148:   PetscInt       **ar = obj->intstarcomposeddata,**new_ar;
149:   PetscInt       *ir = obj->intstarcomposedstate,*new_ir,n = obj->intstar_idmax,new_n,i;

153:   new_n = globalmaxstate;
154:   PetscMalloc(new_n*sizeof(PetscInt*),&new_ar);
155:   PetscMemzero(new_ar,new_n*sizeof(PetscInt*));
156:   PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
157:   PetscMemzero(new_ir,new_n*sizeof(PetscInt));
158:   if (n) {
159:     for (i=0; i<n; i++) {
160:       new_ar[i] = ar[i]; new_ir[i] = ir[i];
161:     }
162:     PetscFree(ar);
163:     PetscFree(ir);
164:   }
165:   obj->intstar_idmax = new_n;
166:   obj->intstarcomposeddata = new_ar; obj->intstarcomposedstate = new_ir;
167:   return(0);
168: }

170: PetscErrorCode PetscObjectIncreaseRealComposedData(PetscObject obj)
171: {
172:   PetscReal      *ar = obj->realcomposeddata,*new_ar;
173:   PetscInt       *ir = obj->realcomposedstate,*new_ir,n = obj->real_idmax,new_n,i;

177:   new_n = globalmaxstate;
178:   PetscMalloc(new_n*sizeof(PetscReal),&new_ar);
179:   PetscMemzero(new_ar,new_n*sizeof(PetscReal));
180:   PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
181:   PetscMemzero(new_ir,new_n*sizeof(PetscInt));
182:   if (n) {
183:     for (i=0; i<n; i++) {
184:       new_ar[i] = ar[i]; new_ir[i] = ir[i];
185:     }
186:     PetscFree(ar);
187:     PetscFree(ir);
188:   }
189:   obj->real_idmax = new_n;
190:   obj->realcomposeddata = new_ar; obj->realcomposedstate = new_ir;
191:   return(0);
192: }

194: PetscErrorCode PetscObjectIncreaseRealstarComposedData(PetscObject obj)
195: {
196:   PetscReal      **ar = obj->realstarcomposeddata,**new_ar;
197:   PetscInt       *ir = obj->realstarcomposedstate,*new_ir,n = obj->realstar_idmax,new_n,i;

201:   new_n = globalmaxstate;
202:   PetscMalloc(new_n*sizeof(PetscReal*),&new_ar);
203:   PetscMemzero(new_ar,new_n*sizeof(PetscReal*));
204:   PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
205:   PetscMemzero(new_ir,new_n*sizeof(PetscInt));
206:   if (n) {
207:     for (i=0; i<n; i++) {
208:       new_ar[i] = ar[i]; new_ir[i] = ir[i];
209:     }
210:     PetscFree(ar);
211:     PetscFree(ir);
212:   }
213:   obj->realstar_idmax = new_n;
214:   obj->realstarcomposeddata = new_ar; obj->realstarcomposedstate = new_ir;
215:   return(0);
216: }

218: PetscErrorCode PetscObjectIncreaseScalarComposedData(PetscObject obj)
219: {
220:   PetscScalar    *ar = obj->scalarcomposeddata,*new_ar;
221:   PetscInt       *ir = obj->scalarcomposedstate,*new_ir,n = obj->scalar_idmax,new_n,i;

225:   new_n = globalmaxstate;
226:   PetscMalloc(new_n*sizeof(PetscScalar),&new_ar);
227:   PetscMemzero(new_ar,new_n*sizeof(PetscScalar));
228:   PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
229:   PetscMemzero(new_ir,new_n*sizeof(PetscInt));
230:   if (n) {
231:     for (i=0; i<n; i++) {
232:       new_ar[i] = ar[i]; new_ir[i] = ir[i];
233:     }
234:     PetscFree(ar);
235:     PetscFree(ir);
236:   }
237:   obj->scalar_idmax = new_n;
238:   obj->scalarcomposeddata = new_ar; obj->scalarcomposedstate = new_ir;
239:   return(0);
240: }

242: PetscErrorCode PetscObjectIncreaseScalarstarComposedData(PetscObject obj)
243: {
244:   PetscScalar    **ar = obj->scalarstarcomposeddata,**new_ar;
245:   PetscInt       *ir = obj->scalarstarcomposedstate,*new_ir,n = obj->scalarstar_idmax,new_n,i;

249:   new_n = globalmaxstate;
250:   PetscMalloc(new_n*sizeof(PetscScalar*),&new_ar);
251:   PetscMemzero(new_ar,new_n*sizeof(PetscScalar*));
252:   PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
253:   PetscMemzero(new_ir,new_n*sizeof(PetscInt));
254:   if (n) {
255:     for (i=0; i<n; i++) {
256:       new_ar[i] = ar[i]; new_ir[i] = ir[i];
257:     }
258:     PetscFree(ar);
259:     PetscFree(ir);
260:   }
261:   obj->scalarstar_idmax = new_n;
262:   obj->scalarstarcomposeddata = new_ar; obj->scalarstarcomposedstate = new_ir;
263:   return(0);
264: }