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_EXTERN PetscErrorCode PetscDSViewFromOptions(PetscDS,const char[],const char[]); 51 PETSC_EXTERN PetscErrorCode PetscDSView(PetscDS,PetscViewer); 52 PETSC_EXTERN PetscErrorCode PetscDSRegister(const char [], PetscErrorCode (*)(PetscDS)); 53 PETSC_EXTERN PetscErrorCode PetscDSRegisterDestroy(void); 54 55 PETSC_EXTERN PetscErrorCode PetscDSGetSpatialDimension(PetscDS, PetscInt *); 56 PETSC_EXTERN PetscErrorCode PetscDSGetNumFields(PetscDS, PetscInt *); 57 PETSC_EXTERN PetscErrorCode PetscDSGetTotalDimension(PetscDS, PetscInt *); 58 PETSC_EXTERN PetscErrorCode PetscDSGetTotalBdDimension(PetscDS, PetscInt *); 59 PETSC_EXTERN PetscErrorCode PetscDSGetTotalComponents(PetscDS, PetscInt *); 60 PETSC_EXTERN PetscErrorCode PetscDSGetFieldOffset(PetscDS, PetscInt, PetscInt *); 61 PETSC_EXTERN PetscErrorCode PetscDSGetBdFieldOffset(PetscDS, PetscInt, PetscInt *); 62 PETSC_EXTERN PetscErrorCode PetscDSGetComponentOffset(PetscDS, PetscInt, PetscInt *); 63 PETSC_EXTERN PetscErrorCode PetscDSGetComponentOffsets(PetscDS, PetscInt *[]); 64 PETSC_EXTERN PetscErrorCode PetscDSGetComponentBdOffsets(PetscDS, PetscInt *[]); 65 PETSC_EXTERN PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS, PetscInt *[]); 66 PETSC_EXTERN PetscErrorCode PetscDSGetComponentBdDerivativeOffsets(PetscDS, PetscInt *[]); 67 68 PETSC_EXTERN PetscErrorCode PetscDSGetDiscretization(PetscDS, PetscInt, PetscObject *); 69 PETSC_EXTERN PetscErrorCode PetscDSSetDiscretization(PetscDS, PetscInt, PetscObject); 70 PETSC_EXTERN PetscErrorCode PetscDSAddDiscretization(PetscDS, PetscObject); 71 PETSC_EXTERN PetscErrorCode PetscDSGetBdDiscretization(PetscDS, PetscInt, PetscObject *); 72 PETSC_EXTERN PetscErrorCode PetscDSSetBdDiscretization(PetscDS, PetscInt, PetscObject); 73 PETSC_EXTERN PetscErrorCode PetscDSAddBdDiscretization(PetscDS, PetscObject); 74 PETSC_EXTERN PetscErrorCode PetscDSGetImplicit(PetscDS, PetscInt, PetscBool*); 75 PETSC_EXTERN PetscErrorCode PetscDSSetImplicit(PetscDS, PetscInt, PetscBool); 76 PETSC_EXTERN PetscErrorCode PetscDSGetAdjacency(PetscDS, PetscInt, PetscBool*, PetscBool*); 77 PETSC_EXTERN PetscErrorCode PetscDSSetAdjacency(PetscDS, PetscInt, PetscBool, PetscBool); 78 PETSC_EXTERN PetscErrorCode PetscDSGetObjective(PetscDS, PetscInt, 79 void (**)(PetscInt, PetscInt, PetscInt, 80 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 81 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 82 PetscReal, const PetscReal[], PetscScalar[])); 83 PETSC_EXTERN PetscErrorCode PetscDSSetObjective(PetscDS, PetscInt, 84 void (*)(PetscInt, PetscInt, PetscInt, 85 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 86 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 87 PetscReal, const PetscReal[], PetscScalar[])); 88 PETSC_EXTERN PetscErrorCode PetscDSGetResidual(PetscDS, PetscInt, 89 void (**)(PetscInt, PetscInt, PetscInt, 90 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 91 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 92 PetscReal, const PetscReal[], PetscScalar[]), 93 void (**)(PetscInt, PetscInt, PetscInt, 94 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 95 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 96 PetscReal, const PetscReal[], PetscScalar[])); 97 PETSC_EXTERN PetscErrorCode PetscDSSetResidual(PetscDS, PetscInt, 98 void (*)(PetscInt, PetscInt, PetscInt, 99 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 100 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 101 PetscReal, const PetscReal[], PetscScalar[]), 102 void (*)(PetscInt, PetscInt, PetscInt, 103 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 104 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 105 PetscReal, const PetscReal[], PetscScalar[])); 106 PETSC_EXTERN PetscErrorCode PetscDSGetJacobian(PetscDS, PetscInt, PetscInt, 107 void (**)(PetscInt, PetscInt, PetscInt, 108 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 109 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 110 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 111 void (**)(PetscInt, PetscInt, PetscInt, 112 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 113 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 114 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 115 void (**)(PetscInt, PetscInt, PetscInt, 116 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 117 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 118 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 119 void (**)(PetscInt, PetscInt, PetscInt, 120 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 121 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 122 PetscReal, PetscReal, const PetscReal[], PetscScalar[])); 123 PETSC_EXTERN PetscErrorCode PetscDSSetJacobian(PetscDS, PetscInt, PetscInt, 124 void (*)(PetscInt, PetscInt, PetscInt, 125 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 126 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 127 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 128 void (*)(PetscInt, PetscInt, PetscInt, 129 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 130 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 131 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 132 void (*)(PetscInt, PetscInt, PetscInt, 133 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 134 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 135 PetscReal, PetscReal, const PetscReal[], PetscScalar[]), 136 void (*)(PetscInt, PetscInt, PetscInt, 137 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 138 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 139 PetscReal, PetscReal, const PetscReal[], PetscScalar[])); 140 PETSC_EXTERN PetscErrorCode PetscDSGetRiemannSolver(PetscDS, PetscInt, 141 void (**)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscScalar[], void *)); 142 PETSC_EXTERN PetscErrorCode PetscDSSetRiemannSolver(PetscDS, PetscInt, 143 void (*)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscScalar[], void *)); 144 PETSC_EXTERN PetscErrorCode PetscDSGetContext(PetscDS, PetscInt, void **); 145 PETSC_EXTERN PetscErrorCode PetscDSSetContext(PetscDS, PetscInt, void *); 146 PETSC_EXTERN PetscErrorCode PetscDSGetBdResidual(PetscDS, PetscInt, 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, const 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, const PetscReal[], const PetscReal[], PetscScalar[])); 155 PETSC_EXTERN PetscErrorCode PetscDSSetBdResidual(PetscDS, PetscInt, 156 void (*)(PetscInt, PetscInt, PetscInt, 157 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 158 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 159 PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]), 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, const PetscReal[], const PetscReal[], PetscScalar[])); 164 PETSC_EXTERN PetscErrorCode PetscDSGetBdJacobian(PetscDS, PetscInt, PetscInt, 165 void (**)(PetscInt, PetscInt, PetscInt, 166 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 167 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 168 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]), 169 void (**)(PetscInt, PetscInt, PetscInt, 170 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 171 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 172 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]), 173 void (**)(PetscInt, PetscInt, PetscInt, 174 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 175 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 176 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]), 177 void (**)(PetscInt, PetscInt, PetscInt, 178 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 179 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 180 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[])); 181 PETSC_EXTERN PetscErrorCode PetscDSSetBdJacobian(PetscDS, PetscInt, PetscInt, 182 void (*)(PetscInt, PetscInt, PetscInt, 183 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 184 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 185 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]), 186 void (*)(PetscInt, PetscInt, PetscInt, 187 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 188 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 189 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]), 190 void (*)(PetscInt, PetscInt, PetscInt, 191 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 192 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 193 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[]), 194 void (*)(PetscInt, PetscInt, PetscInt, 195 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 196 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 197 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscScalar[])); 198 PETSC_EXTERN PetscErrorCode PetscDSGetTabulation(PetscDS, PetscReal ***, PetscReal ***); 199 PETSC_EXTERN PetscErrorCode PetscDSGetBdTabulation(PetscDS, PetscReal ***, PetscReal ***); 200 PETSC_EXTERN PetscErrorCode PetscDSGetEvaluationArrays(PetscDS, PetscScalar **, PetscScalar **, PetscScalar **); 201 PETSC_EXTERN PetscErrorCode PetscDSGetWeakFormArrays(PetscDS, PetscScalar **, PetscScalar **, PetscScalar **, PetscScalar **, PetscScalar **, PetscScalar **); 202 PETSC_EXTERN PetscErrorCode PetscDSGetRefCoordArrays(PetscDS, PetscReal **, PetscScalar **); 203 204 #endif 205