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 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;
80: /*@C
81: PetscObjectComposedDataRegister - Get an available id for
82: composed data
84: Not Collective
86: Output parameter:
87: . id - an identifier under which data can be stored
89: Level: developer
91: seealso: PetscObjectComposedDataSetInt
93: @*/
94: PetscErrorCode PetscObjectComposedDataRegister(PetscInt *id)
95: {
97: *id = globalcurrentstate++;
98: if (globalcurrentstate > globalmaxstate) globalmaxstate += 10;
99: return(0);
100: }
102: PetscErrorCode PetscObjectComposedDataIncreaseInt(PetscObject obj)
103: {
104: PetscInt *ar = obj->intcomposeddata,*new_ar;
105: PetscInt *ir = obj->intcomposedstate,*new_ir,n = obj->int_idmax,new_n,i;
109: new_n = globalmaxstate;
110: PetscMalloc(new_n*sizeof(PetscInt),&new_ar);
111: PetscMemzero(new_ar,new_n*sizeof(PetscInt));
112: PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
113: PetscMemzero(new_ir,new_n*sizeof(PetscInt));
114: if (n) {
115: for (i=0; i<n; i++) {
116: new_ar[i] = ar[i]; new_ir[i] = ir[i];
117: }
118: PetscFree(ar);
119: PetscFree(ir);
120: }
121: obj->int_idmax = new_n;
122: obj->intcomposeddata = new_ar; obj->intcomposedstate = new_ir;
123: return(0);
124: }
125: PetscErrorCode PetscObjectComposedDataIncreaseIntstar(PetscObject obj)
126: {
127: PetscInt **ar = obj->intstarcomposeddata,**new_ar;
128: PetscInt *ir = obj->intstarcomposedstate,*new_ir,n = obj->intstar_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->intstar_idmax = new_n;
145: obj->intstarcomposeddata = new_ar; obj->intstarcomposedstate = new_ir;
146: return(0);
147: }
149: PetscErrorCode PetscObjectComposedDataIncreaseReal(PetscObject obj)
150: {
151: PetscReal *ar = obj->realcomposeddata,*new_ar;
152: PetscInt *ir = obj->realcomposedstate,*new_ir,n = obj->real_idmax,new_n,i;
156: new_n = globalmaxstate;
157: PetscMalloc(new_n*sizeof(PetscReal),&new_ar);
158: PetscMemzero(new_ar,new_n*sizeof(PetscReal));
159: PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
160: PetscMemzero(new_ir,new_n*sizeof(PetscInt));
161: if (n) {
162: for (i=0; i<n; i++) {
163: new_ar[i] = ar[i]; new_ir[i] = ir[i];
164: }
165: PetscFree(ar);
166: PetscFree(ir);
167: }
168: obj->real_idmax = new_n;
169: obj->realcomposeddata = new_ar; obj->realcomposedstate = new_ir;
170: return(0);
171: }
173: PetscErrorCode PetscObjectComposedDataIncreaseRealstar(PetscObject obj)
174: {
175: PetscReal **ar = obj->realstarcomposeddata,**new_ar;
176: PetscInt *ir = obj->realstarcomposedstate,*new_ir,n = obj->realstar_idmax,new_n,i;
180: new_n = globalmaxstate;
181: PetscMalloc(new_n*sizeof(PetscReal*),&new_ar);
182: PetscMemzero(new_ar,new_n*sizeof(PetscReal*));
183: PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
184: PetscMemzero(new_ir,new_n*sizeof(PetscInt));
185: if (n) {
186: for (i=0; i<n; i++) {
187: new_ar[i] = ar[i]; new_ir[i] = ir[i];
188: }
189: PetscFree(ar);
190: PetscFree(ir);
191: }
192: obj->realstar_idmax = new_n;
193: obj->realstarcomposeddata = new_ar; obj->realstarcomposedstate = new_ir;
194: return(0);
195: }
197: PetscErrorCode PetscObjectComposedDataIncreaseScalar(PetscObject obj)
198: {
199: PetscScalar *ar = obj->scalarcomposeddata,*new_ar;
200: PetscInt *ir = obj->scalarcomposedstate,*new_ir,n = obj->scalar_idmax,new_n,i;
204: new_n = globalmaxstate;
205: PetscMalloc(new_n*sizeof(PetscScalar),&new_ar);
206: PetscMemzero(new_ar,new_n*sizeof(PetscScalar));
207: PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
208: PetscMemzero(new_ir,new_n*sizeof(PetscInt));
209: if (n) {
210: for (i=0; i<n; i++) {
211: new_ar[i] = ar[i]; new_ir[i] = ir[i];
212: }
213: PetscFree(ar);
214: PetscFree(ir);
215: }
216: obj->scalar_idmax = new_n;
217: obj->scalarcomposeddata = new_ar; obj->scalarcomposedstate = new_ir;
218: return(0);
219: }
221: PetscErrorCode PetscObjectComposedDataIncreaseScalarstar(PetscObject obj)
222: {
223: PetscScalar **ar = obj->scalarstarcomposeddata,**new_ar;
224: PetscInt *ir = obj->scalarstarcomposedstate,*new_ir,n = obj->scalarstar_idmax,new_n,i;
228: new_n = globalmaxstate;
229: PetscMalloc(new_n*sizeof(PetscScalar*),&new_ar);
230: PetscMemzero(new_ar,new_n*sizeof(PetscScalar*));
231: PetscMalloc(new_n*sizeof(PetscInt),&new_ir);
232: PetscMemzero(new_ir,new_n*sizeof(PetscInt));
233: if (n) {
234: for (i=0; i<n; i++) {
235: new_ar[i] = ar[i]; new_ir[i] = ir[i];
236: }
237: PetscFree(ar);
238: PetscFree(ir);
239: }
240: obj->scalarstar_idmax = new_n;
241: obj->scalarstarcomposeddata = new_ar; obj->scalarstarcomposedstate = new_ir;
242: return(0);
243: }