1 /* 2 Objects which encapsulate discretizations+continuum residuals 3 */ 4 #if !defined(__PETSCDS_H) 5 #define __PETSCDS_H 6 #include <petscfe.h> 7 #include <petscfv.h> 8 #include <petscdstypes.h> 9 10 PETSC_EXTERN PetscErrorCode PetscDSInitializePackage(void); 11 12 PETSC_EXTERN PetscClassId PETSCDS_CLASSID; 13 14 /*J 15 PetscDSType - String with the name of a PETSc discrete system 16 17 Level: beginner 18 19 .seealso: PetscDSSetType(), PetscDS 20 J*/ 21 typedef const char *PetscDSType; 22 #define PETSCDSBASIC "basic" 23 24 typedef void (*PetscPointFunc)(PetscInt, PetscInt, PetscInt, 25 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 26 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 27 PetscReal, const PetscReal[], PetscScalar[]); 28 typedef void (*PetscPointJac)(PetscInt, PetscInt, PetscInt, 29 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 30 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 31 PetscReal, PetscReal, const PetscReal[], PetscScalar[]); 32 typedef void (*PetscBdPointFunc)(PetscInt, PetscInt, PetscInt, 33 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 34 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 35 PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]); 36 typedef void (*PetscBdPointJac)(PetscInt, PetscInt, PetscInt, 37 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 38 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 39 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]); 40 typedef void (*PetscRiemannFunc)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscScalar[], void *); 41 42 43 PETSC_EXTERN PetscFunctionList PetscDSList; 44 PETSC_EXTERN PetscErrorCode PetscDSCreate(MPI_Comm, PetscDS *); 45 PETSC_EXTERN PetscErrorCode PetscDSDestroy(PetscDS *); 46 PETSC_EXTERN PetscErrorCode PetscDSSetType(PetscDS, PetscDSType); 47 PETSC_EXTERN PetscErrorCode PetscDSGetType(PetscDS, PetscDSType *); 48 PETSC_EXTERN PetscErrorCode PetscDSSetUp(PetscDS); 49 PETSC_EXTERN PetscErrorCode PetscDSSetFromOptions(PetscDS); 50 PETSC_STATIC_INLINE PetscErrorCode PetscDSViewFromOptions(PetscDS A,PetscObject B,const char name[]) {return PetscObjectViewFromOptions((PetscObject)A,B,name);} 51 52 PETSC_EXTERN PetscErrorCode PetscDSView(PetscDS,PetscViewer); 53 PETSC_EXTERN PetscErrorCode PetscDSRegister(const char [], PetscErrorCode (*)(PetscDS)); 54 PETSC_EXTERN PetscErrorCode PetscDSRegisterDestroy(void); 55 56 PETSC_EXTERN PetscErrorCode PetscDSGetSpatialDimension(PetscDS, PetscInt *); 57 PETSC_EXTERN PetscErrorCode PetscDSGetNumFields(PetscDS, PetscInt *); 58 PETSC_EXTERN PetscErrorCode PetscDSGetTotalDimension(PetscDS, PetscInt *); 59 PETSC_EXTERN PetscErrorCode PetscDSGetTotalBdDimension(PetscDS, PetscInt *); 60 PETSC_EXTERN PetscErrorCode PetscDSGetTotalComponents(PetscDS, PetscInt *); 61 PETSC_EXTERN PetscErrorCode PetscDSGetFieldOffset(PetscDS, PetscInt, PetscInt *); 62 PETSC_EXTERN PetscErrorCode PetscDSGetBdFieldOffset(PetscDS, PetscInt, PetscInt *); 63 PETSC_EXTERN PetscErrorCode PetscDSGetComponentOffset(PetscDS, PetscInt, PetscInt *); 64 PETSC_EXTERN PetscErrorCode PetscDSGetComponentOffsets(PetscDS, PetscInt *[]); 65 PETSC_EXTERN PetscErrorCode PetscDSGetComponentBdOffsets(PetscDS, PetscInt *[]); 66 PETSC_EXTERN PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS, PetscInt *[]); 67 PETSC_EXTERN PetscErrorCode PetscDSGetComponentBdDerivativeOffsets(PetscDS, PetscInt *[]); 68 69 PETSC_EXTERN PetscErrorCode PetscDSGetDiscretization(PetscDS, PetscInt, PetscObject *); 70 PETSC_EXTERN PetscErrorCode PetscDSSetDiscretization(PetscDS, PetscInt, PetscObject); 71 PETSC_EXTERN PetscErrorCode PetscDSAddDiscretization(PetscDS, PetscObject); 72 PETSC_EXTERN PetscErrorCode PetscDSGetBdDiscretization(PetscDS, PetscInt, PetscObject *); 73 PETSC_EXTERN PetscErrorCode PetscDSSetBdDiscretization(PetscDS, PetscInt, PetscObject); 74 PETSC_EXTERN PetscErrorCode PetscDSAddBdDiscretization(PetscDS, PetscObject); 75 PETSC_EXTERN PetscErrorCode PetscDSGetImplicit(PetscDS, PetscInt, PetscBool*); 76 PETSC_EXTERN PetscErrorCode PetscDSSetImplicit(PetscDS, PetscInt, PetscBool); 77 PETSC_EXTERN PetscErrorCode PetscDSGetAdjacency(PetscDS, PetscInt, PetscBool*, PetscBool*); 78 PETSC_EXTERN PetscErrorCode PetscDSSetAdjacency(PetscDS, PetscInt, PetscBool, PetscBool); 79 PETSC_EXTERN PetscErrorCode PetscDSGetObjective(PetscDS, PetscInt, 80 void (**)(PetscInt, PetscInt, PetscInt, 81 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 82 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 83 PetscReal, const PetscReal[], PetscScalar[])); 84 PETSC_EXTERN PetscErrorCode PetscDSSetObjective(PetscDS, PetscInt, 85 void (*)(PetscInt, PetscInt, PetscInt, 86 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 87 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 88 PetscReal, const PetscReal[], PetscScalar[])); 89 PETSC_EXTERN PetscErrorCode PetscDSGetResidual(PetscDS, PetscInt, 90 void (**)(PetscInt, PetscInt, PetscInt, 91 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 92 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 93 PetscReal, const PetscReal[], PetscScalar[]), 94 void (**)(PetscInt, PetscInt, PetscInt, 95 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 96 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 97 PetscReal, const PetscReal[], PetscScalar[])); 98 PETSC_EXTERN PetscErrorCode PetscDSSetResidual(PetscDS, PetscInt, 99 void (*)(PetscInt, PetscInt, PetscInt, 100 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 101 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 102 PetscReal, const PetscReal[], PetscScalar[]), 103 void (*)(PetscInt, PetscInt, PetscInt, 104 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 105 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 106 PetscReal, const PetscReal[], PetscScalar[])); 107 PETSC_EXTERN PetscErrorCode PetscDSGetJacobian(PetscDS, PetscInt, PetscInt, 108 void (**)(PetscInt, PetscInt, PetscInt, 109 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 110 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 111 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 112 void (**)(PetscInt, PetscInt, PetscInt, 113 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 114 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 115 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 116 void (**)(PetscInt, PetscInt, PetscInt, 117 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 118 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 119 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 120 void (**)(PetscInt, PetscInt, PetscInt, 121 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 122 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 123 PetscReal, PetscReal, const PetscReal[], PetscScalar[])); 124 PETSC_EXTERN PetscErrorCode PetscDSSetJacobian(PetscDS, PetscInt, PetscInt, 125 void (*)(PetscInt, PetscInt, PetscInt, 126 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 127 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 128 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 129 void (*)(PetscInt, PetscInt, PetscInt, 130 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 131 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 132 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 133 void (*)(PetscInt, PetscInt, PetscInt, 134 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 135 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 136 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 137 void (*)(PetscInt, PetscInt, PetscInt, 138 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 139 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 140 PetscReal, PetscReal, const PetscReal[], PetscScalar[])); 141 PETSC_EXTERN PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS, PetscBool *); 142 PETSC_EXTERN PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS, PetscInt, PetscInt, 143 void (**)(PetscInt, PetscInt, PetscInt, 144 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 145 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 146 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 147 void (**)(PetscInt, PetscInt, PetscInt, 148 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 149 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 150 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 151 void (**)(PetscInt, PetscInt, PetscInt, 152 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 153 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 154 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 155 void (**)(PetscInt, PetscInt, PetscInt, 156 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 157 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 158 PetscReal, PetscReal, const PetscReal[], PetscScalar[])); 159 PETSC_EXTERN PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS, PetscInt, PetscInt, 160 void (*)(PetscInt, PetscInt, PetscInt, 161 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 162 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 163 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 164 void (*)(PetscInt, PetscInt, PetscInt, 165 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 166 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 167 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 168 void (*)(PetscInt, PetscInt, PetscInt, 169 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 170 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 171 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 172 void (*)(PetscInt, PetscInt, PetscInt, 173 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 174 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 175 PetscReal, PetscReal, const PetscReal[], PetscScalar[])); 176 PETSC_EXTERN PetscErrorCode PetscDSGetRiemannSolver(PetscDS, PetscInt, 177 void (**)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscScalar[], void *)); 178 PETSC_EXTERN PetscErrorCode PetscDSSetRiemannSolver(PetscDS, PetscInt, 179 void (*)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscScalar[], void *)); 180 PETSC_EXTERN PetscErrorCode PetscDSGetContext(PetscDS, PetscInt, void **); 181 PETSC_EXTERN PetscErrorCode PetscDSSetContext(PetscDS, PetscInt, void *); 182 PETSC_EXTERN PetscErrorCode PetscDSGetBdResidual(PetscDS, PetscInt, 183 void (**)(PetscInt, PetscInt, PetscInt, 184 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 185 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 186 PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]), 187 void (**)(PetscInt, PetscInt, PetscInt, 188 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 189 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 190 PetscReal, const PetscReal[], const PetscReal[], PetscScalar[])); 191 PETSC_EXTERN PetscErrorCode PetscDSSetBdResidual(PetscDS, PetscInt, 192 void (*)(PetscInt, PetscInt, PetscInt, 193 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 194 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 195 PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]), 196 void (*)(PetscInt, PetscInt, PetscInt, 197 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 198 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 199 PetscReal, const PetscReal[], const PetscReal[], PetscScalar[])); 200 PETSC_EXTERN PetscErrorCode PetscDSGetBdJacobian(PetscDS, PetscInt, PetscInt, 201 void (**)(PetscInt, PetscInt, PetscInt, 202 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 203 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 204 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]), 205 void (**)(PetscInt, PetscInt, PetscInt, 206 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 207 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 208 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]), 209 void (**)(PetscInt, PetscInt, PetscInt, 210 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 211 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 212 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]), 213 void (**)(PetscInt, PetscInt, PetscInt, 214 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 215 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 216 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[])); 217 PETSC_EXTERN PetscErrorCode PetscDSSetBdJacobian(PetscDS, PetscInt, PetscInt, 218 void (*)(PetscInt, PetscInt, PetscInt, 219 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 220 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 221 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]), 222 void (*)(PetscInt, PetscInt, PetscInt, 223 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 224 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 225 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]), 226 void (*)(PetscInt, PetscInt, PetscInt, 227 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 228 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 229 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]), 230 void (*)(PetscInt, PetscInt, PetscInt, 231 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 232 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 233 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[])); 234 PETSC_EXTERN PetscErrorCode PetscDSGetTabulation(PetscDS, PetscReal ***, PetscReal ***); 235 PETSC_EXTERN PetscErrorCode PetscDSGetBdTabulation(PetscDS, PetscReal ***, PetscReal ***); 236 PETSC_EXTERN PetscErrorCode PetscDSGetEvaluationArrays(PetscDS, PetscScalar **, PetscScalar **, PetscScalar **); 237 PETSC_EXTERN PetscErrorCode PetscDSGetWeakFormArrays(PetscDS, PetscScalar **, PetscScalar **, PetscScalar **, PetscScalar **, PetscScalar **, PetscScalar **); 238 PETSC_EXTERN PetscErrorCode PetscDSGetRefCoordArrays(PetscDS, PetscReal **, PetscScalar **); 239 PETSC_EXTERN PetscErrorCode PetscDSCopyEquations(PetscDS, PetscDS); 240 241 #endif 242