1# Copyright (c) 2017, Lawrence Livermore National Security, LLC. Produced at 2# the Lawrence Livermore National Laboratory. LLNL-CODE-734707. All Rights 3# reserved. See files LICENSE and NOTICE for details. 4# 5# This file is part of CEED, a collection of benchmarks, miniapps, software 6# libraries and APIs for efficient high-order finite element and spectral 7# element discretizations for exascale applications. For more information and 8# source code availability see http://github.com/ceed. 9# 10# The CEED research is supported by the Exascale Computing Project 17-SC-20-SC, 11# a collaborative effort of two U.S. Department of Energy organizations (Office 12# of Science and the National Nuclear Security Administration) responsible for 13# the planning and preparation of a capable exascale ecosystem, including 14# software, applications, hardware, advanced system engineering and early 15# testbed platforms, in support of the nation's exascale computing imperative. 16 17from _ceed_cffi import ffi, lib 18import sys 19import io 20from abc import ABC 21from .ceed_vector import Vector 22from .ceed_basis import BasisTensorH1, BasisTensorH1Lagrange, BasisH1 23from .ceed_elemrestriction import ElemRestriction, IdentityElemRestriction, BlockedElemRestriction 24from .ceed_qfunction import QFunction, QFunctionByName, IdentityQFunction 25from .ceed_operator import Operator, CompositeOperator 26from .ceed_constants import * 27 28# ------------------------------------------------------------------------------ 29class Ceed(): 30 """Ceed: core components.""" 31 # Attributes 32 _pointer = ffi.NULL 33 34 # Constructor 35 def __init__(self, resource = "/cpu/self"): 36 # libCEED object 37 self._pointer = ffi.new("Ceed *") 38 39 # libCEED call 40 resourceAscii = ffi.new("char[]", resource.encode("ascii")) 41 lib.CeedInit(resourceAscii, self._pointer) 42 43 # Representation 44 def __repr__(self): 45 return "<Ceed instance at " + hex(id(self)) + ">" 46 47 # Get Resource 48 def get_resource(self): 49 """Get the full resource name for a Ceed context. 50 51 Returns: 52 resource: resource name""" 53 54 # libCEED call 55 resource = ffi.new("char **") 56 lib.CeedGetResource(self._pointer[0], resource) 57 58 return ffi.string(resource[0]).decode("UTF-8") 59 60 # Preferred MemType 61 def get_preferred_memtype(self): 62 """Return Ceed preferred memory type. 63 64 Returns: 65 memtype: Ceed preferred memory type""" 66 67 # libCEED call 68 memtype = ffi.new("CeedMemType *", MEM_HOST) 69 lib.CeedGetPreferredMemType(self._pointer[0], memtype) 70 71 return memtype[0] 72 73 # CeedVector 74 def Vector(self, size): 75 """Ceed Vector: storing and manipulating vectors. 76 77 Args: 78 size: length of vector 79 80 Returns: 81 vector: Ceed Vector""" 82 83 return Vector(self, size) 84 85 # CeedElemRestriction 86 def ElemRestriction(self, nelem, elemsize, nnodes, ncomp, indices, 87 memtype=lib.CEED_MEM_HOST, cmode=lib.CEED_COPY_VALUES, 88 imode=lib.CEED_NONINTERLACED): 89 """Ceed ElemRestriction: restriction from local vectors to elements. 90 91 Args: 92 nelem: number of elements described in the indices array 93 elemsize: size (number of nodes) per element 94 nnodes: the number of nodes in the local vector. The input Ceed Vector 95 to which the restriction will be applied is of size 96 (nnodes * ncomp). This size may include data 97 used by other Ceed ElemRestriction objects describing 98 different types of elements. 99 ncomp: number of field components per interpolation node (1 for scalar fields) 100 *indices: Numpy array of shape [nelem, elemsize]. Row i holds the 101 ordered list of the indices (into the input Ceed Vector) 102 for the unknowns corresponding to element i, where 103 0 <= i < nelem. All indices must be in the range 104 [0, nnodes - 1]. 105 **memtype: memory type of the indices array, default CEED_MEM_HOST 106 **cmode: copy mode for the indices array, default CEED_COPY_VALUES 107 **imode: ordering of the ncomp components, i.e. it specifies 108 the ordering of the components of the local vector used 109 by this CeedElemRestriction. CEED_NONINTERLACED indicates 110 the component is the outermost index and CEED_INTERLACED 111 indicates the component is the innermost index in 112 ordering of the local vector, default CEED_NONINTERLACED 113 114 Returns: 115 elemrestriction: Ceed ElemRestiction""" 116 117 return ElemRestriction(self, nelem, elemsize, nnodes, ncomp, indices, 118 memtype=memtype, cmode=cmode, imode=imode) 119 120 def IdentityElemRestriction(self, nelem, elemsize, nnodes, ncomp, 121 imode=lib.CEED_NONINTERLACED): 122 """Ceed Identity ElemRestriction: identity restriction from local vectors to elements. 123 124 Args: 125 nelem: number of elements described in the indices array 126 elemsize: size (number of nodes) per element 127 nnodes: the number of nodes in the local vector. The input Ceed Vector 128 to which the restriction will be applied is of size 129 (nnodes * ncomp). This size may include data 130 used by other Ceed ElemRestriction objects describing 131 different types of elements. 132 ncomp: number of field components per interpolation node 133 (1 for scalar fields) 134 **imode: ordering of the ncomp components, i.e. it specifies 135 the ordering of the components of the local vector used 136 by this CeedElemRestriction. CEED_NONINTERLACED indicates 137 the component is the outermost index and CEED_INTERLACED 138 indicates the component is the innermost index in 139 ordering of the local vector, default CEED_NONINTERLACED 140 141 Returns: 142 elemrestriction: Ceed Identity ElemRestiction""" 143 144 return IdentityElemRestriction(self, nelem, elemsize, nnodes, ncomp, 145 imode=imode) 146 147 def BlockedElemRestriction(self, nelem, elemsize, blksize, nnodes, ncomp, 148 indices, memtype=lib.CEED_MEM_HOST, 149 cmode=lib.CEED_COPY_VALUES, 150 imode=lib.CEED_NONINTERLACED): 151 """Ceed Blocked ElemRestriction: blocked restriction from local vectors to elements. 152 153 Args: 154 nelem: number of elements described in the indices array 155 elemsize: size (number of nodes) per element 156 blksize: number of elements in a block 157 nnodes: the number of nodes in the local vector. The input Ceed Vector 158 to which the restriction will be applied is of size 159 (nnodes * ncomp). This size may include data 160 used by other Ceed ElemRestriction objects describing 161 different types of elements. 162 ncomp: number of field components per interpolation node (1 for scalar fields) 163 *indices: Numpy array of shape [nelem, elemsize]. Row i holds the 164 ordered list of the indices (into the input Ceed Vector) 165 for the unknowns corresponding to element i, where 166 0 <= i < nelem. All indices must be in the range 167 [0, nnodes - 1]. The backend will permute and pad this 168 array to the desired ordering for the blocksize, which is 169 typically given by the backend. The default reordering is 170 to interlace elements. 171 **memtype: memory type of the indices array, default CEED_MEM_HOST 172 **cmode: copy mode for the indices array, default CEED_COPY_VALUES 173 **imode: ordering of the ncomp components, i.e. it specifies 174 the ordering of the components of the local vector used 175 by this CeedElemRestriction. CEED_NONINTERLACED indicates 176 the component is the outermost index and CEED_INTERLACED 177 indicates the component is the innermost index in 178 ordering of the local vector, default CEED_NONINTERLACED 179 180 Returns: 181 elemrestriction: Ceed Blocked ElemRestiction""" 182 183 return BlockedElemRestriction(self, nelem, elemsize, blksize, nnodes, 184 ncomp, indices, memtype=memtype, 185 cmode=cmode, imode=imode) 186 187 # CeedBasis 188 def BasisTensorH1(self, dim, ncomp, P1d, Q1d, interp1d, grad1d, 189 qref1d, qweight1d): 190 """Ceed Tensor H1 Basis: finite element tensor-product basis objects for 191 H^1 discretizations. 192 193 Args: 194 dim: topological dimension 195 ncomp: number of field components (1 for scalar fields) 196 P1d: number of nodes in one dimension 197 Q1d: number of quadrature points in one dimension 198 *interp1d: Numpy array holding the row-major (Q1d * P1d) matrix expressing the 199 values of nodal basis functions at quadrature points 200 *grad1d: Numpy array holding the row-major (Q1d * P1d) matrix expressing the 201 derivatives of nodal basis functions at quadrature points 202 *qref1d: Numpy array of length Q1d holding the locations of quadrature points 203 on the 1D reference element [-1, 1] 204 *qweight1d: Numpy array of length Q1d holding the quadrature weights on the 205 reference element 206 207 Returns: 208 basis: Ceed Basis""" 209 210 return BasisTensorH1(self, dim, ncomp, P1d, Q1d, interp1d, grad1d, 211 qref1d, qweight1d) 212 213 def BasisTensorH1Lagrange(self, dim, ncomp, P, Q, qmode): 214 """Ceed Tensor H1 Lagrange Basis: finite element tensor-product Lagrange basis 215 objects for H^1 discretizations. 216 217 Args: 218 dim: topological dimension 219 ncomp: number of field components (1 for scalar fields) 220 P: number of Gauss-Lobatto nodes in one dimension. The 221 polynomial degree of the resulting Q_k element is k=P-1. 222 Q: number of quadrature points in one dimension 223 qmode: distribution of the Q quadrature points (affects order of 224 accuracy for the quadrature) 225 226 Returns: 227 basis: Ceed Basis""" 228 229 return BasisTensorH1Lagrange(self, dim, ncomp, P, Q, qmode) 230 231 def BasisH1(self, topo, ncomp, nnodes, nqpts, interp, grad, qref, qweight): 232 """Ceed H1 Basis: finite element non tensor-product basis for H^1 discretizations. 233 234 Args: 235 topo: topology of the element, e.g. hypercube, simplex, etc 236 ncomp: number of field components (1 for scalar fields) 237 nnodes: total number of nodes 238 nqpts: total number of quadrature points 239 *interp: Numpy array holding the row-major (nqpts * nnodes) matrix expressing 240 the values of nodal basis functions at quadrature points 241 *grad: Numpy array holding the row-major (nqpts * dim * nnodes) matrix 242 expressing the derivatives of nodal basis functions at 243 quadrature points 244 *qref: Numpy array of length (nqpts * dim) holding the locations of quadrature 245 points on the reference element [-1, 1] 246 *qweight: Numpy array of length nnodes holding the quadrature weights on the 247 reference element 248 249 Returns: 250 basis: Ceed Basis""" 251 252 return BasisH1(self, topo, ncomp, nnodes, nqpts, interp, grad, qref, qweight) 253 254 # CeedQFunction 255 def QFunction(self, vlength, f, source): 256 """Ceed QFunction: point-wise operation at quadrature points for evaluating 257 volumetric terms. 258 259 Args: 260 vlength: vector length. Caller must ensure that number of quadrature 261 points is a multiple of vlength 262 f: ctypes function pointer to evaluate action at quadrature points 263 source: absolute path to source of QFunction, 264 "\\abs_path\\file.h:function_name 265 266 Returns: 267 qfunction: Ceed QFunction""" 268 269 return QFunction(self, vlength, f, source) 270 271 def QFunctionByName(self, name): 272 """Ceed QFunction By Name: point-wise operation at quadrature points 273 from a given gallery, for evaluating volumetric terms. 274 275 Args: 276 name: name of QFunction to use from gallery 277 278 Returns: 279 qfunction: Ceed QFunction By Name""" 280 281 return QFunctionByName(self, name) 282 283 def IdentityQFunction(self, size, inmode, outmode): 284 """Ceed Idenity QFunction: identity qfunction operation. 285 286 Args: 287 size: size of the qfunction fields 288 **inmode: CeedEvalMode for input to Ceed QFunction 289 **outmode: CeedEvalMode for output to Ceed QFunction 290 291 Returns: 292 qfunction: Ceed Identity QFunction""" 293 294 return IdentityQFunction(self, size, inmode, outmode) 295 296 # CeedOperator 297 def Operator(self, qf, dqf=None, qdfT=None): 298 """Ceed Operator: composed FE-type operations on vectors. 299 300 Args: 301 qf: QFunction defining the action of the operator at quadrature points 302 **dqf: QFunction defining the action of the Jacobian, default None 303 **dqfT: QFunction defining the action of the transpose of the Jacobian, 304 default None 305 306 Returns: 307 operator: Ceed Operator""" 308 309 return Operator(self, qf, dqf, qdfT) 310 311 def CompositeOperator(self): 312 """Composite Ceed Operator: composition of multiple CeedOperators. 313 314 Returns: 315 operator: Ceed Composite Operator""" 316 317 return CompositeOperator(self) 318 319 # Destructor 320 def __del__(self): 321 # libCEED call 322 lib.CeedDestroy(self._pointer) 323 324# ------------------------------------------------------------------------------ 325