Actual source code: state.c

  1: /*$Id: gcomm.c,v 1.25 2001/03/23 23:20:38 balay Exp $*/
  2: /*
  3:      Provides utility routines for manulating any type of PETSc object.
  4: */
 5:  #include petsc.h

  9: /*@C
 10:    PetscObjectGetState - 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:          PetscObjectGetState((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: PetscObjectIncreaseState

 32:    Concepts: state

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

 45: /*@C
 46:    PetscObjectIncreaseState - Increases 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:          PetscObjectIncreaseState((PetscObject)mat);

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

 61:    This routine is mostly for internal use by PETSc; a developer need only
 62:    call it after explicit access to an object's internals. Routines such
 63:    as VecSet or MatScale already call this routine.

 65:    Level: developer

 67:    seealso: PetscObjectGetState

 69:    Concepts: state

 71: @*/
 72: int PetscObjectIncreaseState(PetscObject obj)
 73: {
 75:   if (!obj) SETERRQ(PETSC_ERR_ARG_CORRUPT,"Null object");
 76:   obj->state++;
 77:   return(0);
 78: }

 80: int globalcurrentstate = 0, globalmaxstate = 10;
 81: int PetscRegisterComposedData(int *id)
 82: {
 84:   *id = globalcurrentstate++;
 85:   if (globalcurrentstate > globalmaxstate) globalmaxstate += 10;
 86:   return(0);
 87: }

 89: int PetscObjectIncreaseIntComposedData(PetscObject obj)
 90: {
 91:   int *ar = obj->intcomposeddata,*new_ar;
 92:   int *ir = obj->intcomposedstate,*new_ir,
 93:     n = obj->int_idmax,new_n,i,ierr;

 96:   new_n = globalmaxstate;
 97:   PetscMalloc(new_n*sizeof(int),&new_ar);
 98:   PetscMemzero(new_ar,new_n*sizeof(int));
 99:   PetscMalloc(new_n*sizeof(int),&new_ir);
100:   PetscMemzero(new_ir,new_n*sizeof(int));
101:   if (n) {
102:     for (i=0; i<n; i++) {
103:       new_ar[i] = ar[i]; new_ir[i] = ir[i];
104:     }
105:     PetscFree(ar);
106:     PetscFree(ir);
107:   }
108:   obj->int_idmax = new_n;
109:   obj->intcomposeddata = new_ar; obj->intcomposedstate = new_ir;
110:   return(0);
111: }

113: int PetscObjectIncreaseRealComposedData(PetscObject obj)
114: {
115:   PetscReal *ar = obj->realcomposeddata,*new_ar;
116:   int *ir = obj->realcomposedstate,*new_ir,
117:     n = obj->real_idmax,new_n,i,ierr;

120:   new_n = globalmaxstate;
121:   PetscMalloc(new_n*sizeof(PetscReal),&new_ar);
122:   PetscMemzero(new_ar,new_n*sizeof(PetscReal));
123:   PetscMalloc(new_n*sizeof(int),&new_ir);
124:   PetscMemzero(new_ir,new_n*sizeof(int));
125:   if (n) {
126:     for (i=0; i<n; i++) {
127:       new_ar[i] = ar[i]; new_ir[i] = ir[i];
128:     }
129:     PetscFree(ar);
130:     PetscFree(ir);
131:   }
132:   obj->real_idmax = new_n;
133:   obj->realcomposeddata = new_ar; obj->realcomposedstate = new_ir;
134:   return(0);
135: }

137: int PetscObjectIncreaseScalarComposedData(PetscObject obj)
138: {
139:   PetscScalar *ar = obj->scalarcomposeddata,*new_ar;
140:   int *ir = obj->scalarcomposedstate,*new_ir,
141:     n = obj->scalar_idmax,new_n,i,ierr;

144:   new_n = globalmaxstate;
145:   PetscMalloc(new_n*sizeof(PetscScalar),&new_ar);
146:   PetscMemzero(new_ar,new_n*sizeof(PetscScalar));
147:   PetscMalloc(new_n*sizeof(int),&new_ir);
148:   PetscMemzero(new_ir,new_n*sizeof(int));
149:   if (n) {
150:     for (i=0; i<n; i++) {
151:       new_ar[i] = ar[i]; new_ir[i] = ir[i];
152:     }
153:     PetscFree(ar);
154:     PetscFree(ir);
155:   }
156:   obj->scalar_idmax = new_n;
157:   obj->scalarcomposeddata = new_ar; obj->scalarcomposedstate = new_ir;
158:   return(0);
159: }