20 #ifndef GEOS_FUNCTIONS_TABLEFUNCTION_HPP_
21 #define GEOS_FUNCTIONS_TABLEFUNCTION_HPP_
26 #include "LvArray/src/tensorOps.hpp"
85 m_coordinates = std::move( other.m_coordinates );
86 m_values = std::move( other.m_values );
87 m_interpolationMethod = other.m_interpolationMethod;
99 template<
typename IN_ARRAY >
109 template<
typename IN_ARRAY,
typename OUT_ARRAY >
120 void move( LvArray::MemorySpace
const space,
bool const touch )
122 m_coordinates.move( space, touch );
123 m_values.move( space, touch );
150 template<
typename IN_ARRAY >
153 interpolateLinear( IN_ARRAY
const & input )
const;
161 template<
typename IN_ARRAY,
typename OUT_ARRAY >
164 interpolateLinear( IN_ARRAY
const & input, OUT_ARRAY && derivatives )
const;
171 template<
typename IN_ARRAY >
174 interpolateRound( IN_ARRAY
const & input )
const;
183 template<
typename IN_ARRAY >
194 template<
typename IN_ARRAY,
typename OUT_ARRAY >
197 interpolateRound( IN_ARRAY
const & input, OUT_ARRAY && derivatives )
const;
247 FunctionBase::evaluateT< TableFunction, parallelHostPolicy >( group, time,
set, result );
314 return size_t(dim) < m_dimUnits.size() ? m_dimUnits[dim] : units::Unknown;
337 m_dimUnits = dimUnits;
453 template<
typename IN_ARRAY >
459 if( m_interpolationMethod == TableFunction::InterpolationType::Linear )
461 return interpolateLinear( input );
465 return interpolateRound( input );
469 template<
typename IN_ARRAY >
473 TableFunction::KernelWrapper::interpolateLinear( IN_ARRAY
const & input )
const
482 arraySlice1d< real64 const >
const coords = m_coordinates[dim];
483 if( input[dim] <= coords[0] )
488 weights[dim][0] = 0.0;
489 weights[dim][1] = 1.0;
491 else if( input[dim] >= coords[coords.size() - 1] )
494 bounds[dim][0] = coords.size() - 1;
495 bounds[dim][1] = bounds[dim][0];
496 weights[dim][0] = 1.0;
497 weights[dim][1] = 0.0;
508 auto const lower = LvArray::sortedArrayManipulation::find( coords.begin(), coords.size(), input[dim] );
509 bounds[dim][1] = LvArray::integerConversion< localIndex >( lower );
510 bounds[dim][0] = bounds[dim][1] - 1;
512 real64 const dx = coords[bounds[dim][1]] - coords[bounds[dim][0]];
513 weights[dim][0] = 1.0 - ( input[dim] - coords[bounds[dim][0]]) / dx;
514 weights[dim][1] = 1.0 - weights[dim][0];
521 for(
integer point = 0; point < numCorners; ++point )
528 integer const corner = (point >> dim) & 1;
529 tableIndex += bounds[dim][corner] * stride;
530 stride *= m_coordinates.sizeOfArray( dim );
534 real64 cornerValue = m_values[tableIndex];
537 integer const corner = (point >> dim) & 1;
538 cornerValue *= weights[dim][corner];
540 value += cornerValue;
545 template<
typename IN_ARRAY >
549 TableFunction::KernelWrapper::interpolateRound( IN_ARRAY
const & input )
const
558 arraySlice1d< real64 const >
const coords = m_coordinates[dim];
561 if( input[dim] <= coords[0] )
566 else if( input[dim] >= coords[coords.size() - 1] )
569 subIndex = coords.size() - 1;
575 auto const lower = LvArray::sortedArrayManipulation::find( coords.begin(), coords.size(), input[dim] );
576 subIndex = LvArray::integerConversion< localIndex >( lower );
582 if( m_interpolationMethod == TableFunction::InterpolationType::Nearest )
584 if( ( input[dim] - coords[subIndex - 1]) <= ( coords[subIndex] - input[dim]) )
589 else if( m_interpolationMethod == TableFunction::InterpolationType::Lower )
599 tableIndex += subIndex * stride;
600 stride *= coords.size();
604 return m_values[tableIndex];
607 template<
typename IN_ARRAY >
611 TableFunction::KernelWrapper::getCoord( IN_ARRAY
const & input,
localIndex const dim,
InterpolationType interpolationMethod )
const
615 arraySlice1d< real64 const >
const coords = m_coordinates[dim];
617 if( input[dim] <= coords[0] )
622 else if( input[dim] >= coords[coords.size() - 1] )
625 subIndex = coords.size() - 1;
631 auto const lower = LvArray::sortedArrayManipulation::find( coords.begin(), coords.size(), input[dim] );
632 subIndex = LvArray::integerConversion< localIndex >( lower );
638 if( interpolationMethod == TableFunction::InterpolationType::Nearest )
640 if( ( input[dim] - coords[subIndex - 1]) <= ( coords[subIndex] - input[dim]) )
645 else if( interpolationMethod == TableFunction::InterpolationType::Lower )
655 return coords[subIndex];
658 template<
typename IN_ARRAY,
typename OUT_ARRAY >
665 if( m_interpolationMethod == TableFunction::InterpolationType::Linear )
667 return interpolateLinear( input, derivatives );
672 return interpolateRound( input, derivatives );
676 template<
typename IN_ARRAY,
typename OUT_ARRAY >
680 TableFunction::KernelWrapper::interpolateLinear( IN_ARRAY
const & input, OUT_ARRAY && derivatives )
const
691 arraySlice1d< real64 const >
const coords = m_coordinates[dim];
692 if( input[dim] <= coords[0] )
699 dWeights_dInput[dim][0] = 0;
700 dWeights_dInput[dim][1] = 0;
702 else if( input[dim] >= coords[coords.size() - 1] )
705 bounds[dim][0] = coords.size() - 1;
706 bounds[dim][1] = bounds[dim][0];
709 dWeights_dInput[dim][0] = 0;
710 dWeights_dInput[dim][1] = 0;
718 auto lower = LvArray::sortedArrayManipulation::find( coords.begin(), coords.size(), input[dim] );
719 bounds[dim][1] = LvArray::integerConversion< localIndex >( lower );
720 bounds[dim][0] = bounds[dim][1] - 1;
722 real64 const dx = coords[bounds[dim][1]] - coords[bounds[dim][0]];
723 weights[dim][0] = 1.0 - ( input[dim] - coords[bounds[dim][0]]) / dx;
724 weights[dim][1] = 1.0 - weights[dim][0];
725 dWeights_dInput[dim][0] = -1.0 / dx;
726 dWeights_dInput[dim][1] = -dWeights_dInput[dim][0];
734 derivatives[dim] = 0.0;
738 for(
integer point = 0; point < numCorners; ++point )
745 integer const corner = (point >> dim) & 1;
746 tableIndex += bounds[dim][corner] * stride;
747 stride *= m_coordinates.sizeOfArray( dim );
751 real64 cornerValue = m_values[tableIndex];
755 dCornerValue_dInput[dim] = cornerValue;
760 integer const corner = (point >> dim) & 1;
761 cornerValue *= weights[dim][corner];
764 dCornerValue_dInput[kk] *= ( dim == kk ) ? dWeights_dInput[dim][corner] : weights[dim][corner];
770 derivatives[dim] += dCornerValue_dInput[dim];
772 value += cornerValue;
777 template<
typename IN_ARRAY,
typename OUT_ARRAY >
781 TableFunction::KernelWrapper::interpolateRound( IN_ARRAY
const & input, OUT_ARRAY && derivatives )
const
784 GEOS_ERROR(
"Rounding interpolation with derivatives not implemented" );
803 string TableTextFormatter::toString< TableFunction >(
TableFunction const & tableData )
const;
811 string TableCSVFormatter::toString< TableFunction >(
TableFunction const & tableData )
const;
#define GEOS_HOST_DEVICE
Marks a host-device function.
#define GEOS_UNUSED_VAR(...)
Mark an unused variable and silence compiler warnings.
#define GEOS_FORCE_INLINE
Marks a function or lambda for inlining.
#define GEOS_ERROR(msg)
Raise a hard error and terminate the program.
Enumerates the Units that are in use in GEOS and regroups useful conversion and formatting functions.
Unit
Enumerator of available unit types for given physical scales. Units are in SI by default.
Class describing a file Path.
GEOS_HOST_DEVICE real64 compute(IN_ARRAY const &input, OUT_ARRAY &&derivatives) const
Interpolate in the table with derivatives.
void move(LvArray::MemorySpace const space, bool const touch)
Move the KernelWrapper to the given execution space, optionally touching it.
GEOS_HOST_DEVICE real64 compute(IN_ARRAY const &input) const
Interpolate in the table.
string getCoordsDescription(integer dimId, bool shortUnitsToVariables) const
InterpolationType
Enumerator of available interpolation types.
void setTableValues(real64_array values, units::Unit unit=units::Unknown)
Set the table values.
array1d< real64 > & getValues()
Get the table values.
arrayView1d< real64 const > getValues() const
Get the table values.
void checkCoord(real64 coord, localIndex dim) const
Check if the given coordinate is in the bounds of the table coordinates in the specified dimension,...
virtual void evaluate(dataRepository::Group const &group, real64 const time, SortedArrayView< localIndex const > const &set, arrayView1d< real64 > const &result) const override final
Method to evaluate a function on a target object.
ArrayOfArraysView< real64 > getCoordinates()
Get the table axes definitions.
real64 getCoord(real64 const *const input, localIndex dim, InterpolationType interpolationMethod) const
Method to get coordinates.
void setTableCoordinates(array1d< real64_array > const &coordinates, stdVector< units::Unit > const &dimUnits={})
Set the table coordinates.
void setValueUnits(units::Unit unit)
Set the table value units.
string getTableDescription() const
void reInitializeFunction()
Build the maps used to evaluate the table function.
KernelWrapper createKernelWrapper() const
Create an instance of the kernel wrapper.
void setDimUnits(stdVector< units::Unit > const &dimUnits)
Set the units of each dimension.
void initializePostSubGroups() override
Called by Initialize() after to initializing sub-Groups.
TableFunction(const string &name, dataRepository::Group *const parent)
The constructor.
units::Unit getValueUnit() const
void setInterpolationMethod(InterpolationType const method)
Set the interpolation method.
static constexpr integer maxDimensions
maximum dimensions for the coordinates in the table
virtual real64 evaluate(real64 const *const input) const override final
Method to evaluate a function.
static string catalogName()
The catalog name interface.
virtual void initializeFunction() override
Initialize the table function.
ArrayOfArraysView< real64 const > getCoordinates() const
Get the table axes definitions.
integer numDimensions() const
units::Unit getDimUnit(localIndex const dim) const
InterpolationType getInterpolationMethod() const
Get the interpolation method.
string getValuesDescription() const
void outputTableData(OutputOptions const outputOpts) const
Print the table(s) in the log and/or CSV files when requested by the user.
Group & operator=(Group const &)=delete
Deleted copy assignment operator.
ArrayView< T, 1 > arrayView1d
Alias for 1D array view.
LvArray::ArrayOfArraysView< T, INDEX_TYPE const, CONST_SIZES, LvArray::ChaiBuffer > ArrayOfArraysView
View of array of variable-sized arrays. See LvArray::ArrayOfArraysView for details.
array1d< Path > path_array
A 1-dimensional array of geos::Path types.
std::set< T > set
A set of local indices.
double real64
64-bit floating point type.
GEOS_LOCALINDEX_TYPE localIndex
Local index type (for indexing objects within an MPI partition).
array1d< real64 > real64_array
A 1-dimensional array of geos::real64 types.
std::size_t size_t
Unsigned size type.
LvArray::SortedArrayView< T, localIndex, LvArray::ChaiBuffer > SortedArrayView
A sorted array view of local indices.
int integer
Signed integer type.
Array< T, 1 > array1d
Alias for 1D array.
LvArray::ArrayOfArrays< T, INDEX_TYPE, LvArray::ChaiBuffer > ArrayOfArrays
Array of variable-sized arrays. See LvArray::ArrayOfArrays for details.
ENUM_STRINGS(LinearSolverParameters::SolverType, "direct", "cg", "gmres", "fgmres", "bicgstab", "richardson", "preconditioner")
Declare strings associated with enumeration values.
internal::StdVectorWrapper< T, Allocator, USE_STD_CONTAINER_BOUNDS_CHECKING > stdVector
Struct containing output options.
bool writeInLog
Request table output in log.
bool writeCSV
Request table output in CSV file.
Struct containing lookup keys for data repository wrappers.
static constexpr char const * coordinatesString()
static constexpr char const * interpolationString()
static constexpr char const * valuesString()
static constexpr char const * voxelFileString()
static constexpr char const * coordinateFilesString()
static constexpr char const * writeCSVFlagString()