1# -------------------------------------------------------------------- 2 3cdef class Quad(Object): 4 """Quadrature rule for integration.""" 5 def __cinit__(self): 6 self.obj = <PetscObject*> &self.quad 7 self.quad = NULL 8 9 def view(self, Viewer viewer=None) -> None: 10 """View a `Quad` object. 11 12 Collective. 13 14 Parameters 15 ---------- 16 viewer 17 A `Viewer` to display the graph. 18 19 See Also 20 -------- 21 petsc.PetscQuadratureView 22 23 """ 24 cdef PetscViewer vwr = NULL 25 if viewer is not None: vwr = viewer.vwr 26 CHKERR(PetscQuadratureView(self.quad, vwr)) 27 28 def create(self, comm: Comm | None = None) -> Self: 29 """Create a `Quad` object. 30 31 Collective. 32 33 Parameters 34 ---------- 35 comm 36 MPI communicator, defaults to `Sys.getDefaultComm`. 37 38 See Also 39 -------- 40 petsc.PetscQuadratureCreate 41 42 """ 43 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 44 cdef PetscQuadrature newquad = NULL 45 CHKERR(PetscQuadratureCreate(ccomm, &newquad)) 46 CHKERR(PetscCLEAR(self.obj)); self.quad = newquad 47 return self 48 49 def duplicate(self) -> Quad: 50 """Create a deep copy of the `Quad` object. 51 52 Collective. 53 54 See Also 55 -------- 56 petsc.PetscQuadratureDuplicate 57 58 """ 59 cdef Quad newquad = Quad() 60 CHKERR(PetscQuadratureDuplicate(self.quad, &newquad.quad)) 61 return newquad 62 63 def destroy(self) -> Self: 64 """Destroy the `Quad` object. 65 66 Collective. 67 68 See Also 69 -------- 70 petsc.PetscQuadratureDestroy 71 72 """ 73 CHKERR(PetscQuadratureDestroy(&self.quad)) 74 return self 75 76 def getData(self) -> tuple[ArrayReal, ArrayReal]: 77 """Return the data defining the `Quad`. 78 79 Not collective. 80 81 Returns 82 ------- 83 points : ArrayReal 84 The coordinates of the quadrature points. 85 weights : ArrayReal 86 The quadrature weights. 87 88 See Also 89 -------- 90 petsc.PetscQuadratureGetData 91 92 """ 93 cdef PetscInt cdim = 0 94 cdef PetscInt cnc = 0 95 cdef PetscInt cnpoints = 0 96 cdef const PetscReal *cpoints = NULL 97 cdef const PetscReal *cweights = NULL 98 CHKERR(PetscQuadratureGetData(self.quad, &cdim, &cnc, &cnpoints, &cpoints, &cweights)) 99 return array_r(cnpoints*cdim, cpoints), array_r(cnpoints*cnc, cweights) 100 101 # FIXME: 102 # def setData(???) 103 104 def getNumComponents(self) -> int: 105 """Return the number of components for functions to be integrated. 106 107 Not collective. 108 109 See Also 110 -------- 111 setNumComponents, petsc.PetscQuadratureGetNumComponents 112 113 """ 114 cdef PetscInt cnc = 0 115 CHKERR(PetscQuadratureGetNumComponents(self.quad, &cnc)) 116 return toInt(cnc) 117 118 def setNumComponents(self, nc: int) -> None: 119 """Return the number of components for functions to be integrated. 120 121 Not collective. 122 123 Parameters 124 ---------- 125 nc 126 The number of components. 127 128 See Also 129 -------- 130 getNumComponents, petsc.PetscQuadratureSetNumComponents 131 132 """ 133 cdef PetscInt cnc = asInt(nc) 134 CHKERR(PetscQuadratureSetNumComponents(self.quad, cnc)) 135 136 def getOrder(self) -> int: 137 """Return the order of the method in the `Quad`. 138 139 Not collective. 140 141 See Also 142 -------- 143 setOrder, petsc.PetscQuadratureGetOrder 144 145 """ 146 cdef PetscInt corder = 0 147 CHKERR(PetscQuadratureGetOrder(self.quad, &corder)) 148 return toInt(corder) 149 150 def setOrder(self, order: int) -> None: 151 """Set the order of the method in the `Quad`. 152 153 Not collective. 154 155 Parameters 156 ---------- 157 order 158 The order of the quadrature, i.e. the highest degree polynomial 159 that is exactly integrated. 160 161 See Also 162 -------- 163 getOrder, petsc.PetscQuadratureSetOrder 164 165 """ 166 cdef PetscInt corder = asInt(order) 167 CHKERR(PetscQuadratureSetOrder(self.quad, corder)) 168 169 170# -------------------------------------------------------------------- 171