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 20import contextlib 21from .ceed_constants import MEM_HOST, COPY_VALUES, NORM_2 22 23# ------------------------------------------------------------------------------ 24 25 26class Vector(): 27 """Ceed Vector: storing and manipulating vectors.""" 28 29 # Attributes 30 _ceed = ffi.NULL 31 _pointer = ffi.NULL 32 33 # Constructor 34 def __init__(self, ceed, size): 35 # CeedVector object 36 self._pointer = ffi.new("CeedVector *") 37 38 # Reference to Ceed 39 self._ceed = ceed 40 41 # libCEED call 42 lib.CeedVectorCreate(self._ceed._pointer[0], size, self._pointer) 43 44 # Destructor 45 def __del__(self): 46 # libCEED call 47 lib.CeedVectorDestroy(self._pointer) 48 49 # Representation 50 def __repr__(self): 51 return "<CeedVector instance at " + hex(id(self)) + ">" 52 53 # String conversion for print() to stdout 54 def __str__(self): 55 """View a Vector via print().""" 56 57 # libCEED call 58 fmt = ffi.new("char[]", "%f".encode('ascii')) 59 with tempfile.NamedTemporaryFile() as key_file: 60 with open(key_file.name, 'r+') as stream_file: 61 stream = ffi.cast("FILE *", stream_file) 62 63 lib.CeedVectorView(self._pointer[0], fmt, stream) 64 65 stream_file.seek(0) 66 out_string = stream_file.read() 67 68 return out_string 69 70 # Set Vector's data array 71 def set_array(self, array, memtype=MEM_HOST, cmode=COPY_VALUES): 72 """Set the array used by a Vector, freeing any previously allocated 73 array if applicable. 74 75 Args: 76 *array: Numpy or Numba array to be used 77 **memtype: memory type of the array being passed, default CEED_MEM_HOST 78 **cmode: copy mode for the array, default CEED_COPY_VALUES""" 79 80 # Setup the numpy array for the libCEED call 81 array_pointer = ffi.new("CeedScalar *") 82 if memtype == MEM_HOST: 83 array_pointer = ffi.cast( 84 "CeedScalar *", 85 array.__array_interface__['data'][0]) 86 else: 87 array_pointer = ffi.cast( 88 "CeedScalar *", 89 array.__cuda_array_interface__['data'][0]) 90 91 # libCEED call 92 lib.CeedVectorSetArray(self._pointer[0], memtype, cmode, array_pointer) 93 94 # Get Vector's data array 95 def get_array(self, memtype=MEM_HOST): 96 """Get read/write access to a Vector via the specified memory type. 97 98 Args: 99 **memtype: memory type of the array being passed, default CEED_MEM_HOST 100 101 Returns: 102 *array: Numpy or Numba array""" 103 104 # Retrieve the length of the array 105 length_pointer = ffi.new("CeedInt *") 106 lib.CeedVectorGetLength(self._pointer[0], length_pointer) 107 108 # Setup the pointer's pointer 109 array_pointer = ffi.new("CeedScalar **") 110 111 # libCEED call 112 lib.CeedVectorGetArray(self._pointer[0], memtype, array_pointer) 113 114 # Return array created from buffer 115 if memtype == MEM_HOST: 116 # Create buffer object from returned pointer 117 buff = ffi.buffer( 118 array_pointer[0], 119 ffi.sizeof("CeedScalar") * 120 length_pointer[0]) 121 # return Numpy array 122 return np.frombuffer(buff, dtype="float64") 123 else: 124 # CUDA array interface 125 # https://numba.pydata.org/numba-doc/latest/cuda/cuda_array_interface.html 126 import numba.cuda as nbcuda 127 desc = { 128 'shape': (length_pointer[0]), 129 'typestr': '>f8', 130 'data': (int(ffi.cast("intptr_t", array_pointer[0])), False), 131 'version': 2 132 } 133 # return Numba array 134 return nbcuda.from_cuda_array_interface(desc) 135 136 # Get Vector's data array in read-only mode 137 def get_array_read(self, memtype=MEM_HOST): 138 """Get read-only access to a Vector via the specified memory type. 139 140 Args: 141 **memtype: memory type of the array being passed, default CEED_MEM_HOST 142 143 Returns: 144 *array: Numpy or Numba array""" 145 146 # Retrieve the length of the array 147 length_pointer = ffi.new("CeedInt *") 148 lib.CeedVectorGetLength(self._pointer[0], length_pointer) 149 150 # Setup the pointer's pointer 151 array_pointer = ffi.new("CeedScalar **") 152 153 # libCEED call 154 lib.CeedVectorGetArrayRead(self._pointer[0], memtype, array_pointer) 155 156 # Return array created from buffer 157 if memtype == MEM_HOST: 158 # Create buffer object from returned pointer 159 buff = ffi.buffer( 160 array_pointer[0], 161 ffi.sizeof("CeedScalar") * 162 length_pointer[0]) 163 # return read only Numpy array 164 ret = np.frombuffer(buff, dtype="float64") 165 ret.flags['WRITEABLE'] = False 166 return ret 167 else: 168 # CUDA array interface 169 # https://numba.pydata.org/numba-doc/latest/cuda/cuda_array_interface.html 170 import numba.cuda as nbcuda 171 desc = { 172 'shape': (length_pointer[0]), 173 'typestr': '>f8', 174 'data': (int(ffi.cast("intptr_t", array_pointer[0])), False), 175 'version': 2 176 } 177 # return read only Numba array 178 return nbcuda.from_cuda_array_interface(desc) 179 180 # Restore the Vector's data array 181 def restore_array(self): 182 """Restore an array obtained using get_array().""" 183 184 # Setup the pointer's pointer 185 array_pointer = ffi.new("CeedScalar **") 186 187 # libCEED call 188 lib.CeedVectorRestoreArray(self._pointer[0], array_pointer) 189 190 # Restore an array obtained using getArrayRead 191 def restore_array_read(self): 192 """Restore an array obtained using get_array_read().""" 193 194 # Setup the pointer's pointer 195 array_pointer = ffi.new("CeedScalar **") 196 197 # libCEED call 198 lib.CeedVectorRestoreArrayRead(self._pointer[0], array_pointer) 199 200 @contextlib.contextmanager 201 def array(self, *shape, memtype=MEM_HOST): 202 """Context manager for array access. 203 204 Args: 205 shape (tuple): shape of returned numpy.array 206 **memtype: memory type of the array being passed, default CEED_MEM_HOST 207 208 209 Returns: 210 np.array: writable view of vector 211 212 Examples: 213 Constructing the identity inside a libceed.Vector: 214 215 >>> vec = ceed.Vector(16) 216 >>> with vec.array(4, 4) as x: 217 >>> x[...] = np.eye(4) 218 """ 219 x = self.get_array(memtype=memtype) 220 if shape: 221 x = x.reshape(shape) 222 yield x 223 self.restore_array() 224 225 @contextlib.contextmanager 226 def array_read(self, *shape, memtype=MEM_HOST): 227 """Context manager for read-only array access. 228 229 Args: 230 shape (tuple): shape of returned numpy.array 231 **memtype: memory type of the array being passed, default CEED_MEM_HOST 232 233 Returns: 234 np.array: read-only view of vector 235 236 Examples: 237 Constructing the identity inside a libceed.Vector: 238 239 >>> vec = ceed.Vector(6) 240 >>> vec.set_value(1.3) 241 >>> with vec.array_read(2, 3) as x: 242 >>> print(x) 243 """ 244 x = self.get_array_read(memtype=memtype) 245 if shape: 246 x = x.reshape(shape) 247 yield x 248 self.restore_array_read() 249 250 # Get the length of a Vector 251 def get_length(self): 252 """Get the length of a Vector. 253 254 Returns: 255 length: length of the Vector""" 256 257 length_pointer = ffi.new("CeedInt *") 258 259 # libCEED call 260 lib.CeedVectorGetLength(self._pointer[0], length_pointer) 261 262 return length_pointer[0] 263 264 # Get the length of a Vector 265 def __len__(self): 266 """Get the length of a Vector. 267 268 Returns: 269 length: length of the Vector""" 270 271 length_pointer = ffi.new("CeedInt *") 272 273 # libCEED call 274 lib.CeedVectorGetLength(self._pointer[0], length_pointer) 275 276 return length_pointer[0] 277 278 # Set the Vector to a given constant value 279 def set_value(self, value): 280 """Set the Vector to a constant value. 281 282 Args: 283 value: value to be used""" 284 285 # libCEED call 286 lib.CeedVectorSetValue(self._pointer[0], value) 287 288 # Sync the Vector to a specified memtype 289 def sync_array(self, memtype=MEM_HOST): 290 """Sync the Vector to a specified memtype. 291 292 Args: 293 **memtype: memtype to be synced""" 294 295 # libCEED call 296 lib.CeedVectorSyncArray(self._pointer[0], memtype) 297 298 # Compute the norm of a vector 299 def norm(self, normtype=NORM_2): 300 """Get the norm of a Vector. 301 302 Args: 303 **normtype: type of norm to be computed""" 304 305 norm_pointer = ffi.new("CeedScalar *") 306 307 # libCEED call 308 lib.CeedVectorNorm(self._pointer[0], normtype, norm_pointer) 309 310 return norm_pointer[0] 311 312# ------------------------------------------------------------------------------ 313 314 315class _VectorWrap(Vector): 316 """Wrap a CeedVector pointer in a Vector object.""" 317 318 # Constructor 319 def __init__(self, ceed, pointer): 320 # CeedVector object 321 self._pointer = pointer 322 323 # Reference to Ceed 324 self._ceed = ceed 325 326# ------------------------------------------------------------------------------ 327