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