1 #define PETSC_DLL 2 /* 3 Provides utility routines for manulating any type of PETSc object. 4 */ 5 #include "petsc.h" /*I "petsc.h" I*/ 6 7 #undef __FUNCT__ 8 #define __FUNCT__ "PetscObjectStateQuery" 9 /*@C 10 PetscObjectStateQuery - Gets the state of any PetscObject, 11 regardless of the type. 12 13 Not Collective 14 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); 19 20 Output Parameter: 21 . state - the object state 22 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. 27 28 Level: advanced 29 30 seealso: PetscObjectStateIncrease, PetscObjectSetState 31 32 Concepts: state 33 34 @*/ 35 PetscErrorCode PETSC_DLLEXPORT PetscObjectStateQuery(PetscObject obj,PetscInt *state) 36 { 37 PetscFunctionBegin; 38 if (!obj) SETERRQ(PETSC_ERR_ARG_CORRUPT,"Null object"); 39 *state = obj->state; 40 PetscFunctionReturn(0); 41 } 42 43 #undef __FUNCT__ 44 #define __FUNCT__ "PetscObjectSetState" 45 /*@C 46 PetscObjectSetState - Sets the state of any PetscObject, 47 regardless of the type. 48 49 Not Collective 50 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 56 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. 61 62 Level: advanced 63 64 seealso: PetscObjectStateQuery,PetscObjectStateIncrease 65 66 Concepts: state 67 68 @*/ 69 PetscErrorCode PETSC_DLLEXPORT PetscObjectSetState(PetscObject obj,PetscInt state) 70 { 71 PetscFunctionBegin; 72 if (!obj) SETERRQ(PETSC_ERR_ARG_CORRUPT,"Null object"); 73 obj->state = state; 74 PetscFunctionReturn(0); 75 } 76 77 PetscInt PETSC_DLLEXPORT globalcurrentstate = 0; 78 PetscInt PETSC_DLLEXPORT globalmaxstate = 10; 79 /*@C 80 PetscObjectComposedDataRegister - Get an available id for 81 composed data 82 83 Not Collective 84 85 Output parameter: 86 . id - an identifier under which data can be stored 87 88 Level: developer 89 90 seealso: PetscObjectComposedDataSetInt 91 92 @*/ 93 PetscErrorCode PETSC_DLLEXPORT PetscObjectComposedDataRegister(PetscInt *id) 94 { 95 PetscFunctionBegin; 96 *id = globalcurrentstate++; 97 if (globalcurrentstate > globalmaxstate) globalmaxstate += 10; 98 PetscFunctionReturn(0); 99 } 100 101 PetscErrorCode PETSC_DLLEXPORT PetscObjectComposedDataIncreaseInt(PetscObject obj) 102 { 103 PetscInt *ar = obj->intcomposeddata,*new_ar; 104 PetscInt *ir = obj->intcomposedstate,*new_ir,n = obj->int_idmax,new_n,i; 105 PetscErrorCode ierr; 106 107 PetscFunctionBegin; 108 new_n = globalmaxstate; 109 ierr = PetscMalloc(new_n*sizeof(PetscInt),&new_ar);CHKERRQ(ierr); 110 ierr = PetscMemzero(new_ar,new_n*sizeof(PetscInt));CHKERRQ(ierr); 111 ierr = PetscMalloc(new_n*sizeof(PetscInt),&new_ir);CHKERRQ(ierr); 112 ierr = PetscMemzero(new_ir,new_n*sizeof(PetscInt));CHKERRQ(ierr); 113 if (n) { 114 for (i=0; i<n; i++) { 115 new_ar[i] = ar[i]; new_ir[i] = ir[i]; 116 } 117 ierr = PetscFree(ar);CHKERRQ(ierr); 118 ierr = PetscFree(ir);CHKERRQ(ierr); 119 } 120 obj->int_idmax = new_n; 121 obj->intcomposeddata = new_ar; obj->intcomposedstate = new_ir; 122 PetscFunctionReturn(0); 123 } 124 PetscErrorCode PETSC_DLLEXPORT PetscObjectComposedDataIncreaseIntstar(PetscObject obj) 125 { 126 PetscInt **ar = obj->intstarcomposeddata,**new_ar; 127 PetscInt *ir = obj->intstarcomposedstate,*new_ir,n = obj->intstar_idmax,new_n,i; 128 PetscErrorCode ierr; 129 130 PetscFunctionBegin; 131 new_n = globalmaxstate; 132 ierr = PetscMalloc(new_n*sizeof(PetscInt*),&new_ar);CHKERRQ(ierr); 133 ierr = PetscMemzero(new_ar,new_n*sizeof(PetscInt*));CHKERRQ(ierr); 134 ierr = PetscMalloc(new_n*sizeof(PetscInt),&new_ir);CHKERRQ(ierr); 135 ierr = PetscMemzero(new_ir,new_n*sizeof(PetscInt));CHKERRQ(ierr); 136 if (n) { 137 for (i=0; i<n; i++) { 138 new_ar[i] = ar[i]; new_ir[i] = ir[i]; 139 } 140 ierr = PetscFree(ar);CHKERRQ(ierr); 141 ierr = PetscFree(ir);CHKERRQ(ierr); 142 } 143 obj->intstar_idmax = new_n; 144 obj->intstarcomposeddata = new_ar; obj->intstarcomposedstate = new_ir; 145 PetscFunctionReturn(0); 146 } 147 148 PetscErrorCode PETSC_DLLEXPORT PetscObjectComposedDataIncreaseReal(PetscObject obj) 149 { 150 PetscReal *ar = obj->realcomposeddata,*new_ar; 151 PetscInt *ir = obj->realcomposedstate,*new_ir,n = obj->real_idmax,new_n,i; 152 PetscErrorCode ierr; 153 154 PetscFunctionBegin; 155 new_n = globalmaxstate; 156 ierr = PetscMalloc(new_n*sizeof(PetscReal),&new_ar);CHKERRQ(ierr); 157 ierr = PetscMemzero(new_ar,new_n*sizeof(PetscReal));CHKERRQ(ierr); 158 ierr = PetscMalloc(new_n*sizeof(PetscInt),&new_ir);CHKERRQ(ierr); 159 ierr = PetscMemzero(new_ir,new_n*sizeof(PetscInt));CHKERRQ(ierr); 160 if (n) { 161 for (i=0; i<n; i++) { 162 new_ar[i] = ar[i]; new_ir[i] = ir[i]; 163 } 164 ierr = PetscFree(ar);CHKERRQ(ierr); 165 ierr = PetscFree(ir);CHKERRQ(ierr); 166 } 167 obj->real_idmax = new_n; 168 obj->realcomposeddata = new_ar; obj->realcomposedstate = new_ir; 169 PetscFunctionReturn(0); 170 } 171 172 PetscErrorCode PETSC_DLLEXPORT PetscObjectComposedDataIncreaseRealstar(PetscObject obj) 173 { 174 PetscReal **ar = obj->realstarcomposeddata,**new_ar; 175 PetscInt *ir = obj->realstarcomposedstate,*new_ir,n = obj->realstar_idmax,new_n,i; 176 PetscErrorCode ierr; 177 178 PetscFunctionBegin; 179 new_n = globalmaxstate; 180 ierr = PetscMalloc(new_n*sizeof(PetscReal*),&new_ar);CHKERRQ(ierr); 181 ierr = PetscMemzero(new_ar,new_n*sizeof(PetscReal*));CHKERRQ(ierr); 182 ierr = PetscMalloc(new_n*sizeof(PetscInt),&new_ir);CHKERRQ(ierr); 183 ierr = PetscMemzero(new_ir,new_n*sizeof(PetscInt));CHKERRQ(ierr); 184 if (n) { 185 for (i=0; i<n; i++) { 186 new_ar[i] = ar[i]; new_ir[i] = ir[i]; 187 } 188 ierr = PetscFree(ar);CHKERRQ(ierr); 189 ierr = PetscFree(ir);CHKERRQ(ierr); 190 } 191 obj->realstar_idmax = new_n; 192 obj->realstarcomposeddata = new_ar; obj->realstarcomposedstate = new_ir; 193 PetscFunctionReturn(0); 194 } 195 196 PetscErrorCode PETSC_DLLEXPORT PetscObjectComposedDataIncreaseScalar(PetscObject obj) 197 { 198 PetscScalar *ar = obj->scalarcomposeddata,*new_ar; 199 PetscInt *ir = obj->scalarcomposedstate,*new_ir,n = obj->scalar_idmax,new_n,i; 200 PetscErrorCode ierr; 201 202 PetscFunctionBegin; 203 new_n = globalmaxstate; 204 ierr = PetscMalloc(new_n*sizeof(PetscScalar),&new_ar);CHKERRQ(ierr); 205 ierr = PetscMemzero(new_ar,new_n*sizeof(PetscScalar));CHKERRQ(ierr); 206 ierr = PetscMalloc(new_n*sizeof(PetscInt),&new_ir);CHKERRQ(ierr); 207 ierr = PetscMemzero(new_ir,new_n*sizeof(PetscInt));CHKERRQ(ierr); 208 if (n) { 209 for (i=0; i<n; i++) { 210 new_ar[i] = ar[i]; new_ir[i] = ir[i]; 211 } 212 ierr = PetscFree(ar);CHKERRQ(ierr); 213 ierr = PetscFree(ir);CHKERRQ(ierr); 214 } 215 obj->scalar_idmax = new_n; 216 obj->scalarcomposeddata = new_ar; obj->scalarcomposedstate = new_ir; 217 PetscFunctionReturn(0); 218 } 219 220 PetscErrorCode PETSC_DLLEXPORT PetscObjectComposedDataIncreaseScalarstar(PetscObject obj) 221 { 222 PetscScalar **ar = obj->scalarstarcomposeddata,**new_ar; 223 PetscInt *ir = obj->scalarstarcomposedstate,*new_ir,n = obj->scalarstar_idmax,new_n,i; 224 PetscErrorCode ierr; 225 226 PetscFunctionBegin; 227 new_n = globalmaxstate; 228 ierr = PetscMalloc(new_n*sizeof(PetscScalar*),&new_ar);CHKERRQ(ierr); 229 ierr = PetscMemzero(new_ar,new_n*sizeof(PetscScalar*));CHKERRQ(ierr); 230 ierr = PetscMalloc(new_n*sizeof(PetscInt),&new_ir);CHKERRQ(ierr); 231 ierr = PetscMemzero(new_ir,new_n*sizeof(PetscInt));CHKERRQ(ierr); 232 if (n) { 233 for (i=0; i<n; i++) { 234 new_ar[i] = ar[i]; new_ir[i] = ir[i]; 235 } 236 ierr = PetscFree(ar);CHKERRQ(ierr); 237 ierr = PetscFree(ir);CHKERRQ(ierr); 238 } 239 obj->scalarstar_idmax = new_n; 240 obj->scalarstarcomposeddata = new_ar; obj->scalarstarcomposedstate = new_ir; 241 PetscFunctionReturn(0); 242 } 243 244