1c8e769f6Sjeremylt# Copyright (c) 2017, Lawrence Livermore National Security, LLC. Produced at 2c8e769f6Sjeremylt# the Lawrence Livermore National Laboratory. LLNL-CODE-734707. All Rights 3c8e769f6Sjeremylt# reserved. See files LICENSE and NOTICE for details. 4c8e769f6Sjeremylt# 5c8e769f6Sjeremylt# This file is part of CEED, a collection of benchmarks, miniapps, software 6c8e769f6Sjeremylt# libraries and APIs for efficient high-order finite element and spectral 7c8e769f6Sjeremylt# element discretizations for exascale applications. For more information and 8c8e769f6Sjeremylt# source code availability see http://github.com/ceed. 9c8e769f6Sjeremylt# 10c8e769f6Sjeremylt# The CEED research is supported by the Exascale Computing Project 17-SC-20-SC, 11c8e769f6Sjeremylt# a collaborative effort of two U.S. Department of Energy organizations (Office 12c8e769f6Sjeremylt# of Science and the National Nuclear Security Administration) responsible for 13c8e769f6Sjeremylt# the planning and preparation of a capable exascale ecosystem, including 14c8e769f6Sjeremylt# software, applications, hardware, advanced system engineering and early 15c8e769f6Sjeremylt# testbed platforms, in support of the nation's exascale computing imperative. 16c8e769f6Sjeremylt 17c8e769f6Sjeremyltfrom _ceed_cffi import ffi, lib 18c8e769f6Sjeremyltimport tempfile 19c8e769f6Sjeremyltimport numpy as np 20c8e769f6Sjeremyltfrom abc import ABC 2161dbc9d2Sjeremyltfrom .ceed_constants import REQUEST_IMMEDIATE, REQUEST_ORDERED, MEM_HOST, COPY_VALUES, TRANSPOSE, NOTRANSPOSE, INTERLACED, NONINTERLACED 22c8e769f6Sjeremyltfrom .ceed_vector import _VectorWrap 23c8e769f6Sjeremylt 24c8e769f6Sjeremylt# ------------------------------------------------------------------------------ 25c8e769f6Sjeremyltclass _ElemRestrictionBase(ABC): 26c8e769f6Sjeremylt """Ceed ElemRestriction: restriction from local vectors to elements.""" 27c8e769f6Sjeremylt 28c8e769f6Sjeremylt # Attributes 29c8e769f6Sjeremylt _ceed = ffi.NULL 30c8e769f6Sjeremylt _pointer = ffi.NULL 31c8e769f6Sjeremylt 32c8e769f6Sjeremylt # Destructor 33c8e769f6Sjeremylt def __del__(self): 34c8e769f6Sjeremylt # libCEED call 35c8e769f6Sjeremylt lib.CeedElemRestrictionDestroy(self._pointer) 36c8e769f6Sjeremylt 37c8e769f6Sjeremylt # Representation 38c8e769f6Sjeremylt def __repr__(self): 39c8e769f6Sjeremylt return "<CeedElemRestriction instance at " + hex(id(self)) + ">" 40c8e769f6Sjeremylt 41c8e769f6Sjeremylt # String conversion for print() to stdout 42c8e769f6Sjeremylt def __str__(self): 43c8e769f6Sjeremylt """View an ElemRestriction via print().""" 44c8e769f6Sjeremylt 45c8e769f6Sjeremylt # libCEED call 46c8e769f6Sjeremylt with tempfile.NamedTemporaryFile() as key_file: 47c8e769f6Sjeremylt with open(key_file.name, 'r+') as stream_file: 48c8e769f6Sjeremylt stream = ffi.cast("FILE *", stream_file) 49c8e769f6Sjeremylt 50c8e769f6Sjeremylt lib.CeedElemRestrictionView(self._pointer[0], stream) 51c8e769f6Sjeremylt 52c8e769f6Sjeremylt stream_file.seek(0) 53c8e769f6Sjeremylt out_string = stream_file.read() 54c8e769f6Sjeremylt 55c8e769f6Sjeremylt return out_string 56c8e769f6Sjeremylt 57c8e769f6Sjeremylt # Apply CeedElemRestriction 58a8d32208Sjeremylt def apply(self, u, v, tmode=NOTRANSPOSE, request=REQUEST_IMMEDIATE): 59c8e769f6Sjeremylt """Restrict a local vector to an element vector or apply its transpose. 60c8e769f6Sjeremylt 61c8e769f6Sjeremylt Args: 62c8e769f6Sjeremylt u: input vector 63c8e769f6Sjeremylt v: output vector 64c8e769f6Sjeremylt **tmode: apply restriction or transpose, default CEED_NOTRANSPOSE 65c8e769f6Sjeremylt **request: Ceed request, default CEED_REQUEST_IMMEDIATE""" 66c8e769f6Sjeremylt 67c8e769f6Sjeremylt # libCEED call 68a8d32208Sjeremylt lib.CeedElemRestrictionApply(self._pointer[0], tmode, u._pointer[0], 69a8d32208Sjeremylt v._pointer[0], request) 70c8e769f6Sjeremylt 71c8e769f6Sjeremylt # Transpose an ElemRestriction 72c8e769f6Sjeremylt @property 73c8e769f6Sjeremylt def T(self): 74c8e769f6Sjeremylt """Transpose an ElemRestriction.""" 75c8e769f6Sjeremylt 76c8e769f6Sjeremylt return TransposeElemRestriction(self) 77c8e769f6Sjeremylt 78c8e769f6Sjeremylt # Transpose an ElemRestriction 79c8e769f6Sjeremylt @property 80c8e769f6Sjeremylt def transpose(self): 81c8e769f6Sjeremylt """Transpose an ElemRestriction.""" 82c8e769f6Sjeremylt 83c8e769f6Sjeremylt return TransposeElemRestriction(self) 84c8e769f6Sjeremylt 85c8e769f6Sjeremylt # Create restriction vectors 86c8e769f6Sjeremylt def create_vector(self, createLvec = True, createEvec = True): 87c8e769f6Sjeremylt """Create Vectors associated with an ElemRestriction. 88c8e769f6Sjeremylt 89c8e769f6Sjeremylt Args: 90c8e769f6Sjeremylt **createLvec: flag to create local vector, default True 91c8e769f6Sjeremylt **createEvec: flag to create element vector, default True 92c8e769f6Sjeremylt 93c8e769f6Sjeremylt Returns: 94c8e769f6Sjeremylt [lvec, evec]: local vector and element vector, or None if flag set to false""" 95c8e769f6Sjeremylt 96c8e769f6Sjeremylt # Vector pointers 97c8e769f6Sjeremylt lvecPointer = ffi.new("CeedVector *") if createLvec else ffi.NULL 98c8e769f6Sjeremylt evecPointer = ffi.new("CeedVector *") if createEvec else ffi.NULL 99c8e769f6Sjeremylt 100c8e769f6Sjeremylt # libCEED call 101c8e769f6Sjeremylt lib.CeedElemRestrictionCreateVector(self._pointer[0], lvecPointer, 102c8e769f6Sjeremylt evecPointer) 103c8e769f6Sjeremylt 104c8e769f6Sjeremylt # Return vectors 105c8e769f6Sjeremylt lvec = _VectorWrap(self._ceed._pointer, lvecPointer) if createLvec else None 106c8e769f6Sjeremylt evec = _VectorWrap(self._ceed._pointer, evecPointer) if createEvec else None 107c8e769f6Sjeremylt 108c8e769f6Sjeremylt # Return 109c8e769f6Sjeremylt return [lvec, evec] 110c8e769f6Sjeremylt 111c8e769f6Sjeremylt # Get ElemRestriction multiplicity 112a8d32208Sjeremylt def get_multiplicity(self): 113c8e769f6Sjeremylt """Get the multiplicity of nodes in an ElemRestriction. 114c8e769f6Sjeremylt 115c8e769f6Sjeremylt Returns: 116c8e769f6Sjeremylt mult: local vector containing multiplicity of nodes in ElemRestriction""" 117c8e769f6Sjeremylt 118c8e769f6Sjeremylt # Create mult vector 119c8e769f6Sjeremylt [mult, evec] = self.create_vector(createEvec = False) 120c8e769f6Sjeremylt mult.set_value(0) 121c8e769f6Sjeremylt 122c8e769f6Sjeremylt # libCEED call 123a8d32208Sjeremylt lib.CeedElemRestrictionGetMultiplicity(self._pointer[0], mult._pointer[0]) 124c8e769f6Sjeremylt 125c8e769f6Sjeremylt # Return 126c8e769f6Sjeremylt return mult 127c8e769f6Sjeremylt 128c8e769f6Sjeremylt# ------------------------------------------------------------------------------ 129c8e769f6Sjeremyltclass ElemRestriction(_ElemRestrictionBase): 130c8e769f6Sjeremylt """Ceed ElemRestriction: restriction from local vectors to elements.""" 131c8e769f6Sjeremylt 132c8e769f6Sjeremylt # Constructor 133c8e769f6Sjeremylt def __init__(self, ceed, nelem, elemsize, nnodes, ncomp, indices, 13461dbc9d2Sjeremylt memtype=MEM_HOST, cmode=COPY_VALUES, imode=NONINTERLACED): 135c8e769f6Sjeremylt # CeedVector object 136c8e769f6Sjeremylt self._pointer = ffi.new("CeedElemRestriction *") 137c8e769f6Sjeremylt 138c8e769f6Sjeremylt # Reference to Ceed 139c8e769f6Sjeremylt self._ceed = ceed 140c8e769f6Sjeremylt 141c8e769f6Sjeremylt # Setup the numpy array for the libCEED call 142c8e769f6Sjeremylt indices_pointer = ffi.new("const CeedInt *") 143c8e769f6Sjeremylt indices_pointer = ffi.cast("const CeedInt *", 144c8e769f6Sjeremylt indices.__array_interface__['data'][0]) 145c8e769f6Sjeremylt 146c8e769f6Sjeremylt # libCEED call 14761dbc9d2Sjeremylt lib.CeedElemRestrictionCreate(self._ceed._pointer[0], imode, nelem, 148a8d32208Sjeremylt elemsize, nnodes, ncomp, memtype, cmode, 149c8e769f6Sjeremylt indices_pointer, self._pointer) 150c8e769f6Sjeremylt 151c8e769f6Sjeremylt# ------------------------------------------------------------------------------ 152*69a53589Sjeremyltclass StridedElemRestriction(_ElemRestrictionBase): 153*69a53589Sjeremylt """Ceed Strided ElemRestriction: strided restriction from local vectors to elements.""" 154c8e769f6Sjeremylt 155c8e769f6Sjeremylt # Constructor 156*69a53589Sjeremylt def __init__(self, ceed, nelem, elemsize, nnodes, ncomp, strides): 157c8e769f6Sjeremylt # CeedVector object 158c8e769f6Sjeremylt self._pointer = ffi.new("CeedElemRestriction *") 159c8e769f6Sjeremylt 160c8e769f6Sjeremylt # Reference to Ceed 161c8e769f6Sjeremylt self._ceed = ceed 162c8e769f6Sjeremylt 163*69a53589Sjeremylt # Setup the numpy array for the libCEED call 164*69a53589Sjeremylt strides_pointer = ffi.new("const CeedInt *") 165*69a53589Sjeremylt strides_pointer = ffi.cast("const CeedInt *", 166*69a53589Sjeremylt strides.__array_interface__['data'][0]) 167*69a53589Sjeremylt 168c8e769f6Sjeremylt # libCEED call 169*69a53589Sjeremylt lib.CeedElemRestrictionCreateStrided(self._ceed._pointer[0], nelem, 170c8e769f6Sjeremylt elemsize, nnodes, ncomp, 171*69a53589Sjeremylt strides_pointer, self._pointer) 172c8e769f6Sjeremylt 173c8e769f6Sjeremylt# ------------------------------------------------------------------------------ 174c8e769f6Sjeremyltclass BlockedElemRestriction(_ElemRestrictionBase): 175c8e769f6Sjeremylt """Ceed Blocked ElemRestriction: blocked restriction from local vectors to elements.""" 176c8e769f6Sjeremylt 177c8e769f6Sjeremylt # Constructor 178c8e769f6Sjeremylt def __init__(self, ceed, nelem, elemsize, blksize, nnodes, ncomp, indices, 17961dbc9d2Sjeremylt memtype=MEM_HOST, cmode=COPY_VALUES, imode=NONINTERLACED): 180c8e769f6Sjeremylt # CeedVector object 181c8e769f6Sjeremylt self._pointer = ffi.new("CeedElemRestriction *") 182c8e769f6Sjeremylt 183c8e769f6Sjeremylt # Reference to Ceed 184c8e769f6Sjeremylt self._ceed = ceed 185c8e769f6Sjeremylt 186c8e769f6Sjeremylt # Setup the numpy array for the libCEED call 187c8e769f6Sjeremylt indices_pointer = ffi.new("const CeedInt *") 188c8e769f6Sjeremylt indices_pointer = ffi.cast("const CeedInt *", 189c8e769f6Sjeremylt indices.__array_interface__['data'][0]) 190c8e769f6Sjeremylt 191c8e769f6Sjeremylt # libCEED call 19261dbc9d2Sjeremylt lib.CeedElemRestrictionCreateBlocked(self._ceed._pointer[0], imode, nelem, 193c8e769f6Sjeremylt elemsize, blksize, nnodes, ncomp, 194c8e769f6Sjeremylt memtype, cmode, indices_pointer, 195c8e769f6Sjeremylt self._pointer) 196c8e769f6Sjeremylt 197c8e769f6Sjeremylt # Transpose a Blocked ElemRestriction 198c8e769f6Sjeremylt @property 199c8e769f6Sjeremylt def T(self): 200c8e769f6Sjeremylt """Transpose a BlockedElemRestriction.""" 201c8e769f6Sjeremylt 202c8e769f6Sjeremylt return TransposeBlockedElemRestriction(self) 203c8e769f6Sjeremylt 204c8e769f6Sjeremylt # Transpose a Blocked ElemRestriction 205c8e769f6Sjeremylt @property 206c8e769f6Sjeremylt def transpose(self): 207c8e769f6Sjeremylt """Transpose a BlockedElemRestriction.""" 208c8e769f6Sjeremylt 209c8e769f6Sjeremylt return TransposeBlockedElemRestriction(self) 210c8e769f6Sjeremylt 211c8e769f6Sjeremylt # Apply CeedElemRestriction to single block 212a8d32208Sjeremylt def apply_block(self, block, u, v, tmode=NOTRANSPOSE, 213c8e769f6Sjeremylt request=REQUEST_IMMEDIATE): 214c8e769f6Sjeremylt """Restrict a local vector to a block of an element vector or apply its transpose. 215c8e769f6Sjeremylt 216c8e769f6Sjeremylt Args: 217c8e769f6Sjeremylt block: block number to restrict to/from, i.e. block=0 will handle 218c8e769f6Sjeremylt elements [0 : blksize] and block=3 will handle elements 219c8e769f6Sjeremylt [3*blksize : 4*blksize] 220c8e769f6Sjeremylt u: input vector 221c8e769f6Sjeremylt v: output vector 222c8e769f6Sjeremylt **tmode: apply restriction or transpose, default CEED_NOTRANSPOSE 223c8e769f6Sjeremylt **request: Ceed request, default CEED_REQUEST_IMMEDIATE""" 224c8e769f6Sjeremylt 225c8e769f6Sjeremylt # libCEED call 226c8e769f6Sjeremylt lib.CeedElemRestrictionApplyBlock(self._pointer[0], block, tmode, 227a8d32208Sjeremylt u._pointer[0], v._pointer[0], request) 228c8e769f6Sjeremylt 229c8e769f6Sjeremylt# ------------------------------------------------------------------------------ 230*69a53589Sjeremyltclass BlockedStridedElemRestriction(BlockedElemRestriction): 231*69a53589Sjeremylt """Ceed Blocked Strided ElemRestriction: strided restriction from local vectors to elements.""" 232*69a53589Sjeremylt 233*69a53589Sjeremylt # Constructor 234*69a53589Sjeremylt def __init__(self, ceed, nelem, elemsize, blksize, nnodes, ncomp, strides): 235*69a53589Sjeremylt # CeedVector object 236*69a53589Sjeremylt self._pointer = ffi.new("CeedElemRestriction *") 237*69a53589Sjeremylt 238*69a53589Sjeremylt # Reference to Ceed 239*69a53589Sjeremylt self._ceed = ceed 240*69a53589Sjeremylt 241*69a53589Sjeremylt # Setup the numpy array for the libCEED call 242*69a53589Sjeremylt strides_pointer = ffi.new("const CeedInt *") 243*69a53589Sjeremylt strides_pointer = ffi.cast("const CeedInt *", 244*69a53589Sjeremylt strides.__array_interface__['data'][0]) 245*69a53589Sjeremylt 246*69a53589Sjeremylt # libCEED call 247*69a53589Sjeremylt lib.CeedElemRestrictionCreateBlockedStrided(self._ceed._pointer[0], nelem, 248*69a53589Sjeremylt elemsize, blksize, nnodes, 249*69a53589Sjeremylt ncomp, strides_pointer, 250*69a53589Sjeremylt self._pointer) 251*69a53589Sjeremylt 252*69a53589Sjeremylt# ------------------------------------------------------------------------------ 253c8e769f6Sjeremyltclass TransposeElemRestriction(): 254c8e769f6Sjeremylt """Ceed ElemRestriction: transpose restriction from elements to local vectors.""" 255c8e769f6Sjeremylt 256c8e769f6Sjeremylt # Attributes 257c8e769f6Sjeremylt _elemrestriction = None 258c8e769f6Sjeremylt 259c8e769f6Sjeremylt # Constructor 260c8e769f6Sjeremylt def __init__(self, elemrestriction): 261c8e769f6Sjeremylt 262c8e769f6Sjeremylt # Reference elemrestriction 263c8e769f6Sjeremylt self._elemrestriction = elemrestriction 264c8e769f6Sjeremylt 265c8e769f6Sjeremylt # Representation 266c8e769f6Sjeremylt def __repr__(self): 267c8e769f6Sjeremylt return "<Transpose CeedElemRestriction instance at " + hex(id(self)) + ">" 268c8e769f6Sjeremylt 269c8e769f6Sjeremylt 270c8e769f6Sjeremylt # Apply Transpose CeedElemRestriction 271a8d32208Sjeremylt def apply(self, u, v, request=REQUEST_IMMEDIATE): 272c8e769f6Sjeremylt """Restrict an element vector to a local vector. 273c8e769f6Sjeremylt 274c8e769f6Sjeremylt Args: 275c8e769f6Sjeremylt u: input vector 276c8e769f6Sjeremylt v: output vector 277c8e769f6Sjeremylt **request: Ceed request, default CEED_REQUEST_IMMEDIATE""" 278c8e769f6Sjeremylt 279c8e769f6Sjeremylt # libCEED call 280a8d32208Sjeremylt self._elemrestriction.apply(u, v, request=request, tmode=TRANSPOSE) 281c8e769f6Sjeremylt 282c8e769f6Sjeremylt# ------------------------------------------------------------------------------ 283c8e769f6Sjeremyltclass TransposeBlockedElemRestriction(TransposeElemRestriction): 284c8e769f6Sjeremylt """Transpose Ceed Blocked ElemRestriction: blocked transpose restriction from elements 285c8e769f6Sjeremylt to local vectors.""" 286c8e769f6Sjeremylt 287c8e769f6Sjeremylt # Apply Transpose CeedElemRestriction 288a8d32208Sjeremylt def apply_block(self, block, u, v, request=REQUEST_IMMEDIATE): 289c8e769f6Sjeremylt """Restrict a block of an element vector to a local vector. 290c8e769f6Sjeremylt 291c8e769f6Sjeremylt Args: 292c8e769f6Sjeremylt block: block number to restrict to/from, i.e. block=0 will handle 293c8e769f6Sjeremylt elements [0 : blksize] and block=3 will handle elements 294c8e769f6Sjeremylt [3*blksize : 4*blksize] 295c8e769f6Sjeremylt u: input vector 296c8e769f6Sjeremylt v: output vector 297c8e769f6Sjeremylt **request: Ceed request, default CEED_REQUEST_IMMEDIATE""" 298c8e769f6Sjeremylt 299c8e769f6Sjeremylt # libCEED call 300a8d32208Sjeremylt self._elemrestriction.apply_block(block, u, v, request=request, 301c8e769f6Sjeremylt tmode=TRANSPOSE) 302c8e769f6Sjeremylt 303c8e769f6Sjeremylt# ------------------------------------------------------------------------------ 304