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