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