GEOS
wrapperHelpers.hpp
Go to the documentation of this file.
1 /*
2  * ------------------------------------------------------------------------------------------------------------
3  * SPDX-License-Identifier: LGPL-2.1-only
4  *
5  * Copyright (c) 2016-2024 Lawrence Livermore National Security LLC
6  * Copyright (c) 2018-2024 TotalEnergies
7  * Copyright (c) 2018-2024 The Board of Trustees of the Leland Stanford Junior University
8  * Copyright (c) 2023-2024 Chevron
9  * Copyright (c) 2019- GEOS/GEOSX Contributors
10  * All rights reserved
11  *
12  * See top level LICENSE, COPYRIGHT, CONTRIBUTORS, NOTICE, and ACKNOWLEDGEMENTS files for details.
13  * ------------------------------------------------------------------------------------------------------------
14  */
15 
20 #ifndef GEOS_DATAREPOSITORY_WRAPPERHELPERS_HPP_
21 #define GEOS_DATAREPOSITORY_WRAPPERHELPERS_HPP_
22 
23 
25 #define RESTART_TYPE_LOGGING 0
26 
27 // Source includes
28 #include "BufferOps.hpp"
29 #include "BufferOpsDevice.hpp"
30 #include "DefaultValue.hpp"
31 #include "ConduitRestart.hpp"
32 #include "common/DataTypes.hpp"
33 #include "common/GeosxMacros.hpp"
34 #include "common/Span.hpp"
35 #include "codingUtilities/traits.hpp"
36 #include "LvArray/src/system.hpp"
37 
38 #if defined(GEOS_USE_PYGEOSX)
39 #include "LvArray/src/python/python.hpp"
40 #endif
41 
42 // TPL includes
43 #include <conduit.hpp>
44 
45 // System includes
46 #include <cstring>
47 
48 #if RESTART_TYPE_LOGGING
49 #include <unordered_set>
50 #endif
51 
52 namespace geos
53 {
54 namespace dataRepository
55 {
56 namespace wrapperHelpers
57 {
58 namespace internal
59 {
60 
61 inline void logOutputType( string const & typeString, string const & msg )
62 {
63 #if RESTART_TYPE_LOGGING
64  static std::unordered_set< string > m_types;
65 
66  if( !m_types.count( typeString ) )
67  {
68  m_types.insert( typeString );
69  GEOS_LOG( msg << typeString );
70  }
71 #else
72  GEOS_DEBUG_VAR( typeString );
73  GEOS_DEBUG_VAR( msg );
74 #endif
75 }
76 
77 template< typename T, typename ... INDICES >
78 string getIndicesToComponent( T const &, int const component, INDICES const ... existingIndices )
79 {
80  GEOS_ERROR_IF_NE( component, 0 );
81  return LvArray::indexing::getIndexString( existingIndices ... );
82 }
83 
84 template< typename ... INDICES >
85 string getIndicesToComponent( R1Tensor const &, int const component, INDICES const ... existingIndices )
86 { return LvArray::indexing::getIndexString( existingIndices ..., component ); }
87 
88 template< typename T >
89 T const * getPointerToComponent( T const & var, int const component )
90 {
91  GEOS_ERROR_IF_NE( component, 0 );
92  return &var;
93 }
94 
95 inline
96 real64 const * getPointerToComponent( R1Tensor const & var, int const component )
97 {
98  GEOS_ERROR_IF_GE( component, 3 );
99  return &var[ component ];
100 }
101 
102 } // namespace internal
103 
104 template< typename T >
105 class ArrayDimLabels
106 {
107 public:
108 
109  void set( integer const, Span< string const > )
110  {
111  GEOS_ERROR( "Dimension labels are only available in Array wrappers" );
112  }
113 
114  Span< string const > get( integer const ) const
115  {
116  GEOS_ERROR( "Dimension labels are only available in Array wrappers" );
117  return {};
118  }
119 };
120 
121 template< typename T, int NDIM, typename PERM >
122 class ArrayDimLabels< Array< T, NDIM, PERM > >
123 {
124 public:
125 
126  void set( integer const dim, Span< string const > labels )
127  {
128  GEOS_ERROR_IF_LT( dim, 0 );
129  GEOS_ERROR_IF_GE( dim, NDIM );
130  m_values[dim].resize( labels.size() );
131  std::copy( labels.begin(), labels.end(), m_values[dim].begin() );
132  }
133 
134  Span< string const > get( integer const dim ) const
135  {
136  GEOS_ERROR_IF_LT( dim, 0 );
137  GEOS_ERROR_IF_GE( dim, NDIM );
138  return { m_values[dim].data(), m_values[dim].size() };
139  }
140 
141 private:
142 
143  string_array m_values[NDIM]{};
144 };
145 
146 template< typename T >
147 inline std::enable_if_t< traits::HasMemberFunction_size< T >, size_t >
148 size( T const & value )
149 { return value.size(); }
150 
151 template< typename T >
152 inline std::enable_if_t< !traits::HasMemberFunction_size< T >, size_t >
153 size( T const & GEOS_UNUSED_PARAM( value ) )
154 { return 1; }
155 
156 
157 inline char *
158 dataPtr( string & var )
159 { return const_cast< char * >( var.data() ); }
160 
161 inline char *
162 dataPtr( Path & var )
163 { return const_cast< char * >( var.data() ); }
164 
165 template< typename T >
166 inline std::enable_if_t< traits::HasMemberFunction_data< T >, typename traits::Pointer< T > >
167 dataPtr( T & value )
168 { return value.data(); }
169 
170 template< typename T >
171 inline std::enable_if_t< !traits::HasMemberFunction_data< T >, typename traits::Pointer< T > >
172 dataPtr( T & value )
173 { return &value; }
174 
175 template< class T >
176 inline typename traits::ConstPointer< T >
177 dataPtr( T const & value )
178 { return dataPtr( const_cast< T & >( value ) ); }
179 
180 
181 template< typename T >
182 inline std::enable_if_t< traits::HasMemberFunction_resize< T > >
183 resize( T & value, localIndex const newSize )
184 { value.resize( newSize ); }
185 
186 template< typename T >
187 inline std::enable_if_t< !traits::HasMemberFunction_resize< T > >
188 resize( T & GEOS_UNUSED_PARAM( value ),
189  localIndex const GEOS_UNUSED_PARAM( newSize ) )
190 {}
191 
192 template< typename T >
193 inline void
194 resizeDefault( T & value,
195  localIndex const newSize,
196  DefaultValue< T > const & defaultValue,
197  string const & )
198 { value.resizeDefault( newSize, defaultValue.value ); }
199 
200 
201 
202 template< typename T, int NDIM, typename PERMUTATION >
203 inline void
204 resizeDimensions( Array< T, NDIM, PERMUTATION > & value, int num_dims, localIndex const * const dims )
205 { value.resize( num_dims, dims ); }
206 
207 template< typename T >
208 inline void
209 resizeDimensions( T & value, int num_dims, localIndex const * const dims )
210 {
211  if( num_dims != 1 )
212  {
213  GEOS_ERROR( "Data is not multidimensional" );
214  return;
215  }
216  resize( value, dims[ 0 ] );
217 }
218 
219 
220 template< typename T >
221 inline localIndex
222 byteSizeOfElement()
223 { return sizeof( *dataPtr( std::declval< T >() ) ); }
224 
225 
226 template< typename T >
227 inline size_t
228 byteSize( T const & value )
229 { return wrapperHelpers::size( value ) * byteSizeOfElement< T >(); }
230 
231 
232 template< typename T >
233 inline localIndex
234 numElementsFromByteSize( localIndex const byteSize )
235 {
236  GEOS_ERROR_IF_NE( byteSize % byteSizeOfElement< T >(), 0 );
237  return byteSize / byteSizeOfElement< T >();
238 }
239 
240 
242 namespace detail
243 {
244 template< typename T >
245 std::enable_if_t< traits::HasMemberFunction_reserveValues< T const > >
246 reserveValuesIfAvailable( T & value, localIndex const newCapacity )
247 {
248  localIndex const oldCapacity = value.size();
249  double const oldValueCapacity = value.valueCapacity();
250  localIndex const newValueCapacity = oldValueCapacity * newCapacity / oldCapacity;
251  value.reserveValues( newValueCapacity );
252 }
253 
254 template< typename T >
255 std::enable_if_t< !traits::HasMemberFunction_reserveValues< T const > >
256 reserveValuesIfAvailable( T &, localIndex const )
257 {}
258 }
260 
261 template< typename T >
262 std::enable_if_t< traits::HasMemberFunction_reserve< T > >
263 reserve( T & value, localIndex const newCapacity )
264 {
265  value.reserve( newCapacity );
266  detail::reserveValuesIfAvailable( value, newCapacity );
267 }
268 
269 template< typename T,
270  int NDIM,
271  typename PERMUTATION >
272 void reserve( Array< T, NDIM, PERMUTATION > & value, localIndex newCapacity )
273 {
274  localIndex const * const dims = value.dims();
275  for( int i=1; i<NDIM; ++i )
276  {
277  newCapacity *= dims[ i ];
278  }
279  value.reserve( newCapacity );
280 }
281 
282 template< typename T >
283 std::enable_if_t< !traits::HasMemberFunction_reserve< T > >
284 reserve( T & GEOS_UNUSED_PARAM( value ), localIndex const GEOS_UNUSED_PARAM( newCapacity ) )
285 {}
286 
287 
288 template< typename T >
289 std::enable_if_t< traits::HasMemberFunction_capacity< T const >, localIndex >
290 capacity( T const & value )
291 { return value.capacity(); }
292 
293 template< typename T >
294 std::enable_if_t< !traits::HasMemberFunction_capacity< T const >, localIndex >
295 capacity( T const & value )
296 { return wrapperHelpers::size( value ); }
297 
298 
299 
300 template< typename T >
301 std::enable_if_t< traits::HasMemberFunction_setName< T > >
302 setName( T & value, string const & name )
303 { value.setName( name ); }
304 
305 template< typename T >
306 std::enable_if_t< !traits::HasMemberFunction_setName< T > >
307 setName( T & GEOS_UNUSED_PARAM( value ), string const & GEOS_UNUSED_PARAM( name ) )
308 {}
309 
310 template< typename T >
311 std::enable_if_t< traits::HasMemberFunction_move< T > >
312 move( T & value, LvArray::MemorySpace const space, bool const touch )
313 { value.move( space, touch ); }
314 
315 template< typename T >
316 std::enable_if_t< !traits::HasMemberFunction_move< T > >
317 move( T & GEOS_UNUSED_PARAM( value ),
318  LvArray::MemorySpace const GEOS_UNUSED_PARAM( space ),
319  bool const GEOS_UNUSED_PARAM( touch ) )
320 {}
321 
322 // This is for an object that needs to be packed.
323 template< typename T >
324 std::enable_if_t< !bufferOps::can_memcpy< typename traits::Pointer< T > > >
325 pushDataToConduitNode( T const & var, conduit::Node & node )
326 {
327  internal::logOutputType( LvArray::system::demangleType( var ), "Packing for output: " );
328 
329  // Get the number of bytes in the packed object.
330  localIndex const byteSize = bufferOps::PackSize( var );
331 
332  // Create a conduit data type that describes the array.
333  conduit::DataType const dtype( conduitTypeInfo< buffer_unit_type >::id, byteSize );
334 
335  // Allocate the array in the "__values__" child.
336  conduit::Node & valuesNode = node[ "__values__" ];
337  valuesNode.set( dtype );
338 
339  // Get the pointer to the array and pack the object into it.
340  buffer_unit_type * buffer = valuesNode.value();
341  bufferOps::Pack< true >( buffer, var );
342 }
343 
344 // This is for an object that needs to be packed.
345 template< typename T >
346 std::enable_if_t< !bufferOps::can_memcpy< typename traits::Pointer< T > > >
347 pullDataFromConduitNode( T & var, conduit::Node const & node )
348 {
349  conduit::Node const & valuesNode = node.fetch_existing( "__values__" );
350 
351  // Get the number of bytes in the array and a pointer to the array.
352  localIndex const byteSize = valuesNode.dtype().number_of_elements();
353  buffer_unit_type const * buffer = valuesNode.value();
354 
355  // Unpack the object from the array.
356  localIndex const bytesRead = bufferOps::Unpack( buffer, var );
357  GEOS_ERROR_IF_NE( bytesRead, byteSize );
358 }
359 
360 // This is for an string since the type of char is different on different platforms :(.
361 inline
362 void
363 pushDataToConduitNode( string const & var, conduit::Node & node )
364 {
365  internal::logOutputType( LvArray::system::demangleType( var ), "Output via external pointer: " );
366 
367  constexpr int conduitTypeID = conduitTypeInfo< signed char >::id;
368  conduit::DataType const dtype( conduitTypeID, var.size() );
369 
370  signed char * const ptr = const_cast< signed char * >( reinterpret_cast< signed char const * >( var.data() ) );
371  node[ "__values__" ].set_external( dtype, ptr );
372 }
373 
374 // This is for Path since it derives from string. See overload for string.
375 inline
376 void
377 pushDataToConduitNode( Path const & var, conduit::Node & node )
378 {
379  pushDataToConduitNode( static_cast< string const & >(var), node );
380 }
381 
382 // This is for an object that doesn't need to be packed but isn't an LvArray.
383 template< typename T >
384 std::enable_if_t< bufferOps::can_memcpy< typename traits::Pointer< T > > >
385 pushDataToConduitNode( T const & var, conduit::Node & node )
386 {
387  internal::logOutputType( LvArray::system::demangleType( var ), "Output via external pointer: " );
388 
389  constexpr int conduitTypeID = conduitTypeInfo< typename traits::Pointer< T > >::id;
390  constexpr int sizeofConduitType = conduitTypeInfo< typename traits::Pointer< T > >::sizeOfConduitType;
391  localIndex const numBytes = byteSize( var );
392  conduit::DataType const dtype( conduitTypeID, numBytes / sizeofConduitType );
393 
394  void * const ptr = const_cast< void * >( static_cast< void const * >( dataPtr( var ) ) );
395  node[ "__values__" ].set_external( dtype, ptr );
396 }
397 
398 // This is for an object that doesn't need to be packed but isn't an LvArray or a SortedArray.
399 template< typename T >
400 std::enable_if_t< bufferOps::can_memcpy< typename traits::Pointer< T > > >
401 pullDataFromConduitNode( T & var, conduit::Node const & node )
402 {
403  conduit::Node const & valuesNode = node.fetch_existing( "__values__" );
404 
405  localIndex const byteSize = LvArray::integerConversion< localIndex >( valuesNode.dtype().strided_bytes() );
406  localIndex const numElements = numElementsFromByteSize< T >( byteSize );
407 
408  resize( var, numElements );
409 
410  std::memcpy( dataPtr( var ), valuesNode.data_ptr(), byteSize );
411 }
412 
413 // This is for a SortedArray that doesn't need to be packed.
414 template< typename T >
415 std::enable_if_t< bufferOps::can_memcpy< T > >
416 pullDataFromConduitNode( SortedArray< T > & var, conduit::Node const & node )
417 {
418  conduit::Node const & valuesNode = node.fetch_existing( "__values__" );
419 
420  localIndex const byteSize = LvArray::integerConversion< localIndex >( valuesNode.dtype().strided_bytes() );
421  localIndex const numElements = numElementsFromByteSize< T >( byteSize );
422 
423  T const * const values = reinterpret_cast< T const * >( valuesNode.data_ptr() );
424  var.insert( values, values + numElements );
425 }
426 
427 
428 // This is an LvArray that doesn't need to be packed.
429 template< typename T, int NDIM, typename PERMUTATION >
430 std::enable_if_t< bufferOps::can_memcpy< T > >
431 pushDataToConduitNode( Array< T, NDIM, PERMUTATION > const & var,
432  conduit::Node & node )
433 {
434  internal::logOutputType( LvArray::system::demangleType( var ), "Output array via external pointer: " );
435 
436  // Push the data into conduit
437  constexpr int conduitTypeID = conduitTypeInfo< T >::id;
438  constexpr int sizeofConduitType = conduitTypeInfo< T >::sizeOfConduitType;
439  conduit::DataType const dtype( conduitTypeID, var.size() * sizeof( T ) / sizeofConduitType );
440  void * const ptr = const_cast< void * >( static_cast< void const * >( var.data() ) );
441  node[ "__values__" ].set_external( dtype, ptr );
442 
443  // Create a copy of the dimensions
444  camp::idx_t temp[ NDIM + 1 ];
445  for( int i = 0; i < NDIM; ++i )
446  {
447  temp[ i ] = var.size( i );
448  }
449 
450  // If T is something like a Tensor than there is an extra implicit dimension.
451  constexpr int const implicitDimensionLength = conduitTypeInfo< T >::numConduitValues;
452  constexpr bool const hasImplicitDimension = implicitDimensionLength != 1;
453  constexpr int totalNumDimensions = NDIM + hasImplicitDimension;
454  if( hasImplicitDimension )
455  {
456  temp[ NDIM ] = implicitDimensionLength;
457  }
458 
459  // push the dimensions into the node
460  conduit::DataType const dimensionType( conduitTypeInfo< camp::idx_t >::id, totalNumDimensions );
461  node[ "__dimensions__" ].set( dimensionType, temp );
462 
463  // Create a copy of the permutation
464  constexpr std::array< camp::idx_t, NDIM > const perm = to_stdArray( RAJA::as_array< PERMUTATION >::get());
465  for( int i = 0; i < NDIM; ++i )
466  {
467  temp[ i ] = perm[ i ];
468  }
469 
470  if( hasImplicitDimension )
471  {
472  temp[ NDIM ] = NDIM;
473  }
474 
475  node[ "__permutation__" ].set( dimensionType, temp );
476 }
477 
478 // This is an LvArray that doesn't need to be packed.
479 template< typename T, int NDIM, typename PERMUTATION >
480 std::enable_if_t< bufferOps::can_memcpy< T > >
481 pullDataFromConduitNode( Array< T, NDIM, PERMUTATION > & var,
482  conduit::Node const & node )
483 {
484  // Get the number of dimensions written out, accounting for an implicit dimension and the permutation.
485  constexpr int const implicitDimensionLength = conduitTypeInfo< T >::numConduitValues;
486  constexpr bool const hasImplicitDimension = implicitDimensionLength != 1;
487  constexpr int totalNumDimensions = NDIM + hasImplicitDimension;
488 
489  // Check that the permutations match.
490  conduit::Node const & permutationNode = node.fetch_existing( "__permutation__" );
491  GEOS_ERROR_IF_NE( permutationNode.dtype().number_of_elements(), totalNumDimensions );
492 
493  constexpr std::array< camp::idx_t, NDIM > const perm = to_stdArray( RAJA::as_array< PERMUTATION >::get());
494  camp::idx_t const * const permFromConduit = permutationNode.value();
495  for( int i = 0; i < NDIM; ++i )
496  {
497  GEOS_ERROR_IF_NE_MSG( permFromConduit[ i ], perm[ i ],
498  "The permutation of the data in conduit and the provided Array don't match." );
499  }
500 
501  if( hasImplicitDimension )
502  {
503  GEOS_ERROR_IF_NE_MSG( permFromConduit[ NDIM ], NDIM,
504  "The permutation of the data in conduit and the provided Array don't match." );
505  }
506 
507  // Now pull out the dimensions and resize the array.
508  conduit::Node const & dimensionNode = node.fetch_existing( "__dimensions__" );
509  GEOS_ERROR_IF_NE( dimensionNode.dtype().number_of_elements(), totalNumDimensions );
510  camp::idx_t const * const dims = dimensionNode.value();
511 
512  if( hasImplicitDimension )
513  {
514  GEOS_ERROR_IF_NE( dims[ NDIM ], implicitDimensionLength );
515  }
516 
517  var.resize( NDIM, dims );
518 
519  // Finally memcpy
520  conduit::Node const & valuesNode = node.fetch_existing( "__values__" );
521  localIndex numBytesFromArray = var.size() * sizeof( T );
522  GEOS_ERROR_IF_NE( numBytesFromArray, valuesNode.dtype().strided_bytes() );
523  if( numBytesFromArray > 0 )
524  {
525  std::memcpy( var.data(), valuesNode.data_ptr(), numBytesFromArray );
526  }
527 }
528 
529 
530 
531 template< typename T, typename INDEX_TYPE >
532 std::enable_if_t< bufferOps::can_memcpy< T > >
533 pushDataToConduitNode( ArrayOfArrays< T, INDEX_TYPE > const & var2,
534  conduit::Node & node )
535 {
536  ArrayOfArraysView< T const, INDEX_TYPE > const & var = var2.toViewConst();
537  internal::logOutputType( LvArray::system::demangleType( var ), "Output array via external pointer: " );
538 
539  // ArrayOfArrays::m_numArrays
540  INDEX_TYPE const numArrays = var.size();
541  conduit::DataType const numArraysType( conduitTypeInfo< INDEX_TYPE >::id, 1 );
542  node[ "__numberOfArrays__" ].set( numArraysType, const_cast< void * >( static_cast< void const * >(&numArrays) ) );
543 
544  // ArrayOfArrays::m_offsets
545  INDEX_TYPE const * const offsets = var.getOffsets();
546  conduit::DataType const offsetsType( conduitTypeInfo< INDEX_TYPE >::id, numArrays+1 );
547  node[ "__offsets__" ].set_external( offsetsType, const_cast< void * >( static_cast< void const * >( offsets ) ) );
548 
549  // ArrayOfArrays::m_sizes
550  INDEX_TYPE const * const sizes = var.getSizes();
551  conduit::DataType const sizesType( conduitTypeInfo< INDEX_TYPE >::id, numArrays );
552  node[ "__sizes__" ].set_external( sizesType, const_cast< void * >( static_cast< void const * >( sizes ) ) );
553 
554  // **** WARNING: alters the uninitialized values in the ArrayOfArrays ****
555  T * const values = const_cast< T * >(var.getValues());
556  for( INDEX_TYPE i = 0; i < numArrays; ++i )
557  {
558  INDEX_TYPE const curOffset = offsets[ i ];
559  INDEX_TYPE const nextOffset = offsets[ i + 1 ];
560  for( INDEX_TYPE j = curOffset + var.sizeOfArray( i ); j < nextOffset; ++j )
561  {
562  if constexpr ( std::is_arithmetic< T >::value )
563  {
564  values[ j ] = 0;
565  }
566  else
567  {
568  values[ j ] = T();
569  }
570  }
571  }
572 
573  constexpr int conduitTypeID = conduitTypeInfo< T >::id;
574  constexpr int sizeofConduitType = conduitTypeInfo< T >::sizeOfConduitType;
575  conduit::DataType const dtype( conduitTypeID, offsets[numArrays] * sizeof( T ) / sizeofConduitType );
576 
577  // Push the data into conduit
578  node[ "__values__" ].set_external( dtype, values );
579 }
580 
581 template< typename T, typename INDEX_TYPE >
582 std::enable_if_t< bufferOps::can_memcpy< T > >
583 pullDataFromConduitNode( ArrayOfArrays< T, INDEX_TYPE > & var,
584  conduit::Node const & node )
585 {
586 
587  // numArrays node
588  conduit::Node const & numArraysNode = node.fetch_existing( "__numberOfArrays__" );
589  INDEX_TYPE const * const numArrays = numArraysNode.value();
590 
591  // offsets node
592  conduit::Node const & offsetsNode = node.fetch_existing( "__offsets__" );
593  conduit::DataType const & offsetsDataType = offsetsNode.dtype();
594  INDEX_TYPE const * const offsets = offsetsNode.value();
595  INDEX_TYPE const sizeOffsets = offsetsDataType.number_of_elements();
596 
597  // sizes node
598  conduit::Node const & sizesNode = node.fetch_existing( "__sizes__" );
599  conduit::DataType const & sizesDataType = sizesNode.dtype();
600  INDEX_TYPE const * const sizes = sizesNode.value();
601  INDEX_TYPE const sizeSizes = sizesDataType.number_of_elements();
602 
603  // Check that the numArrays, sizes and offsets are consistent.
604  GEOS_ERROR_IF_NE( *numArrays, sizeSizes );
605  GEOS_ERROR_IF_NE( *numArrays+1, sizeOffsets );
606 
607  // values node
608  conduit::Node const & valuesNode = node.fetch_existing( "__values__" );
609  conduit::DataType const & valuesDataType = valuesNode.dtype();
610  const INDEX_TYPE valuesSize = valuesDataType.number_of_elements();
611 
612  // should preallocate var.m_values with estimated sizes
613  INDEX_TYPE const arraySizeEstimate = (*numArrays)==0 ? 0 : valuesSize / (*numArrays);
614  var.resize( *numArrays, arraySizeEstimate );
615  var.reserveValues( valuesSize );
616 
617  // correctly set the sizes and capacities of each sub-array
618  localIndex allocatedSize = 0;
619  for( INDEX_TYPE i = 0; i < *numArrays; ++i )
620  {
621  INDEX_TYPE const arrayAllocation = offsets[i+1] - offsets[i];
622  var.setCapacityOfArray( i, arrayAllocation );
623  var.resizeArray( i, sizes[ i ] );
624  allocatedSize += arrayAllocation;
625  }
626 
627  // make sure that the allocated size is the same as the number of values read
628  GEOS_ERROR_IF_NE( valuesSize, allocatedSize );
629 
630  // make sure the allocatedSize is consistent wit the last offset
631  GEOS_ERROR_IF_NE( allocatedSize, offsets[sizeOffsets-1] );
632 
633  // get a view because the ArrayOfArraysView data accessors are protected
634  ArrayOfArraysView< T const, INDEX_TYPE > const & varView = var.toViewConst();
635  INDEX_TYPE const * const varOffsets = varView.getOffsets();
636  INDEX_TYPE const * const varSizes = varView.getSizes();
637 
638  // check that the offsets that are read are the same as the ones that were allocated
639  GEOS_ERROR_IF_NE( varOffsets[0], offsets[0] );
640 
641  // check each subarray has the identical capacity and size
642  for( INDEX_TYPE i = 0; i<*numArrays; ++i )
643  {
644  GEOS_ERROR_IF_NE( varOffsets[i+1], offsets[i+1] );
645  GEOS_ERROR_IF_NE( varSizes[i], sizes[i] );
646  }
647 
648  // copy the values
649  localIndex numBytesFromArray = allocatedSize * sizeof( T );
650  GEOS_ERROR_IF_NE( numBytesFromArray, valuesDataType.strided_bytes() );
651  if( numBytesFromArray > 0 )
652  {
653  std::memcpy( const_cast< T * >(varView.getValues()), valuesNode.data_ptr(), numBytesFromArray );
654  }
655 }
656 
657 
658 
659 template< typename T >
660 void pushDataToConduitNode( InterObjectRelation< T > const & var,
661  conduit::Node & node )
662 {return pushDataToConduitNode( var.base(), node ); }
663 
664 template< typename T >
665 void pullDataFromConduitNode( InterObjectRelation< T > & var,
666  conduit::Node const & node )
667 { return pullDataFromConduitNode( var.base(), node ); }
668 
669 
671 template< typename T, int NDIM, int USD >
672 std::enable_if_t< std::is_arithmetic< T >::value || traits::is_tensorT< T > >
673 addBlueprintField( ArrayView< T const, NDIM, USD > const & var,
674  conduit::Node & fields,
675  string const & fieldName,
676  string const & topology,
677  stdVector< string > const & componentNames )
678 {
679  GEOS_ERROR_IF_LE( var.size(), 0 );
680 
681  using ConduitType = typename conduitTypeInfo< T >::type;
682  constexpr int conduitTypeID = conduitTypeInfo< T >::id;
683  constexpr int numComponentsPerValue = conduitTypeInfo< T >::numConduitValues;
684 
685  localIndex const totalNumberOfComponents = numComponentsPerValue * var.size() / var.size( 0 );
686  if( !componentNames.empty() )
687  {
688  GEOS_ERROR_IF_NE( localIndex( componentNames.size() ), totalNumberOfComponents );
689  }
690 
691  var.move( hostMemorySpace, false );
692 
693  conduit::DataType dtype( conduitTypeID, var.size( 0 ) );
694  dtype.set_stride( sizeof( ConduitType ) * numComponentsPerValue * var.strides()[ 0 ] );
695 
696  localIndex curComponent = 0;
697  LvArray::forValuesInSliceWithIndices( var[ 0 ], [&fields, &fieldName, &topology, &componentNames, totalNumberOfComponents, &dtype, &curComponent]
698  ( T const & val, auto const ... indices )
699  {
700  for( int i = 0; i < numComponentsPerValue; ++i )
701  {
702  string name;
703  if( totalNumberOfComponents == 1 )
704  {
705  name = fieldName;
706  }
707  else if( componentNames.empty() )
708  {
709  string indexString = internal::getIndicesToComponent( val, i, indices ... );
710  indexString.erase( indexString.begin() );
711  indexString.pop_back();
712  indexString.pop_back();
713  name = fieldName + indexString;
714  }
715  else
716  {
717  name = componentNames[ curComponent++ ];
718  }
719 
720  conduit::Node & field = fields[ name ];
721  field[ "association" ] = "element";
722  field[ "volume_dependent" ] = "false";
723  field[ "topology" ] = topology;
724 
725  void const * pointer = internal::getPointerToComponent( val, i );
726  field[ "values" ].set_external( dtype, const_cast< void * >( pointer ) );
727  }
728  } );
729 }
730 
731 template< typename T >
732 void addBlueprintField( T const &,
733  conduit::Node & fields,
734  string const &,
735  string const &,
736  stdVector< string > const & )
737 {
738  GEOS_ERROR( GEOS_FMT( "Cannot create a mcarray out of {}\nWas trying to write it to {}",
739  LvArray::system::demangleType< T >(),
740  fields.path() ) );
741  GEOS_UNUSED_VAR( fields );
742 }
743 
744 template< typename T, int NDIM, int USD >
745 std::enable_if_t< std::is_arithmetic< T >::value || traits::is_tensorT< T > >
746 populateMCArray( ArrayView< T const, NDIM, USD > const & var,
747  conduit::Node & node,
748  stdVector< string > const & componentNames )
749 {
750  GEOS_ERROR_IF_LE( var.size(), 0 );
751 
752  using ConduitType = typename conduitTypeInfo< T >::type;
753  constexpr int conduitTypeID = conduitTypeInfo< T >::id;
754  constexpr int numComponentsPerValue = conduitTypeInfo< T >::numConduitValues;
755 
756  if( !componentNames.empty() )
757  {
758  GEOS_ERROR_IF_NE( localIndex( componentNames.size() ), numComponentsPerValue * var.size() / var.size( 0 ) );
759  }
760 
761  var.move( hostMemorySpace, false );
762 
763  conduit::DataType dtype( conduitTypeID, var.size( 0 ) );
764  dtype.set_stride( sizeof( ConduitType ) * numComponentsPerValue * var.strides()[ 0 ] );
765 
766  localIndex curComponent = 0;
767  LvArray::forValuesInSliceWithIndices( var[ 0 ], [&componentNames, &node, &dtype, &curComponent]
768  ( T const & val, auto const ... indices )
769  {
770  for( int i = 0; i < numComponentsPerValue; ++i )
771  {
772  string const name = componentNames.empty() ? internal::getIndicesToComponent( val, i, indices ... ) :
773  componentNames[ curComponent++ ];
774 
775  void const * pointer = internal::getPointerToComponent( val, i );
776  node[ name ].set_external( dtype, const_cast< void * >( pointer ) );
777  }
778  } );
779 }
780 
781 template< typename T >
782 void populateMCArray( T const &,
783  conduit::Node & node,
784  stdVector< string > const & )
785 {
786  GEOS_ERROR( GEOS_FMT( "Cannot create a mcarray out of {}\nWas trying to write it to {}",
787  LvArray::system::demangleType< T >(),
788  node.path() ) );
789  GEOS_UNUSED_VAR( node );
790 }
791 
792 template< typename T >
793 std::enable_if_t< std::is_arithmetic< T >::value, std::unique_ptr< Array< T, 1 > > >
794 averageOverSecondDim( ArrayView< T const, 1, 0 > const & var )
795 {
796  std::unique_ptr< Array< T, 1 > > ret = std::make_unique< Array< T, 1 > >();
797 
798  ret->resize( var.size() );
799  ret->template setValues< serialPolicy >( var );
800 
801  return ret;
802 }
803 
804 template< typename T, int NDIM, int USD >
805 std::enable_if_t< std::is_arithmetic< T >::value, std::unique_ptr< Array< T, NDIM - 1 > > >
806 averageOverSecondDim( ArrayView< T const, NDIM, USD > const & var )
807 {
808  std::unique_ptr< Array< T, NDIM - 1 > > ret = std::make_unique< Array< T, NDIM - 1 > >();
809 
810  localIndex newDims[ NDIM - 1 ];
811  newDims[ 0 ] = var.size( 0 );
812  for( int i = 2; i < NDIM; ++i )
813  {
814  newDims[ i - 1 ] = var.size( i );
815  }
816 
817  ret->resize( NDIM - 1, newDims );
818 
819  ArrayView< T, NDIM - 1 > const & output = *ret;
820 
821  localIndex const numSamples = var.size( 1 );
822  forAll< serialPolicy >( var.size( 0 ), [var, numSamples, &output] ( localIndex const i )
823  {
824  LvArray::sumOverFirstDimension( var[ i ], output[ i ] );
825 
826  LvArray::forValuesInSlice( output[ i ], [numSamples] ( T & val )
827  {
828  val /= numSamples;
829  } );
830  } );
831 
832  return ret;
833 }
834 
835 template< typename T >
836 std::unique_ptr< int > averageOverSecondDim( T const & )
837 {
838  GEOS_ERROR( GEOS_FMT( "Cannot average over the second dimension of {}", LvArray::system::demangleType< T >() ) );
839  return std::unique_ptr< int >( nullptr );
840 }
841 
842 template< typename T, int NDIM, int USD >
843 int numArrayDims( ArrayView< T const, NDIM, USD > const & GEOS_UNUSED_PARAM( var ) )
844 {
845  return NDIM;
846 }
847 
848 template< typename T >
849 int numArrayDims( stdVector< T > const & GEOS_UNUSED_PARAM( var ) )
850 {
851  return 1;
852 }
853 
854 template< typename T >
855 int numArrayDims( T const & GEOS_UNUSED_PARAM( var ) )
856 {
857  return 0;
858 }
859 
860 template< typename T, int NDIM, int USD >
861 localIndex numArrayComp( ArrayView< T const, NDIM, USD > const & var )
862 {
863  return LvArray::indexing::multiplyAll< NDIM - 1 >( var.dims() + 1 );
864 }
865 
866 template< typename T >
867 localIndex numArrayComp( ArrayView< T const, 1, 0 > const & GEOS_UNUSED_PARAM( var ) )
868 {
869  return 1;
870 }
871 
872 template< typename T >
873 localIndex numArrayComp( T const & GEOS_UNUSED_PARAM( var ) )
874 {
875  return 0;
876 }
877 
878 template< bool DO_PACKING, typename T, typename IDX >
879 inline std::enable_if_t< bufferOps::is_packable_by_index< T >, localIndex >
880 PackByIndex( buffer_unit_type * & buffer, T & var, IDX & idx )
881 { return bufferOps::PackByIndex< DO_PACKING >( buffer, var, idx ); }
882 
883 template< bool DO_PACKING, typename T, typename IDX >
884 inline std::enable_if_t< !bufferOps::is_packable_by_index< T >, localIndex >
885 PackByIndex( buffer_unit_type * &, T &, IDX & )
886 {
887  GEOS_ERROR( GEOS_FMT( "Trying to pack data type ({}) by index. Operation not supported.",
888  LvArray::system::demangleType< T >() ) );
889  return 0;
890 }
891 
892 template< typename T, typename IDX >
893 inline std::enable_if_t< bufferOps::is_packable_by_index< T >, localIndex >
894 UnpackByIndex( buffer_unit_type const * & buffer, T & var, IDX & idx )
895 { return bufferOps::UnpackByIndex( buffer, var, idx ); }
896 
897 template< typename T, typename IDX >
898 inline std::enable_if_t< !bufferOps::is_packable_by_index< T >, localIndex >
899 UnpackByIndex( buffer_unit_type const * &, T &, IDX & )
900 {
901  GEOS_ERROR( GEOS_FMT( "Trying to unpack data type ({}) by index. Operation not supported.",
902  LvArray::system::demangleType< T >() ) );
903  return 0;
904 }
905 
906 
907 template< bool DO_PACKING, typename T >
908 inline std::enable_if_t< bufferOps::is_container< T > || bufferOps::can_memcpy< T >, localIndex >
909 PackDevice( buffer_unit_type * & buffer, T const & var, parallelDeviceEvents & events )
910 { return bufferOps::PackDevice< DO_PACKING >( buffer, var, events ); }
911 
912 
913 template< bool DO_PACKING, typename T >
914 inline std::enable_if_t< !bufferOps::is_container< T > && !bufferOps::can_memcpy< T >, localIndex >
915 PackDevice( buffer_unit_type * &, T const &, parallelDeviceEvents & )
916 {
917  GEOS_ERROR( GEOS_FMT( "Trying to pack data type ({}) on device. Operation not supported.",
918  LvArray::system::demangleType< T >() ) );
919  return 0;
920 }
921 
922 template< bool DO_PACKING, typename T, typename IDX >
923 inline std::enable_if_t< bufferOps::is_container< T >, localIndex >
924 PackByIndexDevice( buffer_unit_type * & buffer, T const & var, IDX & idx, parallelDeviceEvents & events )
925 { return bufferOps::PackByIndexDevice< DO_PACKING >( buffer, var, idx, events ); }
926 
927 template< bool DO_PACKING, typename T, typename IDX >
928 inline std::enable_if_t< !bufferOps::is_container< T >, localIndex >
929 PackByIndexDevice( buffer_unit_type * &, T const &, IDX &, parallelDeviceEvents & )
930 {
931  GEOS_ERROR( GEOS_FMT( "Trying to pack data type ({}) by index on device. Operation not supported.",
932  LvArray::system::demangleType< T >() ) );
933  return 0;
934 }
935 
936 template< typename T >
937 inline std::enable_if_t< bufferOps::is_container< T >, localIndex >
938 UnpackDevice( buffer_unit_type const * & buffer, T const & var, parallelDeviceEvents & events )
939 { return bufferOps::UnpackDevice( buffer, var, events ); }
940 
941 template< typename T >
942 inline std::enable_if_t< !bufferOps::is_container< T >, localIndex >
943 UnpackDevice( buffer_unit_type const * &, T const &, parallelDeviceEvents & )
944 {
945  GEOS_ERROR( GEOS_FMT( "Trying to unpack data type ({}) on device. Operation not supported.",
946  LvArray::system::demangleType< T >() ) );
947  return 0;
948 }
949 
950 template< typename T, typename IDX >
951 inline std::enable_if_t< bufferOps::is_container< T >, localIndex >
952 UnpackByIndexDevice( buffer_unit_type const * & buffer, T const & var, IDX & idx, parallelDeviceEvents & events, MPI_Op op=MPI_REPLACE )
953 { return bufferOps::UnpackByIndexDevice( buffer, var, idx, events, op ); }
954 
955 template< typename T, typename IDX >
956 inline std::enable_if_t< !bufferOps::is_container< T >, localIndex >
957 UnpackByIndexDevice( buffer_unit_type const * &, T &, IDX &, parallelDeviceEvents &, MPI_Op )
958 {
959  GEOS_ERROR( GEOS_FMT( "Trying to unpack data type ({}) by index on device. Operation not supported.",
960  LvArray::system::demangleType< T >() ) );
961  return 0;
962 }
963 
964 
965 template< bool DO_PACKING, typename T >
967 PackDataDevice( buffer_unit_type * & buffer, T const & var, parallelDeviceEvents & events )
968 { return bufferOps::PackDataDevice< DO_PACKING >( buffer, var, events ); }
969 
970 template< bool DO_PACKING, typename T, typename IDX >
971 inline std::enable_if_t< bufferOps::is_container< T >, localIndex >
972 PackDataByIndexDevice( buffer_unit_type * & buffer, T const & var, IDX & idx, parallelDeviceEvents & events )
973 { return bufferOps::PackDataByIndexDevice< DO_PACKING >( buffer, var, idx, events ); }
974 
975 template< bool DO_PACKING, typename T, typename IDX >
976 inline std::enable_if_t< !bufferOps::is_container< T >, localIndex >
977 PackDataByIndexDevice( buffer_unit_type * &, T const &, IDX &, parallelDeviceEvents & )
978 {
979  GEOS_ERROR( GEOS_FMT( "Trying to pack data type ({}) by index on device. Operation not supported.",
980  LvArray::system::demangleType< T >() ) );
981  return 0;
982 }
983 
984 template< typename T >
985 inline std::enable_if_t< bufferOps::is_container< T >, localIndex >
986 UnpackDataDevice( buffer_unit_type const * & buffer, T const & var, parallelDeviceEvents & events )
987 { return bufferOps::UnpackDataDevice( buffer, var, events ); }
988 
989 template< typename T >
990 inline std::enable_if_t< !bufferOps::is_container< T >, localIndex >
991 UnpackDataDevice( buffer_unit_type const * &, T const &, parallelDeviceEvents & )
992 {
993  GEOS_ERROR( GEOS_FMT( "Trying to unpack data type ({}) on device. Operation not supported.",
994  LvArray::system::demangleType< T >() ) );
995  return 0;
996 }
997 
998 template< typename T, typename IDX >
999 inline std::enable_if_t< bufferOps::is_container< T >, localIndex >
1000 UnpackDataByIndexDevice( buffer_unit_type const * & buffer, T const & var, IDX & idx, parallelDeviceEvents & events, MPI_Op op )
1001 { return bufferOps::UnpackDataByIndexDevice( buffer, var, idx, events, op ); }
1002 
1003 template< typename T, typename IDX >
1004 inline std::enable_if_t< !bufferOps::is_container< T >, localIndex >
1005 UnpackDataByIndexDevice( buffer_unit_type const * &, T const &, IDX &, parallelDeviceEvents &, MPI_Op )
1006 {
1007  GEOS_ERROR( GEOS_FMT( "Trying to unpack data type ({}) by index on device. Operation not supported.",
1008  LvArray::system::demangleType< T >() ) );
1009  return 0;
1010 }
1011 
1012 #if defined(GEOS_USE_PYGEOSX)
1013 
1014 template< typename T >
1015 inline std::enable_if_t< LvArray::python::CanCreate< T >, PyObject * >
1016 createPythonObject( T & object )
1017 { return LvArray::python::create( object ); }
1018 
1019 template< typename T >
1020 inline std::enable_if_t< !LvArray::python::CanCreate< T >, PyObject * >
1021 createPythonObject( T & )
1022 { return nullptr; }
1023 
1024 #endif
1025 
1026 } // namespace wrapperHelpers
1027 } // namespace dataRepository
1028 } // namespace geos
1029 
1030 #undef RESTART_TYPE_LOGGING
1031 
1032 #endif // GEOS_DATAREPOSITORY_WRAPPERHELPERS_HPP_
#define GEOS_UNUSED_VAR(...)
Mark an unused variable and silence compiler warnings.
#define GEOS_DEBUG_VAR(...)
Mark a debug variable and silence compiler warnings.
#define GEOS_UNUSED_PARAM(X)
Mark an unused argument and silence compiler warnings.
Definition: GeosxMacros.hpp:97
#define GEOS_ERROR_IF_LE(lhs, rhs)
Raise a hard error if one value compares less than or equal to the other.
Definition: Logger.hpp:545
#define GEOS_ERROR(...)
Raise a hard error and terminate the program.
Definition: Logger.hpp:225
#define GEOS_ERROR_IF_LT(lhs, rhs)
Raise a hard error if one value compares less than the other.
Definition: Logger.hpp:527
#define GEOS_ERROR_IF_GE(lhs, rhs)
Raise a hard error if one value compares greater than or equal to the other.
Definition: Logger.hpp:510
#define GEOS_LOG(...)
Log a message on screen.
Definition: Logger.hpp:38
#define GEOS_ERROR_IF_NE_MSG(lhs, rhs,...)
Raise a hard error if two values are not equal.
Definition: Logger.hpp:469
#define GEOS_ERROR_IF_NE(lhs, rhs)
Raise a hard error if two values are not equal.
Definition: Logger.hpp:476
stdVector< string > string_array
A 1-dimensional array of geos::string types.
Definition: DataTypes.hpp:361
LvArray::Array< T, NDIM, PERMUTATION, localIndex, LvArray::ChaiBuffer > Array
Multidimensional array type. See LvArray:Array for details.
Definition: DataTypes.hpp:141
std::set< T > set
A set of local indices.
Definition: DataTypes.hpp:262
double real64
64-bit floating point type.
Definition: DataTypes.hpp:98
GEOS_LOCALINDEX_TYPE localIndex
Local index type (for indexing objects within an MPI partition).
Definition: DataTypes.hpp:84
constexpr stdArray< T, N > to_stdArray(std::array< T, N > const &arr)
Convert an std::array to an stdArray.
Tensor< real64, 3 > R1Tensor
Alias for a local (stack-based) rank-1 tensor type.
Definition: DataTypes.hpp:165
signed char buffer_unit_type
Type stored in communication buffers.
Definition: DataTypes.hpp:108
int integer
Signed integer type.
Definition: DataTypes.hpp:81
LvArray::ArrayView< T, NDIM, USD, localIndex, LvArray::ChaiBuffer > ArrayView
Multidimensional array view type. See LvArray:ArrayView for details.
Definition: DataTypes.hpp:147