preCICE v3.1.2
Loading...
Searching...
No Matches
Public Member Functions | Private Types | Private Member Functions | Private Attributes | List of all members
precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T > Class Template Reference

Mapping with radial basis functions using the Petsc library to solve the resulting system. More...

#include <PetRadialBasisFctMapping.hpp>

Inheritance diagram for precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >:
[legend]
Collaboration diagram for precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >:
[legend]

Public Member Functions

 PetRadialBasisFctMapping (Mapping::Constraint constraint, int dimensions, const RADIAL_BASIS_FUNCTION_T &function, std::array< bool, 3 > deadAxis, double solverRtol=1e-9, Polynomial polynomial=Polynomial::SEPARATE)
 Constructor.
 
 ~PetRadialBasisFctMapping () override
 Deletes the PETSc objects and the _deadAxis array.
 
void computeMapping () final override
 Computes the mapping coefficients from the in- and output mesh.
 
void clear () final override
 Removes a computed mapping.
 
std::string getName () const final override
 name of the rbf mapping
 
- Public Member Functions inherited from precice::mapping::RadialBasisFctBaseMapping< RADIAL_BASIS_FUNCTION_T >
 RadialBasisFctBaseMapping (Constraint constraint, int dimensions, const RADIAL_BASIS_FUNCTION_T &function, std::array< bool, 3 > deadAxis, InitialGuessRequirement mappingType)
 Constructor.
 
virtual ~RadialBasisFctBaseMapping ()=default
 
virtual void tagMeshFirstRound () final
 Method used by partition. Tags vertices that could be owned by this rank.
 
virtual void tagMeshSecondRound () final
 Method used by partition. Tags vertices that can be filtered out.
 
- Public Member Functions inherited from precice::mapping::Mapping
 Mapping (Constraint constraint, int dimensions, bool requiresGradientData, InitialGuessRequirement initialGuessRequirement)
 Constructor, takes mapping constraint.
 
Mappingoperator= (Mapping &&)=delete
 
virtual ~Mapping ()=default
 Destructor, empty.
 
void setMeshes (const mesh::PtrMesh &input, const mesh::PtrMesh &output)
 Sets input and output meshes carrying data to be mapped.
 
const mesh::PtrMeshgetInputMesh () const
 
const mesh::PtrMeshgetOutputMesh () const
 
Constraint getConstraint () const
 Returns the constraint (consistent/conservative) of the mapping.
 
MeshRequirement getInputRequirement () const
 Returns the requirement on the input mesh.
 
MeshRequirement getOutputRequirement () const
 Returns the requirement on the output mesh.
 
bool hasComputedMapping () const
 Returns true, if the mapping has been computed.
 
virtual bool hasConstraint (const Constraint &constraint) const
 Checks whether the mapping has the given constraint or not.
 
bool isScaledConsistent () const
 Returns true if mapping is a form of scaled consistent mapping.
 
bool requiresInitialGuess () const
 Return true if the mapping requires an initial guess.
 
const Eigen::VectorXd & initialGuess () const
 Return the provided initial guess of a mapping using an initialGuess.
 
Eigen::VectorXd & initialGuess ()
 
bool hasInitialGuess () const
 True if initialGuess().size() == 0.
 
void map (int inputDataID, int outputDataID)
 
void map (int inputDataID, int outputDataID, Eigen::VectorXd &initialGuess)
 
void map (const time::Sample &input, Eigen::VectorXd &output)
 Maps an input Sample to output data from input mesh to output mesh.
 
void map (const time::Sample &input, Eigen::VectorXd &output, Eigen::VectorXd &initialGuess)
 Maps an input Sample to output data from input mesh to output mesh, given an initialGuess.
 
virtual void scaleConsistentMapping (const Eigen::VectorXd &input, Eigen::VectorXd &output, Constraint type) const
 Scales the consistently mapped output data such that the surface integral of the values on input mesh and output mesh are equal.
 
bool requiresGradientData () const
 Returns whether the mapping requires gradient data.
 

Private Types

using VertexData = std::vector<std::vector<std::pair<int, double>>>
 Stores col -> value for each row. Used to return the already computed values from the preconditioning.
 

Private Member Functions

void mapConservative (const time::Sample &inData, Eigen::VectorXd &outData) final override
 Maps data using a conservative constraint.
 
void mapConsistent (const time::Sample &inData, Eigen::VectorXd &outData) final override
 Maps data using a consistent constraint.
 
void printMappingInfo (int dim) const
 Prints an INFO about the current mapping.
 
VertexData bgPreallocationMatrixC (const mesh::PtrMesh inMesh)
 Preallocate matrix C and saves the coefficients using a boost::geometry spatial tree for neighbor search.
 
VertexData bgPreallocationMatrixA (const mesh::PtrMesh inMesh, const mesh::PtrMesh outMesh)
 
void loadInitialGuessForDim (int dimension, int allDimensions, petsc::Vector &destination)
 
void storeInitialGuessForDim (int dimension, int allDimensions, petsc::Vector &source)
 

Private Attributes

logging::Logger _log {"mapping::PetRadialBasisFctMapping"}
 
petsc::Matrix _matrixC
 Interpolation system matrix. Evaluated basis function on the input mesh.
 
petsc::Matrix _matrixQ
 Vandermonde Matrix for linear polynomial, constructed from vertices of the input mesh.
 
petsc::Matrix _matrixA
 Interpolation evaluation matrix. Evaluated basis function on the output mesh.
 
petsc::Matrix _matrixV
 Coordinates of the output mesh to evaluate the separated polynomial.
 
petsc::Vector oneInterpolant
 Interpolant of g(x) = 1 evaluated at output sites, used for rescaling.
 
petsc::KSPSolver _solver
 Used to solve matrixC for the RBF weighting factors.
 
petsc::KSPSolver _QRsolver
 Used to solve the under-determined system for the separated polynomial.
 
AO _AOmapping
 Maps globalIndex to local index.
 
const double _solverRtol
 
Polynomial _polynomial
 Toggles the use of the additional polynomial.
 
bool useRescaling = true
 Toggles use of rescaled basis functions, only active when Polynomial == SEPARATE.
 
size_t polyparams
 Number of coefficients for the integrated polynomial. Depends on dimension and number of dead dimensions.
 
size_t sepPolyparams
 Number of coefficients for the separated polynomial. Depends on dimension and number of dead dimensions.
 
size_t localPolyparams
 Equals to polyparams if rank == 0. Is 0 everywhere else.
 
utils::Parallel::CommStatePtr _commState
 The CommState used to communicate.
 

Additional Inherited Members

- Public Types inherited from precice::mapping::Mapping
enum  Constraint { CONSISTENT , CONSERVATIVE , SCALED_CONSISTENT_SURFACE , SCALED_CONSISTENT_VOLUME }
 Specifies additional constraints for a mapping. More...
 
enum class  MeshRequirement { UNDEFINED = 0 , VERTEX = 1 , FULL = 2 }
 Specifies requirements for the input and output meshes of a mapping. More...
 
enum class  InitialGuessRequirement : bool { Required = true , None = false }
 Specifies whether the mapping requires an initial guess. More...
 
- Protected Member Functions inherited from precice::mapping::RadialBasisFctBaseMapping< RADIAL_BASIS_FUNCTION_T >
int getPolynomialParameters () const
 Computes the number of polynomial degrees of freedom based on the problem dimension and the dead axis.
 
- Protected Member Functions inherited from precice::mapping::Mapping
mesh::PtrMesh input () const
 Returns pointer to input mesh.
 
mesh::PtrMesh output () const
 Returns pointer to output mesh.
 
void setInputRequirement (MeshRequirement requirement)
 Sets the mesh requirement for the input mesh.
 
void setOutputRequirement (MeshRequirement requirement)
 Sets the mesh requirement for the output mesh.
 
int getDimensions () const
 
- Protected Attributes inherited from precice::mapping::RadialBasisFctBaseMapping< RADIAL_BASIS_FUNCTION_T >
RADIAL_BASIS_FUNCTION_T _basisFunction
 Radial basis function type used in interpolation.
 
std::vector< bool > _deadAxis
 true if the mapping along some axis should be ignored
 
- Protected Attributes inherited from precice::mapping::Mapping
bool _hasComputedMapping = false
 Flag to indicate whether computeMapping() has been called.
 
bool _requiresGradientData
 Flag if gradient data is required for the mapping.
 

Detailed Description

template<typename RADIAL_BASIS_FUNCTION_T>
class precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >

Mapping with radial basis functions using the Petsc library to solve the resulting system.

With help of the input data points and values an interpolant is constructed. The interpolant is formed by a weighted sum of conditionally positive radial basis functions and a (low order) polynomial, and evaluated at the output data points.

The radial basis function type has to be given as template parameter.

Definition at line 52 of file PetRadialBasisFctMapping.hpp.

Member Typedef Documentation

◆ VertexData

template<typename RADIAL_BASIS_FUNCTION_T >
using precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::VertexData = std::vector<std::vector<std::pair<int, double>>>
private

Stores col -> value for each row. Used to return the already computed values from the preconditioning.

Definition at line 94 of file PetRadialBasisFctMapping.hpp.

Constructor & Destructor Documentation

◆ PetRadialBasisFctMapping()

template<typename RADIAL_BASIS_FUNCTION_T >
precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::PetRadialBasisFctMapping ( Mapping::Constraint constraint,
int dimensions,
const RADIAL_BASIS_FUNCTION_T & function,
std::array< bool, 3 > deadAxis,
double solverRtol = 1e-9,
Polynomial polynomial = Polynomial::SEPARATE )

Constructor.

Parameters
[in]constraintSpecifies mapping to be consistent or conservative.
[in]dimensionsDimensionality of the meshes
[in]functionRadial basis function used for mapping.
[in]xDead,yDead,zDeadDeactivates mapping along an axis
[in]solverRtolRelative tolerance for the linear solver.
[in]polynomialType of polynomial augmentation

For description on convergence testing and meaning of solverRtol see http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/KSP/KSPConvergedDefault.html#KSPConvergedDefault

Definition at line 170 of file PetRadialBasisFctMapping.hpp.

Here is the call graph for this function:

◆ ~PetRadialBasisFctMapping()

template<typename RADIAL_BASIS_FUNCTION_T >
precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::~PetRadialBasisFctMapping ( )
override

Deletes the PETSc objects and the _deadAxis array.

Definition at line 195 of file PetRadialBasisFctMapping.hpp.

Here is the call graph for this function:

Member Function Documentation

◆ bgPreallocationMatrixA()

template<typename RADIAL_BASIS_FUNCTION_T >
PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::VertexData precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::bgPreallocationMatrixA ( const mesh::PtrMesh inMesh,
const mesh::PtrMesh outMesh )
private

Definition at line 922 of file PetRadialBasisFctMapping.hpp.

Here is the call graph for this function:

◆ bgPreallocationMatrixC()

template<typename RADIAL_BASIS_FUNCTION_T >
PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::VertexData precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::bgPreallocationMatrixC ( const mesh::PtrMesh inMesh)
private

Preallocate matrix C and saves the coefficients using a boost::geometry spatial tree for neighbor search.

Definition at line 842 of file PetRadialBasisFctMapping.hpp.

Here is the call graph for this function:

◆ clear()

template<typename RADIAL_BASIS_FUNCTION_T >
void precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::clear ( )
finaloverridevirtual

Removes a computed mapping.

Implements precice::mapping::RadialBasisFctBaseMapping< RADIAL_BASIS_FUNCTION_T >.

Definition at line 505 of file PetRadialBasisFctMapping.hpp.

Here is the call graph for this function:

◆ computeMapping()

template<typename RADIAL_BASIS_FUNCTION_T >
void precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::computeMapping ( )
finaloverridevirtual

Computes the mapping coefficients from the in- and output mesh.

Implements precice::mapping::RadialBasisFctBaseMapping< RADIAL_BASIS_FUNCTION_T >.

Definition at line 201 of file PetRadialBasisFctMapping.hpp.

Here is the call graph for this function:

◆ getName()

template<typename RADIAL_BASIS_FUNCTION_T >
std::string precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::getName ( ) const
finaloverridevirtual

name of the rbf mapping

Implements precice::mapping::Mapping.

Definition at line 521 of file PetRadialBasisFctMapping.hpp.

◆ loadInitialGuessForDim()

template<typename RADIAL_BASIS_FUNCTION_T >
void precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::loadInitialGuessForDim ( int dimension,
int allDimensions,
petsc::Vector & destination )
private

load the initialGuess for a given dimension or allocates the storage for the first iteration

The initialGuess passed to Mapping::map contains one guess for each data dimension as the mapping is evaluated for each dimension. As an examples, a 3D vector thus contains 3 initialGuesses, one for each solver. This extracts it from the overall initialGuess.

Definition at line 527 of file PetRadialBasisFctMapping.hpp.

Here is the call graph for this function:

◆ mapConservative()

template<typename RADIAL_BASIS_FUNCTION_T >
void precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::mapConservative ( const time::Sample & inData,
Eigen::VectorXd & outData )
finaloverrideprivatevirtual

Maps data using a conservative constraint.

Parameters
[in]inputSample to map data from
[in]outputValues to map to

If requiresInitialGuess(), then the initial guess is available via initialGuess(). Provide a new initial guess by overwriting it. The mapping has full control over its size.

See also
For mappings requiring an initialGuess: initialGuess() hasInitialGuess()

Implements precice::mapping::Mapping.

Definition at line 691 of file PetRadialBasisFctMapping.hpp.

Here is the call graph for this function:

◆ mapConsistent()

template<typename RADIAL_BASIS_FUNCTION_T >
void precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::mapConsistent ( const time::Sample & inData,
Eigen::VectorXd & outData )
finaloverrideprivatevirtual

Maps data using a consistent constraint.

Parameters
[in]inputSample to map data from
[in]outputValues to map to

If requiresInitialGuess(), then the initial guess is available via initialGuess(). Provide a new initial guess by overwriting it. The mapping has full control over its size.

See also
For mappings requiring an initialGuess: initialGuess() hasInitialGuess()

Implements precice::mapping::Mapping.

Definition at line 564 of file PetRadialBasisFctMapping.hpp.

Here is the call graph for this function:

◆ printMappingInfo()

template<typename RADIAL_BASIS_FUNCTION_T >
void precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::printMappingInfo ( int dim) const
private

Prints an INFO about the current mapping.

Definition at line 821 of file PetRadialBasisFctMapping.hpp.

◆ storeInitialGuessForDim()

template<typename RADIAL_BASIS_FUNCTION_T >
void precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::storeInitialGuessForDim ( int dimension,
int allDimensions,
petsc::Vector & source )
private

stores the initialGuess for a given dimension

The initialGuess passed to Mapping::map contains one guess for each data dimension as the mapping is evaluated for each dimension. As an examples, a 3D vector thus contains 3 initialGuesses, one for each solver. This stores the initialGuess of a given dimension into the overall initialGuess.

Definition at line 548 of file PetRadialBasisFctMapping.hpp.

Here is the call graph for this function:

Member Data Documentation

◆ _AOmapping

template<typename RADIAL_BASIS_FUNCTION_T >
AO precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::_AOmapping
private

Maps globalIndex to local index.

Definition at line 120 of file PetRadialBasisFctMapping.hpp.

◆ _commState

template<typename RADIAL_BASIS_FUNCTION_T >
utils::Parallel::CommStatePtr precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::_commState
private

The CommState used to communicate.

Definition at line 143 of file PetRadialBasisFctMapping.hpp.

◆ _log

template<typename RADIAL_BASIS_FUNCTION_T >
logging::Logger precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::_log {"mapping::PetRadialBasisFctMapping"}
mutableprivate

Definition at line 96 of file PetRadialBasisFctMapping.hpp.

◆ _matrixA

template<typename RADIAL_BASIS_FUNCTION_T >
petsc::Matrix precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::_matrixA
private

Interpolation evaluation matrix. Evaluated basis function on the output mesh.

Definition at line 105 of file PetRadialBasisFctMapping.hpp.

◆ _matrixC

template<typename RADIAL_BASIS_FUNCTION_T >
petsc::Matrix precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::_matrixC
private

Interpolation system matrix. Evaluated basis function on the input mesh.

Definition at line 99 of file PetRadialBasisFctMapping.hpp.

◆ _matrixQ

template<typename RADIAL_BASIS_FUNCTION_T >
petsc::Matrix precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::_matrixQ
private

Vandermonde Matrix for linear polynomial, constructed from vertices of the input mesh.

Definition at line 102 of file PetRadialBasisFctMapping.hpp.

◆ _matrixV

template<typename RADIAL_BASIS_FUNCTION_T >
petsc::Matrix precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::_matrixV
private

Coordinates of the output mesh to evaluate the separated polynomial.

Definition at line 108 of file PetRadialBasisFctMapping.hpp.

◆ _polynomial

template<typename RADIAL_BASIS_FUNCTION_T >
Polynomial precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::_polynomial
private

Toggles the use of the additional polynomial.

Definition at line 125 of file PetRadialBasisFctMapping.hpp.

◆ _QRsolver

template<typename RADIAL_BASIS_FUNCTION_T >
petsc::KSPSolver precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::_QRsolver
private

Used to solve the under-determined system for the separated polynomial.

Definition at line 117 of file PetRadialBasisFctMapping.hpp.

◆ _solver

template<typename RADIAL_BASIS_FUNCTION_T >
petsc::KSPSolver precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::_solver
private

Used to solve matrixC for the RBF weighting factors.

Definition at line 114 of file PetRadialBasisFctMapping.hpp.

◆ _solverRtol

template<typename RADIAL_BASIS_FUNCTION_T >
const double precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::_solverRtol
private

Definition at line 122 of file PetRadialBasisFctMapping.hpp.

◆ localPolyparams

template<typename RADIAL_BASIS_FUNCTION_T >
size_t precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::localPolyparams
private

Equals to polyparams if rank == 0. Is 0 everywhere else.

Definition at line 137 of file PetRadialBasisFctMapping.hpp.

◆ oneInterpolant

template<typename RADIAL_BASIS_FUNCTION_T >
petsc::Vector precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::oneInterpolant
private

Interpolant of g(x) = 1 evaluated at output sites, used for rescaling.

Definition at line 111 of file PetRadialBasisFctMapping.hpp.

◆ polyparams

template<typename RADIAL_BASIS_FUNCTION_T >
size_t precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::polyparams
private

Number of coefficients for the integrated polynomial. Depends on dimension and number of dead dimensions.

Definition at line 131 of file PetRadialBasisFctMapping.hpp.

◆ sepPolyparams

template<typename RADIAL_BASIS_FUNCTION_T >
size_t precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::sepPolyparams
private

Number of coefficients for the separated polynomial. Depends on dimension and number of dead dimensions.

Definition at line 134 of file PetRadialBasisFctMapping.hpp.

◆ useRescaling

template<typename RADIAL_BASIS_FUNCTION_T >
bool precice::mapping::PetRadialBasisFctMapping< RADIAL_BASIS_FUNCTION_T >::useRescaling = true
private

Toggles use of rescaled basis functions, only active when Polynomial == SEPARATE.

Definition at line 128 of file PetRadialBasisFctMapping.hpp.


The documentation for this class was generated from the following file: