1# -------------------------------------------------------------------- 2 3cdef extern from * nogil: 4 5 ctypedef enum PetscDMStagStencilType"DMStagStencilType": 6 DMSTAG_STENCIL_STAR 7 DMSTAG_STENCIL_BOX 8 DMSTAG_STENCIL_NONE 9 10 ctypedef enum PetscDMStagStencilLocation"DMStagStencilLocation": 11 DMSTAG_NULL_LOCATION 12 DMSTAG_BACK_DOWN_LEFT 13 DMSTAG_BACK_DOWN 14 DMSTAG_BACK_DOWN_RIGHT 15 DMSTAG_BACK_LEFT 16 DMSTAG_BACK 17 DMSTAG_BACK_RIGHT 18 DMSTAG_BACK_UP_LEFT 19 DMSTAG_BACK_UP 20 DMSTAG_BACK_UP_RIGHT 21 DMSTAG_DOWN_LEFT 22 DMSTAG_DOWN 23 DMSTAG_DOWN_RIGHT 24 DMSTAG_LEFT 25 DMSTAG_ELEMENT 26 DMSTAG_RIGHT 27 DMSTAG_UP_LEFT 28 DMSTAG_UP 29 DMSTAG_UP_RIGHT 30 DMSTAG_FRONT_DOWN_LEFT 31 DMSTAG_FRONT_DOWN 32 DMSTAG_FRONT_DOWN_RIGHT 33 DMSTAG_FRONT_LEFT 34 DMSTAG_FRONT 35 DMSTAG_FRONT_RIGHT 36 DMSTAG_FRONT_UP_LEFT 37 DMSTAG_FRONT_UP 38 DMSTAG_FRONT_UP_RIGHT 39 40 PetscErrorCode DMStagCreate1d(MPI_Comm, PetscDMBoundaryType, PetscInt, PetscInt, PetscInt, PetscDMStagStencilType, PetscInt, const PetscInt[], PetscDM*) 41 PetscErrorCode DMStagCreate2d(MPI_Comm, PetscDMBoundaryType, PetscDMBoundaryType, PetscInt, PetscInt, PetscInt, PetscInt, PetscInt, PetscInt, PetscInt, PetscDMStagStencilType, PetscInt, const PetscInt[], const PetscInt[], PetscDM*) 42 PetscErrorCode DMStagCreate3d(MPI_Comm, PetscDMBoundaryType, PetscDMBoundaryType, PetscDMBoundaryType, PetscInt, PetscInt, PetscInt, PetscInt, PetscInt, PetscInt, PetscInt, PetscInt, PetscInt, PetscInt, PetscDMStagStencilType, PetscInt, const PetscInt[], const PetscInt[], const PetscInt[], PetscDM*) 43 44 PetscErrorCode DMStagGetCorners(PetscDM, PetscInt*, PetscInt*, PetscInt*, PetscInt*, PetscInt*, PetscInt*, PetscInt*, PetscInt*, PetscInt*) 45 PetscErrorCode DMStagGetGhostCorners(PetscDM, PetscInt*, PetscInt*, PetscInt*, PetscInt*, PetscInt*, PetscInt*) 46 PetscErrorCode DMStagGetLocalSizes(PetscDM, PetscInt*, PetscInt*, PetscInt*) 47 PetscErrorCode DMStagGetEntriesPerElement(PetscDM, PetscInt*) 48 49 PetscErrorCode DMStagGetDOF(PetscDM, PetscInt*, PetscInt*, PetscInt*, PetscInt*) 50 PetscErrorCode DMStagGetNumRanks(PetscDM, PetscInt*, PetscInt*, PetscInt*) 51 PetscErrorCode DMStagGetGlobalSizes(PetscDM, PetscInt*, PetscInt*, PetscInt*) 52 PetscErrorCode DMStagGetBoundaryTypes(PetscDM, PetscDMBoundaryType*, PetscDMBoundaryType*, PetscDMBoundaryType*) 53 PetscErrorCode DMStagGetStencilWidth(PetscDM, PetscInt*) 54 PetscErrorCode DMStagGetStencilType(PetscDM, PetscDMStagStencilType*) 55 PetscErrorCode DMStagGetOwnershipRanges(PetscDM, const PetscInt*[], const PetscInt*[], const PetscInt*[]) 56 57 PetscErrorCode DMStagSetDOF(PetscDM, PetscInt, PetscInt, PetscInt, PetscInt) 58 PetscErrorCode DMStagSetNumRanks(PetscDM, PetscInt, PetscInt, PetscInt) 59 PetscErrorCode DMStagSetGlobalSizes(PetscDM, PetscInt, PetscInt, PetscInt) 60 PetscErrorCode DMStagSetBoundaryTypes(PetscDM, PetscDMBoundaryType, PetscDMBoundaryType, PetscDMBoundaryType) 61 PetscErrorCode DMStagSetStencilWidth(PetscDM, PetscInt) 62 PetscErrorCode DMStagSetStencilType(PetscDM, PetscDMStagStencilType) 63 PetscErrorCode DMStagSetOwnershipRanges(PetscDM, const PetscInt[], const PetscInt[], const PetscInt[]) 64 65 PetscErrorCode DMStagGetLocationSlot(PetscDM, PetscDMStagStencilLocation, PetscInt, PetscInt*) 66 PetscErrorCode DMStagGetLocationDOF(PetscDM, PetscDMStagStencilLocation, PetscInt*) 67 PetscErrorCode DMStagGetProductCoordinateLocationSlot(PetscDM, PetscDMStagStencilLocation, PetscInt*) 68 69 PetscErrorCode DMStagGetIsFirstRank(PetscDM, PetscBool*, PetscBool*, PetscBool*) 70 PetscErrorCode DMStagGetIsLastRank(PetscDM, PetscBool*, PetscBool*, PetscBool*) 71 72 PetscErrorCode DMStagSetUniformCoordinatesExplicit(PetscDM, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal) 73 PetscErrorCode DMStagSetUniformCoordinatesProduct(PetscDM, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal) 74 PetscErrorCode DMStagSetCoordinateDMType(PetscDM, PetscDMType) 75 PetscErrorCode DMStagSetUniformCoordinates(PetscDM, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal, PetscReal) 76 77 PetscErrorCode DMStagCreateCompatibleDMStag(PetscDM, PetscInt, PetscInt, PetscInt, PetscInt, PetscDM*) 78 PetscErrorCode DMStagVecSplitToDMDA(PetscDM, PetscVec, PetscDMStagStencilLocation, PetscInt, PetscDM*, PetscVec*) 79 PetscErrorCode DMStagMigrateVec(PetscDM, PetscVec, PetscDM, PetscVec) 80 81# -------------------------------------------------------------------- 82 83cdef inline PetscDMStagStencilType asStagStencil(object stencil) \ 84 except <PetscDMStagStencilType>(-1): 85 if isinstance(stencil, str): 86 if stencil == "star": return DMSTAG_STENCIL_STAR 87 elif stencil == "box": return DMSTAG_STENCIL_BOX 88 elif stencil == "none": return DMSTAG_STENCIL_NONE 89 else: raise ValueError("unknown stencil type: %s" % stencil) 90 return stencil 91 92cdef inline object toStagStencil(PetscDMStagStencilType stype): 93 if stype == DMSTAG_STENCIL_STAR: return "star" 94 elif stype == DMSTAG_STENCIL_BOX: return "box" 95 elif stype == DMSTAG_STENCIL_NONE: return "none" 96 97cdef inline PetscDMStagStencilLocation asStagStencilLocation(object stencil_location) \ 98 except <PetscDMStagStencilLocation>(-1): 99 if isinstance(stencil_location, str): 100 if stencil_location == "null": return DMSTAG_NULL_LOCATION 101 elif stencil_location == "back_down_left": return DMSTAG_BACK_DOWN_LEFT 102 elif stencil_location == "back_down": return DMSTAG_BACK_DOWN 103 elif stencil_location == "back_down_right": return DMSTAG_BACK_DOWN_RIGHT 104 elif stencil_location == "back_left": return DMSTAG_BACK_LEFT 105 elif stencil_location == "back": return DMSTAG_BACK 106 elif stencil_location == "back_right": return DMSTAG_BACK_RIGHT 107 elif stencil_location == "back_up_left": return DMSTAG_BACK_UP_LEFT 108 elif stencil_location == "back_up": return DMSTAG_BACK_UP 109 elif stencil_location == "back_up_right": return DMSTAG_BACK_UP_RIGHT 110 elif stencil_location == "down_left": return DMSTAG_DOWN_LEFT 111 elif stencil_location == "down": return DMSTAG_DOWN 112 elif stencil_location == "down_right": return DMSTAG_DOWN_RIGHT 113 elif stencil_location == "left": return DMSTAG_LEFT 114 elif stencil_location == "element": return DMSTAG_ELEMENT 115 elif stencil_location == "right": return DMSTAG_RIGHT 116 elif stencil_location == "up_left": return DMSTAG_UP_LEFT 117 elif stencil_location == "up": return DMSTAG_UP 118 elif stencil_location == "up_right": return DMSTAG_UP_RIGHT 119 elif stencil_location == "front_down_left": return DMSTAG_FRONT_DOWN_LEFT 120 elif stencil_location == "front_down": return DMSTAG_FRONT_DOWN 121 elif stencil_location == "front_down_right": return DMSTAG_FRONT_DOWN_RIGHT 122 elif stencil_location == "front_left": return DMSTAG_FRONT_LEFT 123 elif stencil_location == "front": return DMSTAG_FRONT 124 elif stencil_location == "front_right": return DMSTAG_FRONT_RIGHT 125 elif stencil_location == "front_up_left": return DMSTAG_FRONT_UP_LEFT 126 elif stencil_location == "front_up": return DMSTAG_FRONT_UP 127 elif stencil_location == "front_up_right": return DMSTAG_FRONT_UP_RIGHT 128 else: raise ValueError("unknown stencil location type: %s" % stencil_location) 129 return stencil_location 130 131 132cdef inline PetscInt asStagDims(dims, 133 PetscInt *_M, 134 PetscInt *_N, 135 PetscInt *_P) except? -1: 136 cdef PetscInt dim = PETSC_DECIDE 137 cdef object M=None, N=None, P=None 138 dims = tuple(dims) 139 dim = <PetscInt>len(dims) 140 if dim == 0: pass 141 elif dim == 1: M, = dims 142 elif dim == 2: M, N = dims 143 elif dim == 3: M, N, P = dims 144 if dim >= 1: _M[0] = asInt(M) 145 if dim >= 2: _N[0] = asInt(N) 146 if dim >= 3: _P[0] = asInt(P) 147 return dim 148 149cdef inline tuple toStagDims(PetscInt dim, 150 PetscInt M, 151 PetscInt N, 152 PetscInt P): 153 if dim == 0: return () 154 elif dim == 1: return (toInt(M),) 155 elif dim == 2: return (toInt(M), toInt(N)) 156 elif dim == 3: return (toInt(M), toInt(N), toInt(P)) 157 158cdef inline PetscInt asDofs(dofs, 159 PetscInt *_dof0, 160 PetscInt *_dof1, 161 PetscInt *_dof2, 162 PetscInt *_dof3) except? -1: 163 cdef PetscInt ndofs = PETSC_DECIDE 164 cdef object dof0=None, dof1=None, dof2=None, dof3=None 165 dofs = tuple(dofs) 166 ndofs = <PetscInt>len(dofs) 167 if ndofs == 2: dof0, dof1 = dofs 168 elif ndofs == 3: dof0, dof1, dof2 = dofs 169 elif ndofs == 4: dof0, dof1, dof2, dof3 = dofs 170 if ndofs >= 2: _dof0[0] = asInt(dof0) 171 if ndofs >= 2: _dof1[0] = asInt(dof1) 172 if ndofs >= 3: _dof2[0] = asInt(dof2) 173 if ndofs >= 4: _dof3[0] = asInt(dof3) 174 return ndofs 175 176cdef inline tuple toDofs(PetscInt ndofs, 177 PetscInt dof0, 178 PetscInt dof1, 179 PetscInt dof2, 180 PetscInt dof3): 181 if ndofs == 2: return (toInt(dof0), toInt(dof1)) 182 elif ndofs == 3: return (toInt(dof0), toInt(dof1), toInt(dof2)) 183 elif ndofs == 4: return (toInt(dof0), toInt(dof1), toInt(dof2), toInt(dof3)) 184 185cdef inline tuple asStagOwnershipRanges(object ownership_ranges, 186 PetscInt dim, 187 PetscInt *m, PetscInt *n, PetscInt *p, 188 PetscInt **_x, 189 PetscInt **_y, 190 PetscInt **_z): 191 cdef object ranges = list(ownership_ranges) 192 cdef PetscInt rdim = <PetscInt>len(ranges) 193 cdef PetscInt nlx=0, nly=0, nlz=0 194 if dim == PETSC_DECIDE: dim = rdim 195 elif dim != rdim: raise ValueError( 196 "number of dimensions %d and number ownership ranges %d" % 197 (toInt(dim), toInt(rdim))) 198 if dim >= 1: 199 ranges[0] = iarray_i(ranges[0], &nlx, _x) 200 if m[0] == PETSC_DECIDE: m[0] = nlx 201 elif m[0] != nlx: raise ValueError( 202 "ownership range size %d and number or processors %d" % 203 (toInt(nlx), toInt(m[0]))) 204 if dim >= 2: 205 ranges[1] = iarray_i(ranges[1], &nly, _y) 206 if n[0] == PETSC_DECIDE: n[0] = nly 207 elif n[0] != nly: raise ValueError( 208 "ownership range size %d and number or processors %d" % 209 (toInt(nly), toInt(n[0]))) 210 if dim >= 3: 211 ranges[2] = iarray_i(ranges[2], &nlz, _z) 212 if p[0] == PETSC_DECIDE: p[0] = nlz 213 elif p[0] != nlz: raise ValueError( 214 "ownership range size %d and number or processors %d" % 215 (toInt(nlz), toInt(p[0]))) 216 return tuple(ranges) 217 218 219cdef inline tuple toStagOwnershipRanges(PetscInt dim, 220 PetscInt m, PetscInt n, PetscInt p, 221 const PetscInt *lx, 222 const PetscInt *ly, 223 const PetscInt *lz): 224 # Returns tuple of arrays containing ownership ranges as Python arrays 225 ranges = [array_i(m, lx)] 226 if dim > 1: 227 ranges.append(array_i(n, ly)) 228 if dim > 2: 229 ranges.append(array_i(p, lz)) 230 return tuple(ranges) 231 232cdef inline object toStagBoundary(PetscDMBoundaryType btype): 233 if btype == DM_BOUNDARY_NONE: return "none" 234 elif btype == DM_BOUNDARY_PERIODIC: return "periodic" 235 elif btype == DM_BOUNDARY_GHOSTED: return "ghosted" 236 237cdef inline tuple toStagBoundaryTypes(PetscInt dim, PetscDMBoundaryType btx, PetscDMBoundaryType bty, PetscDMBoundaryType btz): 238 if dim == 1: return (toStagBoundary(btx),) 239 if dim == 2: return (toStagBoundary(btx), toStagBoundary(bty)) 240 if dim == 3: return (toStagBoundary(btx), toStagBoundary(bty), toStagBoundary(btz)) 241 242# -------------------------------------------------------------------- 243