Actual source code: state.c

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

  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 PETSC_DLLEXPORT PetscObjectStateQuery(PetscObject obj,PetscInt *state)
 36: {
 38:   if (!obj) SETERRQ(PETSC_ERR_ARG_CORRUPT,"Null object");
 39:   *state = obj->state;
 40:   return(0);
 41: }

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

 49:    Not Collective

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

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

 62:    Level: advanced

 64:    seealso: PetscObjectStateQuery,PetscObjectStateIncrease

 66:    Concepts: state

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

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

 83:    Not Collective

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

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

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

100:    Level: developer

102:    seealso: PetscObjectStateQuery

104:    Concepts: state

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

115: PetscInt PETSC_DLLEXPORT globalcurrentstate = 0;
116: PetscInt PETSC_DLLEXPORT globalmaxstate = 10;
117: PetscErrorCode PETSC_DLLEXPORT PetscObjectComposedDataRegister(PetscInt *id)
118: {
120:   *id = globalcurrentstate++;
121:   if (globalcurrentstate > globalmaxstate) globalmaxstate += 10;
122:   return(0);
123: }

125: PetscErrorCode PETSC_DLLEXPORT PetscObjectComposedDataIncreaseInt(PetscObject obj)
126: {
127:   PetscInt       *ar = obj->intcomposeddata,*new_ar;
128:   PetscInt       *ir = obj->intcomposedstate,*new_ir,n = obj->int_idmax,new_n,i;

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

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

172: PetscErrorCode PETSC_DLLEXPORT PetscObjectComposedDataIncreaseReal(PetscObject obj)
173: {
174:   PetscReal      *ar = obj->realcomposeddata,*new_ar;
175:   PetscInt       *ir = obj->realcomposedstate,*new_ir,n = obj->real_idmax,new_n,i;

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

196: PetscErrorCode PETSC_DLLEXPORT PetscObjectComposedDataIncreaseRealstar(PetscObject obj)
197: {
198:   PetscReal      **ar = obj->realstarcomposeddata,**new_ar;
199:   PetscInt       *ir = obj->realstarcomposedstate,*new_ir,n = obj->realstar_idmax,new_n,i;

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

220: PetscErrorCode PETSC_DLLEXPORT PetscObjectComposedDataIncreaseScalar(PetscObject obj)
221: {
222:   PetscScalar    *ar = obj->scalarcomposeddata,*new_ar;
223:   PetscInt       *ir = obj->scalarcomposedstate,*new_ir,n = obj->scalar_idmax,new_n,i;

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

244: PetscErrorCode PETSC_DLLEXPORT PetscObjectComposedDataIncreaseScalarstar(PetscObject obj)
245: {
246:   PetscScalar    **ar = obj->scalarstarcomposeddata,**new_ar;
247:   PetscInt       *ir = obj->scalarstarcomposedstate,*new_ir,n = obj->scalarstar_idmax,new_n,i;

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