1f2d2bf5dSjeremylt# Copyright (c) 2017, Lawrence Livermore National Security, LLC. Produced at 2f2d2bf5dSjeremylt# the Lawrence Livermore National Laboratory. LLNL-CODE-734707. All Rights 3f2d2bf5dSjeremylt# reserved. See files LICENSE and NOTICE for details. 4f2d2bf5dSjeremylt# 5f2d2bf5dSjeremylt# This file is part of CEED, a collection of benchmarks, miniapps, software 6f2d2bf5dSjeremylt# libraries and APIs for efficient high-order finite element and spectral 7f2d2bf5dSjeremylt# element discretizations for exascale applications. For more information and 8f2d2bf5dSjeremylt# source code availability see http://github.com/ceed. 9f2d2bf5dSjeremylt# 10f2d2bf5dSjeremylt# The CEED research is supported by the Exascale Computing Project 17-SC-20-SC, 11f2d2bf5dSjeremylt# a collaborative effort of two U.S. Department of Energy organizations (Office 12f2d2bf5dSjeremylt# of Science and the National Nuclear Security Administration) responsible for 13f2d2bf5dSjeremylt# the planning and preparation of a capable exascale ecosystem, including 14f2d2bf5dSjeremylt# software, applications, hardware, advanced system engineering and early 15f2d2bf5dSjeremylt# testbed platforms, in support of the nation's exascale computing imperative. 16f2d2bf5dSjeremylt 17f2d2bf5dSjeremyltfrom _ceed_cffi import ffi, lib 18f2d2bf5dSjeremyltimport tempfile 19f2d2bf5dSjeremyltfrom abc import ABC 20f2d2bf5dSjeremyltfrom .ceed_constants import REQUEST_IMMEDIATE, REQUEST_ORDERED, NOTRANSPOSE 21f2d2bf5dSjeremylt 22f2d2bf5dSjeremylt# ------------------------------------------------------------------------------ 237a7b0fa3SJed Brown 247a7b0fa3SJed Brown 25f2d2bf5dSjeremyltclass _OperatorBase(ABC): 26f2d2bf5dSjeremylt """Ceed Operator: composed FE-type operations on vectors.""" 27f2d2bf5dSjeremylt 28f2d2bf5dSjeremylt # Attributes 29f2d2bf5dSjeremylt _ceed = ffi.NULL 30f2d2bf5dSjeremylt _pointer = ffi.NULL 31f2d2bf5dSjeremylt 32f2d2bf5dSjeremylt # Destructor 33f2d2bf5dSjeremylt def __del__(self): 34f2d2bf5dSjeremylt # libCEED call 35f2d2bf5dSjeremylt lib.CeedOperatorDestroy(self._pointer) 36f2d2bf5dSjeremylt 37f2d2bf5dSjeremylt # Representation 38f2d2bf5dSjeremylt def __repr__(self): 39f2d2bf5dSjeremylt return "<CeedOperator instance at " + hex(id(self)) + ">" 40f2d2bf5dSjeremylt 41f2d2bf5dSjeremylt # String conversion for print() to stdout 42f2d2bf5dSjeremylt def __str__(self): 43f2d2bf5dSjeremylt """View an Operator via print().""" 44f2d2bf5dSjeremylt 45f2d2bf5dSjeremylt # libCEED call 46f2d2bf5dSjeremylt with tempfile.NamedTemporaryFile() as key_file: 47f2d2bf5dSjeremylt with open(key_file.name, 'r+') as stream_file: 48f2d2bf5dSjeremylt stream = ffi.cast("FILE *", stream_file) 49f2d2bf5dSjeremylt 50f2d2bf5dSjeremylt lib.CeedOperatorView(self._pointer[0], stream) 51f2d2bf5dSjeremylt 52f2d2bf5dSjeremylt stream_file.seek(0) 53f2d2bf5dSjeremylt out_string = stream_file.read() 54f2d2bf5dSjeremylt 55f2d2bf5dSjeremylt return out_string 56f2d2bf5dSjeremylt 57*6397b31bSJeremy L Thompson # Assemble linear diagonal 58*6397b31bSJeremy L Thompson def linear_assemble_diagonal(self, d, request=REQUEST_IMMEDIATE): 59*6397b31bSJeremy L Thompson """Assemble the diagonal of a square linear Operator 60*6397b31bSJeremy L Thompson 61*6397b31bSJeremy L Thompson Args: 62*6397b31bSJeremy L Thompson d: Vector to store assembled Operator diagonal 63*6397b31bSJeremy L Thompson **request: Ceed request, default CEED_REQUEST_IMMEDIATE""" 64*6397b31bSJeremy L Thompson 65*6397b31bSJeremy L Thompson # libCEED call 66*6397b31bSJeremy L Thompson lib.CeedOperatorLinearAssembleDiagonal(self._pointer[0], 67*6397b31bSJeremy L Thompson d._pointer[0], request) 68*6397b31bSJeremy L Thompson 69*6397b31bSJeremy L Thompson # Assemble add linear diagonal 70*6397b31bSJeremy L Thompson def linear_assemble_add_diagonal(self, d, request=REQUEST_IMMEDIATE): 71*6397b31bSJeremy L Thompson """Sum the diagonal of a square linear Operator into a Vector 72*6397b31bSJeremy L Thompson 73*6397b31bSJeremy L Thompson Args: 74*6397b31bSJeremy L Thompson d: Vector to store assembled Operator diagonal 75*6397b31bSJeremy L Thompson **request: Ceed request, default CEED_REQUEST_IMMEDIATE""" 76*6397b31bSJeremy L Thompson 77*6397b31bSJeremy L Thompson # libCEED call 78*6397b31bSJeremy L Thompson lib.CeedOperatorLinearAssembleAddDiagonal(self._pointer[0], 79*6397b31bSJeremy L Thompson d._pointer[0], request) 80*6397b31bSJeremy L Thompson 81*6397b31bSJeremy L Thompson # Assemble linear point block diagonal 82*6397b31bSJeremy L Thompson def linear_assemble_point_block_diagonal(self, d, request=REQUEST_IMMEDIATE): 83*6397b31bSJeremy L Thompson """Assemble the diagonal of a square linear Operator 84*6397b31bSJeremy L Thompson 85*6397b31bSJeremy L Thompson Args: 86*6397b31bSJeremy L Thompson d: Vector to store assembled Operator point block diagonal, 87*6397b31bSJeremy L Thompson provided in row-major form with an ncomp*ncomp block 88*6397b31bSJeremy L Thompson at each node 89*6397b31bSJeremy L Thompson **request: Ceed request, default CEED_REQUEST_IMMEDIATE""" 90*6397b31bSJeremy L Thompson 91*6397b31bSJeremy L Thompson # libCEED call 92*6397b31bSJeremy L Thompson lib.CeedOperatorLinearAssemblePointBlockDiagonal(self._pointer[0], 93*6397b31bSJeremy L Thompson d._pointer[0], request) 94*6397b31bSJeremy L Thompson 95*6397b31bSJeremy L Thompson # Assemble linear point block diagonal 96*6397b31bSJeremy L Thompson def linear_assemble_add_point_block_diagonal(self, d, request=REQUEST_IMMEDIATE): 97*6397b31bSJeremy L Thompson """Sum the diagonal of a square linear Operator into a Vector 98*6397b31bSJeremy L Thompson 99*6397b31bSJeremy L Thompson Args: 100*6397b31bSJeremy L Thompson d: Vector to store assembled Operator point block diagonal, 101*6397b31bSJeremy L Thompson provided in row-major form with an ncomp*ncomp block 102*6397b31bSJeremy L Thompson at each node 103*6397b31bSJeremy L Thompson **request: Ceed request, default CEED_REQUEST_IMMEDIATE""" 104*6397b31bSJeremy L Thompson 105*6397b31bSJeremy L Thompson # libCEED call 106*6397b31bSJeremy L Thompson lib.CeedOperatorLinearAssembleAddPointBlockDiagonal(self._pointer[0], 107*6397b31bSJeremy L Thompson d._pointer[0], request) 108*6397b31bSJeremy L Thompson 109f2d2bf5dSjeremylt # Apply CeedOperator 110f2d2bf5dSjeremylt def apply(self, u, v, request=REQUEST_IMMEDIATE): 111f2d2bf5dSjeremylt """Apply Operator to a vector. 112f2d2bf5dSjeremylt 113f2d2bf5dSjeremylt Args: 114f2d2bf5dSjeremylt u: Vector containing input state or CEED_VECTOR_NONE if there are no 115f2d2bf5dSjeremylt active inputs 116f2d2bf5dSjeremylt v: Vector to store result of applying operator (must be distinct from u) 117f2d2bf5dSjeremylt or CEED_VECTOR_NONE if there are no active outputs 118f2d2bf5dSjeremylt **request: Ceed request, default CEED_REQUEST_IMMEDIATE""" 119f2d2bf5dSjeremylt 120f2d2bf5dSjeremylt # libCEED call 121f2d2bf5dSjeremylt lib.CeedOperatorApply(self._pointer[0], u._pointer[0], v._pointer[0], 122f2d2bf5dSjeremylt request) 123f2d2bf5dSjeremylt 124f2d2bf5dSjeremylt # Apply CeedOperator 125f2d2bf5dSjeremylt def apply_add(self, u, v, request=REQUEST_IMMEDIATE): 126f2d2bf5dSjeremylt """Apply Operator to a vector and add result to output vector. 127f2d2bf5dSjeremylt 128f2d2bf5dSjeremylt Args: 129f2d2bf5dSjeremylt u: Vector containing input state or CEED_VECTOR_NONE if there are no 130f2d2bf5dSjeremylt active inputs 131f2d2bf5dSjeremylt v: Vector to sum in result of applying operator (must be distinct from u) 132f2d2bf5dSjeremylt or CEED_VECTOR_NONE if there are no active outputs 133f2d2bf5dSjeremylt **request: Ceed request, default CEED_REQUEST_IMMEDIATE""" 134f2d2bf5dSjeremylt 135f2d2bf5dSjeremylt # libCEED call 136f2d2bf5dSjeremylt lib.CeedOperatorApplyAdd(self._pointer[0], u._pointer[0], v._pointer[0], 137f2d2bf5dSjeremylt request) 138f2d2bf5dSjeremylt 139f2d2bf5dSjeremylt# ------------------------------------------------------------------------------ 1407a7b0fa3SJed Brown 1417a7b0fa3SJed Brown 142f2d2bf5dSjeremyltclass Operator(_OperatorBase): 143f2d2bf5dSjeremylt """Ceed Operator: composed FE-type operations on vectors.""" 144f2d2bf5dSjeremylt 145f2d2bf5dSjeremylt # Constructor 146f2d2bf5dSjeremylt def __init__(self, ceed, qf, dqf=None, dqfT=None): 147f2d2bf5dSjeremylt # CeedOperator object 148f2d2bf5dSjeremylt self._pointer = ffi.new("CeedOperator *") 149f2d2bf5dSjeremylt 150f2d2bf5dSjeremylt # Reference to Ceed 151f2d2bf5dSjeremylt self._ceed = ceed 152f2d2bf5dSjeremylt 153f2d2bf5dSjeremylt # libCEED call 154f2d2bf5dSjeremylt lib.CeedOperatorCreate(self._ceed._pointer[0], qf._pointer[0], 155f2d2bf5dSjeremylt dqf._pointer[0] if dqf else ffi.NULL, 156f2d2bf5dSjeremylt dqfT._pointer[0] if dqfT else ffi.NULL, 157f2d2bf5dSjeremylt self._pointer) 158f2d2bf5dSjeremylt 159f2d2bf5dSjeremylt # Add field to CeedOperator 160a8d32208Sjeremylt def set_field(self, fieldname, restriction, basis, vector): 161f2d2bf5dSjeremylt """Provide a field to an Operator for use by its QFunction. 162f2d2bf5dSjeremylt 163f2d2bf5dSjeremylt Args: 164f2d2bf5dSjeremylt fieldname: name of the field (to be matched with the same name used 165f2d2bf5dSjeremylt by QFunction) 166f2d2bf5dSjeremylt restriction: ElemRestriction 167f2d2bf5dSjeremylt basis: Basis in which the field resides or CEED_BASIS_COLLOCATED 168f2d2bf5dSjeremylt if collocated with quadrature points 169f2d2bf5dSjeremylt vector: Vector to be used by Operator or CEED_VECTOR_ACTIVE 170f2d2bf5dSjeremylt if field is active or CEED_VECTOR_NONE if using 171a8d32208Sjeremylt CEED_EVAL_WEIGHT in the QFunction""" 172f2d2bf5dSjeremylt 173f2d2bf5dSjeremylt # libCEED call 174f2d2bf5dSjeremylt fieldnameAscii = ffi.new("char[]", fieldname.encode('ascii')) 175f2d2bf5dSjeremylt lib.CeedOperatorSetField(self._pointer[0], fieldnameAscii, 176a8d32208Sjeremylt restriction._pointer[0], basis._pointer[0], 177f2d2bf5dSjeremylt vector._pointer[0]) 178f2d2bf5dSjeremylt 179f2d2bf5dSjeremylt# ------------------------------------------------------------------------------ 1807a7b0fa3SJed Brown 1817a7b0fa3SJed Brown 182f2d2bf5dSjeremyltclass CompositeOperator(_OperatorBase): 183f2d2bf5dSjeremylt """Ceed Composite Operator: composition of multiple Operators.""" 184f2d2bf5dSjeremylt 185f2d2bf5dSjeremylt # Constructor 186f2d2bf5dSjeremylt def __init__(self, ceed): 187f2d2bf5dSjeremylt # CeedOperator object 188f2d2bf5dSjeremylt self._pointer = ffi.new("CeedOperator *") 189f2d2bf5dSjeremylt 190f2d2bf5dSjeremylt # Reference to Ceed 191f2d2bf5dSjeremylt self._ceed = ceed 192f2d2bf5dSjeremylt # libCEED call 193f2d2bf5dSjeremylt lib.CeedCompositeOperatorCreate(self._ceed._pointer[0], self._pointer) 194f2d2bf5dSjeremylt 195f2d2bf5dSjeremylt # Add sub operators 196f2d2bf5dSjeremylt def add_sub(self, subop): 197f2d2bf5dSjeremylt """Add a sub-operator to a composite CeedOperator. 198f2d2bf5dSjeremylt 199f2d2bf5dSjeremylt Args: 200f2d2bf5dSjeremylt subop: sub-operator Operator""" 201f2d2bf5dSjeremylt 202f2d2bf5dSjeremylt # libCEED call 203f2d2bf5dSjeremylt lib.CeedCompositeOperatorAddSub(self._pointer[0], subop._pointer[0]) 204f2d2bf5dSjeremylt 205f2d2bf5dSjeremylt# ------------------------------------------------------------------------------ 206