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 """Ceed ElemRestriction: restriction from local vectors to elements. 89 90 Args: 91 nelem: number of elements described in the indices array 92 elemsize: size (number of nodes) per element 93 nnodes: the number of nodes in the local vector. The input Ceed Vector 94 to which the restriction will be applied is of size 95 (nnodes * ncomp). This size may include data 96 used by other Ceed ElemRestriction objects describing 97 different types of elements. 98 ncomp: number of field components per interpolation node (1 for scalar fields) 99 *indices: Numpy array of shape [nelem, elemsize]. Row i holds the 100 ordered list of the indices (into the input Ceed Vector) 101 for the unknowns corresponding to element i, where 102 0 <= i < nelem. All indices must be in the range 103 [0, nnodes - 1]. 104 **memtype: memory type of the indices array, default CEED_MEM_HOST 105 **cmode: copy mode for the indices array, default CEED_COPY_VALUES 106 107 Returns: 108 elemrestriction: Ceed ElemRestiction""" 109 110 return ElemRestriction(self, nelem, elemsize, nnodes, ncomp, indices, 111 memtype=memtype, cmode=cmode) 112 113 def IdentityElemRestriction(self, nelem, elemsize, nnodes, ncomp): 114 """Ceed Identity ElemRestriction: identity restriction from local vectors to elements. 115 116 Args: 117 nelem: number of elements described in the indices array 118 elemsize: size (number of nodes) per element 119 nnodes: the number of nodes in the local vector. The input Ceed Vector 120 to which the restriction will be applied is of size 121 (nnodes * ncomp). This size may include data 122 used by other Ceed ElemRestriction objects describing 123 different types of elements. 124 ncomp: number of field components per interpolation node (1 for scalar fields) 125 126 Returns: 127 elemrestriction: Ceed Identity ElemRestiction""" 128 129 return IdentityElemRestriction(self, nelem, elemsize, nnodes, ncomp) 130 131 def BlockedElemRestriction(self, nelem, elemsize, blksize, nnodes, ncomp, 132 indices, memtype=lib.CEED_MEM_HOST, 133 cmode=lib.CEED_COPY_VALUES): 134 """Ceed Blocked ElemRestriction: blocked restriction from local vectors to elements. 135 136 Args: 137 nelem: number of elements described in the indices array 138 elemsize: size (number of nodes) per element 139 blksize: number of elements in a block 140 nnodes: the number of nodes in the local vector. The input Ceed Vector 141 to which the restriction will be applied is of size 142 (nnodes * ncomp). This size may include data 143 used by other Ceed ElemRestriction objects describing 144 different types of elements. 145 ncomp: number of field components per interpolation node (1 for scalar fields) 146 *indices: Numpy array of shape [nelem, elemsize]. Row i holds the 147 ordered list of the indices (into the input Ceed Vector) 148 for the unknowns corresponding to element i, where 149 0 <= i < nelem. All indices must be in the range 150 [0, nnodes - 1]. The backend will permute and pad this 151 array to the desired ordering for the blocksize, which is 152 typically given by the backend. The default reordering is 153 to interlace elements. 154 **memtype: memory type of the indices array, default CEED_MEM_HOST 155 **cmode: copy mode for the indices array, default CEED_COPY_VALUES 156 157 Returns: 158 elemrestriction: Ceed Blocked ElemRestiction""" 159 160 return BlockedElemRestriction(self, nelem, elemsize, blksize, nnodes, 161 ncomp, indices, memtype=memtype, 162 cmode=cmode) 163 164 # CeedBasis 165 def BasisTensorH1(self, dim, ncomp, P1d, Q1d, interp1d, grad1d, 166 qref1d, qweight1d): 167 """Ceed Tensor H1 Basis: finite element tensor-product basis objects for 168 H^1 discretizations. 169 170 Args: 171 dim: topological dimension 172 ncomp: number of field components (1 for scalar fields) 173 P1d: number of nodes in one dimension 174 Q1d: number of quadrature points in one dimension 175 *interp1d: Numpy array holding the row-major (Q1d * P1d) matrix expressing the 176 values of nodal basis functions at quadrature points 177 *grad1d: Numpy array holding the row-major (Q1d * P1d) matrix expressing the 178 derivatives of nodal basis functions at quadrature points 179 *qref1d: Numpy array of length Q1d holding the locations of quadrature points 180 on the 1D reference element [-1, 1] 181 *qweight1d: Numpy array of length Q1d holding the quadrature weights on the 182 reference element 183 184 Returns: 185 basis: Ceed Basis""" 186 187 return BasisTensorH1(self, dim, ncomp, P1d, Q1d, interp1d, grad1d, 188 qref1d, qweight1d) 189 190 def BasisTensorH1Lagrange(self, dim, ncomp, P, Q, qmode): 191 """Ceed Tensor H1 Lagrange Basis: finite element tensor-product Lagrange basis 192 objects for H^1 discretizations. 193 194 Args: 195 dim: topological dimension 196 ncomp: number of field components (1 for scalar fields) 197 P: number of Gauss-Lobatto nodes in one dimension. The 198 polynomial degree of the resulting Q_k element is k=P-1. 199 Q: number of quadrature points in one dimension 200 qmode: distribution of the Q quadrature points (affects order of 201 accuracy for the quadrature) 202 203 Returns: 204 basis: Ceed Basis""" 205 206 return BasisTensorH1Lagrange(self, dim, ncomp, P, Q, qmode) 207 208 def BasisH1(self, topo, ncomp, nnodes, nqpts, interp, grad, qref, qweight): 209 """Ceed H1 Basis: finite element non tensor-product basis for H^1 discretizations. 210 211 Args: 212 topo: topology of the element, e.g. hypercube, simplex, etc 213 ncomp: number of field components (1 for scalar fields) 214 nnodes: total number of nodes 215 nqpts: total number of quadrature points 216 *interp: Numpy array holding the row-major (nqpts * nnodes) matrix expressing 217 the values of nodal basis functions at quadrature points 218 *grad: Numpy array holding the row-major (nqpts * dim * nnodes) matrix 219 expressing the derivatives of nodal basis functions at 220 quadrature points 221 *qref: Numpy array of length (nqpts * dim) holding the locations of quadrature 222 points on the reference element [-1, 1] 223 *qweight: Numpy array of length nnodes holding the quadrature weights on the 224 reference element 225 226 Returns: 227 basis: Ceed Basis""" 228 229 return BasisH1(self, topo, ncomp, nnodes, nqpts, interp, grad, qref, qweight) 230 231 # CeedQFunction 232 def QFunction(self, vlength, f, source): 233 """Ceed QFunction: point-wise operation at quadrature points for evaluating 234 volumetric terms. 235 236 Args: 237 vlength: vector length. Caller must ensure that number of quadrature 238 points is a multiple of vlength 239 f: ctypes function pointer to evaluate action at quadrature points 240 source: absolute path to source of QFunction, 241 "\\abs_path\\file.h:function_name 242 243 Returns: 244 qfunction: Ceed QFunction""" 245 246 return QFunction(self, vlength, f, source) 247 248 def QFunctionByName(self, name): 249 """Ceed QFunction By Name: point-wise operation at quadrature points 250 from a given gallery, for evaluating volumetric terms. 251 252 Args: 253 name: name of QFunction to use from gallery 254 255 Returns: 256 qfunction: Ceed QFunction By Name""" 257 258 return QFunctionByName(self, name) 259 260 def IdentityQFunction(self, size, inmode, outmode): 261 """Ceed Idenity QFunction: identity qfunction operation. 262 263 Args: 264 size: size of the qfunction fields 265 **inmode: CeedEvalMode for input to Ceed QFunction 266 **outmode: CeedEvalMode for output to Ceed QFunction 267 268 Returns: 269 qfunction: Ceed Identity QFunction""" 270 271 return IdentityQFunction(self, size, inmode, outmode) 272 273 # CeedOperator 274 def Operator(self, qf, dqf=None, qdfT=None): 275 """Ceed Operator: composed FE-type operations on vectors. 276 277 Args: 278 qf: QFunction defining the action of the operator at quadrature points 279 **dqf: QFunction defining the action of the Jacobian, default None 280 **dqfT: QFunction defining the action of the transpose of the Jacobian, 281 default None 282 283 Returns: 284 operator: Ceed Operator""" 285 286 return Operator(self, qf, dqf, qdfT) 287 288 def CompositeOperator(self): 289 """Composite Ceed Operator: composition of multiple CeedOperators. 290 291 Returns: 292 operator: Ceed Composite Operator""" 293 294 return CompositeOperator(self) 295 296 # Destructor 297 def __del__(self): 298 # libCEED call 299 lib.CeedDestroy(self._pointer) 300 301# ------------------------------------------------------------------------------ 302