xref: /libCEED/python/ceed_elemrestriction.py (revision ae718e2f2476441c07fb6c037655a83f28d5f429)
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 tempfile
19import numpy as np
20from abc import ABC
21from .ceed_constants import REQUEST_IMMEDIATE, REQUEST_ORDERED, MEM_HOST, USE_POINTER, COPY_VALUES, TRANSPOSE, NOTRANSPOSE
22from .ceed_vector import _VectorWrap
23
24# ------------------------------------------------------------------------------
25
26
27class _ElemRestrictionBase(ABC):
28    """Ceed ElemRestriction: restriction from local vectors to elements."""
29
30    # Destructor
31    def __del__(self):
32        # libCEED call
33        err_code = lib.CeedElemRestrictionDestroy(self._pointer)
34        self._ceed._check_error(err_code)
35
36    # Representation
37    def __repr__(self):
38        return "<CeedElemRestriction instance at " + hex(id(self)) + ">"
39
40    # String conversion for print() to stdout
41    def __str__(self):
42        """View an ElemRestriction via print()."""
43
44        # libCEED call
45        with tempfile.NamedTemporaryFile() as key_file:
46            with open(key_file.name, 'r+') as stream_file:
47                stream = ffi.cast("FILE *", stream_file)
48
49                err_code = lib.CeedElemRestrictionView(self._pointer[0], stream)
50                self._ceed._check_error(err_code)
51
52                stream_file.seek(0)
53                out_string = stream_file.read()
54
55        return out_string
56
57    # Apply CeedElemRestriction
58    def apply(self, u, v, tmode=NOTRANSPOSE, request=REQUEST_IMMEDIATE):
59        """Restrict a local vector to an element vector or apply its transpose.
60
61           Args:
62             u: input vector
63             v: output vector
64             **tmode: apply restriction or transpose, default CEED_NOTRANSPOSE
65             **request: Ceed request, default CEED_REQUEST_IMMEDIATE"""
66
67        # libCEED call
68        err_code = lib.CeedElemRestrictionApply(self._pointer[0], tmode, u._pointer[0],
69                                                v._pointer[0], request)
70        self._ceed._check_error(err_code)
71
72    # Transpose an ElemRestriction
73    @property
74    def T(self):
75        """Transpose an ElemRestriction."""
76
77        return TransposeElemRestriction(self)
78
79    # Transpose an ElemRestriction
80    @property
81    def transpose(self):
82        """Transpose an ElemRestriction."""
83
84        return TransposeElemRestriction(self)
85
86    # Create restriction vectors
87    def create_vector(self, createLvec=True, createEvec=True):
88        """Create Vectors associated with an ElemRestriction.
89
90           Args:
91             **createLvec: flag to create local vector, default True
92             **createEvec: flag to create element vector, default True
93
94           Returns:
95             [lvec, evec]: local vector and element vector, or None if flag set to false"""
96
97        # Vector pointers
98        lvecPointer = ffi.new("CeedVector *") if createLvec else ffi.NULL
99        evecPointer = ffi.new("CeedVector *") if createEvec else ffi.NULL
100
101        # libCEED call
102        err_code = lib.CeedElemRestrictionCreateVector(self._pointer[0], lvecPointer,
103                                                       evecPointer)
104        self._ceed._check_error(err_code)
105
106        # Return vectors
107        lvec = _VectorWrap(
108            self._ceed, lvecPointer) if createLvec else None
109        evec = _VectorWrap(
110            self._ceed, evecPointer) if createEvec else None
111
112        # Return
113        return [lvec, evec]
114
115    # Get ElemRestriction multiplicity
116    def get_multiplicity(self):
117        """Get the multiplicity of nodes in an ElemRestriction.
118
119           Returns:
120             mult: local vector containing multiplicity of nodes in ElemRestriction"""
121
122        # Create mult vector
123        [mult, evec] = self.create_vector(createEvec=False)
124        mult.set_value(0)
125
126        # libCEED call
127        err_code = lib.CeedElemRestrictionGetMultiplicity(
128            self._pointer[0], mult._pointer[0])
129        self._ceed._check_error(err_code)
130
131        # Return
132        return mult
133
134    # Get ElemRestrition Layout
135    def get_layout(self):
136        """Get the element vector layout of an ElemRestriction.
137
138           Returns:
139             layout: Vector containing layout array, stored as [nodes, components, elements].
140                     The data for node i, component j, element k in the element
141                     vector is given by i*layout[0] + j*layout[1] + k*layout[2]."""
142
143        # Create output array
144        layout = np.zeros(3, dtype="int32")
145        array_pointer = ffi.cast(
146            "CeedInt *",
147            layout.__array_interface__['data'][0])
148
149        # libCEED call
150        err_code = lib.CeedElemRestrictionGetELayout(
151            self._pointer[0], array_pointer)
152        self._ceed._check_error(err_code)
153
154        # Return
155        return layout
156
157# ------------------------------------------------------------------------------
158
159
160class ElemRestriction(_ElemRestrictionBase):
161    """Ceed ElemRestriction: restriction from local vectors to elements."""
162
163    # Constructor
164    def __init__(self, ceed, nelem, elemsize, ncomp, compstride, lsize, offsets,
165                 memtype=MEM_HOST, cmode=COPY_VALUES):
166        # CeedVector object
167        self._pointer = ffi.new("CeedElemRestriction *")
168
169        # Reference to Ceed
170        self._ceed = ceed
171
172        # Store array reference if needed
173        if cmode == USE_POINTER:
174            self._array_reference = offsets
175        else:
176            self._array_reference = None
177
178        # Setup the numpy array for the libCEED call
179        offsets_pointer = ffi.new("const CeedInt *")
180        offsets_pointer = ffi.cast("const CeedInt *",
181                                   offsets.__array_interface__['data'][0])
182
183        # libCEED call
184        err_code = lib.CeedElemRestrictionCreate(self._ceed._pointer[0], nelem,
185                                                 elemsize, ncomp, compstride,
186                                                 lsize, memtype, cmode,
187                                                 offsets_pointer, self._pointer)
188        self._ceed._check_error(err_code)
189
190# ------------------------------------------------------------------------------
191
192
193class StridedElemRestriction(_ElemRestrictionBase):
194    """Ceed Strided ElemRestriction: strided restriction from local vectors to elements."""
195
196    # Constructor
197    def __init__(self, ceed, nelem, elemsize, ncomp, lsize, strides):
198        # CeedVector object
199        self._pointer = ffi.new("CeedElemRestriction *")
200
201        # Reference to Ceed
202        self._ceed = ceed
203
204        # Setup the numpy array for the libCEED call
205        strides_pointer = ffi.new("const CeedInt *")
206        strides_pointer = ffi.cast("const CeedInt *",
207                                   strides.__array_interface__['data'][0])
208
209        # libCEED call
210        err_code = lib.CeedElemRestrictionCreateStrided(self._ceed._pointer[0],
211                                                        nelem, elemsize, ncomp,
212                                                        lsize, strides_pointer,
213                                                        self._pointer)
214        self._ceed._check_error(err_code)
215
216# ------------------------------------------------------------------------------
217
218
219class BlockedElemRestriction(_ElemRestrictionBase):
220    """Ceed Blocked ElemRestriction: blocked restriction from local vectors to elements."""
221
222    # Constructor
223    def __init__(self, ceed, nelem, elemsize, blksize, ncomp, compstride, lsize,
224                 offsets, memtype=MEM_HOST, cmode=COPY_VALUES):
225        # CeedVector object
226        self._pointer = ffi.new("CeedElemRestriction *")
227
228        # Reference to Ceed
229        self._ceed = ceed
230
231        # Setup the numpy array for the libCEED call
232        offsets_pointer = ffi.new("const CeedInt *")
233        offsets_pointer = ffi.cast("const CeedInt *",
234                                   offsets.__array_interface__['data'][0])
235
236        # libCEED call
237        err_code = lib.CeedElemRestrictionCreateBlocked(self._ceed._pointer[0], nelem,
238                                                        elemsize, blksize, ncomp,
239                                                        compstride, lsize, memtype, cmode,
240                                                        offsets_pointer, self._pointer)
241        self._ceed._check_error(err_code)
242
243    # Transpose a Blocked ElemRestriction
244    @property
245    def T(self):
246        """Transpose a BlockedElemRestriction."""
247
248        return TransposeBlockedElemRestriction(self)
249
250    # Transpose a Blocked ElemRestriction
251    @property
252    def transpose(self):
253        """Transpose a BlockedElemRestriction."""
254
255        return TransposeBlockedElemRestriction(self)
256
257    # Apply CeedElemRestriction to single block
258    def apply_block(self, block, u, v, tmode=NOTRANSPOSE,
259                    request=REQUEST_IMMEDIATE):
260        """Restrict a local vector to a block of an element vector or apply its transpose.
261
262           Args:
263             block: block number to restrict to/from, i.e. block=0 will handle
264                      elements [0 : blksize] and block=3 will handle elements
265                      [3*blksize : 4*blksize]
266             u: input vector
267             v: output vector
268             **tmode: apply restriction or transpose, default CEED_NOTRANSPOSE
269             **request: Ceed request, default CEED_REQUEST_IMMEDIATE"""
270
271        # libCEED call
272        err_code = lib.CeedElemRestrictionApplyBlock(self._pointer[0], block, tmode,
273                                                     u._pointer[0], v._pointer[0],
274                                                     request)
275        self._ceed._check_error(err_code)
276
277# ------------------------------------------------------------------------------
278
279
280class BlockedStridedElemRestriction(BlockedElemRestriction):
281    """Ceed Blocked Strided ElemRestriction: strided restriction from local vectors to elements."""
282
283    # Constructor
284    def __init__(self, ceed, nelem, elemsize, blksize, ncomp, lsize, strides):
285        # CeedVector object
286        self._pointer = ffi.new("CeedElemRestriction *")
287
288        # Reference to Ceed
289        self._ceed = ceed
290
291        # Setup the numpy array for the libCEED call
292        strides_pointer = ffi.new("const CeedInt *")
293        strides_pointer = ffi.cast("const CeedInt *",
294                                   strides.__array_interface__['data'][0])
295
296        # libCEED call
297        err_code = lib.CeedElemRestrictionCreateBlockedStrided(self._ceed._pointer[0], nelem,
298                                                               elemsize, blksize, ncomp,
299                                                               lsize, strides_pointer,
300                                                               self._pointer)
301        self._ceed._check_error(err_code)
302
303# ------------------------------------------------------------------------------
304
305
306class TransposeElemRestriction():
307    """Ceed ElemRestriction: transpose restriction from elements to local vectors."""
308
309    # Attributes
310    _elemrestriction = None
311
312    # Constructor
313    def __init__(self, elemrestriction):
314
315        # Reference elemrestriction
316        self._elemrestriction = elemrestriction
317
318    # Representation
319    def __repr__(self):
320        return "<Transpose CeedElemRestriction instance at " + \
321            hex(id(self)) + ">"
322
323    # Apply Transpose CeedElemRestriction
324
325    def apply(self, u, v, request=REQUEST_IMMEDIATE):
326        """Restrict an element vector to a local vector.
327
328           Args:
329             u: input vector
330             v: output vector
331             **request: Ceed request, default CEED_REQUEST_IMMEDIATE"""
332
333        # libCEED call
334        self._elemrestriction.apply(u, v, request=request, tmode=TRANSPOSE)
335
336# ------------------------------------------------------------------------------
337
338
339class TransposeBlockedElemRestriction(TransposeElemRestriction):
340    """Transpose Ceed Blocked ElemRestriction: blocked transpose restriction from elements
341         to local vectors."""
342
343    # Apply Transpose CeedElemRestriction
344    def apply_block(self, block, u, v, request=REQUEST_IMMEDIATE):
345        """Restrict a block of an element vector to a local vector.
346
347           Args:
348             block: block number to restrict to/from, i.e. block=0 will handle
349                      elements [0 : blksize] and block=3 will handle elements
350                      [3*blksize : 4*blksize]
351             u: input vector
352             v: output vector
353             **request: Ceed request, default CEED_REQUEST_IMMEDIATE"""
354
355        # libCEED call
356        self._elemrestriction.apply_block(block, u, v, request=request,
357                                          tmode=TRANSPOSE)
358
359# ------------------------------------------------------------------------------
360