Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys > Class Template Referenceabstract

CoreDriverBase - Base class for core driver that is agnostic to the Entity_kind. More...

#include <coredriver.h>

Inheritance diagram for Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >:
Portage::CoreDriver< D, ONWHAT, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >

Public Member Functions

 CoreDriverBase ()=default
 
virtual ~CoreDriverBase ()=default
 
virtual Entity_kind onwhat ()=0
 
template<Entity_kind ONWHAT, template< int, Entity_kind, class, class > class Search>
Portage::vector< std::vector< int > > search ()
 search for candidate source entities whose control volumes (cells, dual cells) overlap the control volumes of target cells More...
 
template<Entity_kind ONWHAT, template< Entity_kind, class, class, class, template< class, int, class, class > class, class, class > class Intersect>
Portage::vector< std::vector< Portage::Weights_t > > intersect_meshes (Portage::vector< std::vector< int >> const &intersection_candidates)
 intersect target entities with candidate source entities More...
 
template<template< Entity_kind, class, class, class, template< class, int, class, class > class, class, class > class Intersect>
std::vector< Portage::vector< std::vector< Portage::Weights_t > > > intersect_materials (Portage::vector< std::vector< int >> const &intersection_candidates)
 
template<Entity_kind ONWHAT>
Portage::vector< Vector< D > > compute_source_gradient (std::string const field_name, Limiter_type limiter_type=NOLIMITER, Boundary_Limiter_type boundary_limiter_type=BND_NOLIMITER, int material_id=0, const Part< SourceMesh, SourceState > *source_part=nullptr)
 Compute the gradient field of the given variable on source mesh. More...
 
template<typename T = double, Entity_kind ONWHAT, template< int, Entity_kind, class, class, class, class, class, template< class, int, class, class > class, class, class, class > class Interpolate>
void interpolate_mesh_var (std::string srcvarname, std::string trgvarname, Portage::vector< std::vector< Weights_t >> const &sources_and_weights, Portage::vector< Vector< D >> *gradients=nullptr)
 
template<typename T = double, Entity_kind ONWHAT, template< int, Entity_kind, class, class, class, class, class, template< class, int, class, class > class, class, class, class > class Interpolate>
std::enable_if< ONWHAT==CELL, void > interpolate_mesh_var (std::string srcvarname, std::string trgvarname, Portage::vector< std::vector< Weights_t >> const &sources_and_weights, const PartPair< D, SourceMesh, SourceState, TargetMesh, TargetState > *parts_pair, Portage::vector< Vector< D >> *gradients=nullptr)
 
template<typename T = double, template< int, Entity_kind, class, class, class, class, class, template< class, int, class, class > class, class, class, class > class Interpolate>
void interpolate_mat_var (std::string srcvarname, std::string trgvarname, std::vector< Portage::vector< std::vector< Weights_t >>> const &sources_and_weights_by_mat, std::vector< Portage::vector< Vector< D >>> *gradients=nullptr)
 
template<Entity_kind ONWHAT>
bool check_mismatch (Portage::vector< std::vector< Weights_t >> const &source_weights)
 Check if meshes are mismatched (don't cover identical portions of space) More...
 
template<Entity_kind ONWHAT>
bool has_mismatch ()
 Return if meshes are mismatched (check_mismatch must already have been called) More...
 
template<Entity_kind ONWHAT>
bool fix_mismatch (std::string const &src_var_name, std::string const &trg_var_name, double global_lower_bound=-std::numeric_limits< double >::max(), double global_upper_bound=std::numeric_limits< double >::max(), double conservation_tol=1e2 *std::numeric_limits< double >::epsilon(), int maxiter=5, Partial_fixup_type partial_fixup_type=Partial_fixup_type::SHIFTED_CONSERVATIVE, Empty_fixup_type empty_fixup_type=Empty_fixup_type::EXTRAPOLATE)
 Return if meshes are mismatched (check_mismatch must already have been called) More...
 
template<Entity_kind ONWHAT>
void set_num_tols (const double min_absolute_distance, const double min_absolute_volume)
 Set numerical tolerances for small distances and volumes. More...
 
template<Entity_kind ONWHAT>
void set_num_tols (const NumericTolerances_t &num_tols)
 Set numerical tolerances for small volumes, distances, etc. More...
 

Detailed Description

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
class Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >

CoreDriverBase - Base class for core driver that is agnostic to the Entity_kind.

Template Parameters
SourceMeshA lightweight wrapper to a specific input mesh implementation that provides certain functionality.
SourceStateA lightweight wrapper to a specific input state manager implementation that provides certain functionality.
TargetMeshA lightweight wrapper to a specific target mesh implementation that provides certain functionality.
TargetStateA lightweight wrapper to a specific target state manager implementation that provides certain functionality.
InterfaceReconstructorTypeAn interface reconstruction class that takes the raw interface reconstruction method, the dimension of the problem and the source mesh class as template parameters
Matpoly_SplitterA polygon/polyhedron splitting class (returns both pieces of the polygon)
Matpoly_ClipperA polygon/polyhedron clipping class (returns only the piece below/behind the clipping plane)
CoordinateSystem

Constructor & Destructor Documentation

◆ CoreDriverBase()

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >::CoreDriverBase ( )
default

◆ ~CoreDriverBase()

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
virtual Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >::~CoreDriverBase ( )
virtualdefault

Member Function Documentation

◆ check_mismatch()

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
template<Entity_kind ONWHAT>
bool Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >::check_mismatch ( Portage::vector< std::vector< Weights_t >> const &  source_weights)
inline

Check if meshes are mismatched (don't cover identical portions of space)

Template Parameters
Entity_kindWhat kind of entity are we performing intersection of
Parameters
[in]sources_weightsIntersection sources and moments (vols, centroids)
Returns
Whether the meshes are mismatched

◆ compute_source_gradient()

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
template<Entity_kind ONWHAT>
Portage::vector<Vector<D> > Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >::compute_source_gradient ( std::string const  field_name,
Limiter_type  limiter_type = NOLIMITER,
Boundary_Limiter_type  boundary_limiter_type = BND_NOLIMITER,
int  material_id = 0,
const Part< SourceMesh, SourceState > *  source_part = nullptr 
)
inline

Compute the gradient field of the given variable on source mesh.

Template Parameters
ONWHATentity kind (cell or node).
Parameters
field_namethe variable name.
limiter_typegradient limiter to use on internal regions.
boundary_limiter_typegradient limiter to use on boundary.
source_partthe source mesh part to consider if any.

◆ fix_mismatch()

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
template<Entity_kind ONWHAT>
bool Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >::fix_mismatch ( std::string const &  src_var_name,
std::string const &  trg_var_name,
double  global_lower_bound = -std::numeric_limits<double>::max(),
double  global_upper_bound = std::numeric_limits<double>::max(),
double  conservation_tol = 1e2*std::numeric_limits<double>::epsilon(),
int  maxiter = 5,
Partial_fixup_type  partial_fixup_type = Partial_fixup_type::SHIFTED_CONSERVATIVE,
Empty_fixup_type  empty_fixup_type = Empty_fixup_type::EXTRAPOLATE 
)
inline

Return if meshes are mismatched (check_mismatch must already have been called)

Returns
Whether the meshes are mismatched

◆ has_mismatch()

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
template<Entity_kind ONWHAT>
bool Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >::has_mismatch ( )
inline

Return if meshes are mismatched (check_mismatch must already have been called)

Returns
Whether the meshes are mismatched

◆ interpolate_mat_var()

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
template<typename T = double, template< int, Entity_kind, class, class, class, class, class, template< class, int, class, class > class, class, class, class > class Interpolate>
void Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >::interpolate_mat_var ( std::string  srcvarname,
std::string  trgvarname,
std::vector< Portage::vector< std::vector< Weights_t >>> const &  sources_and_weights_by_mat,
std::vector< Portage::vector< Vector< D >>> *  gradients = nullptr 
)
inline

Interpolate a (multi-)material variable of type T residing on CELLs

Parameters
[in]srcvarnameVariable name on source mesh
[in]trgvarnameVariable name on target mesh
[in]gradientsGradients of variable on source mesh (can be nullptr for 1st order remap)

◆ interpolate_mesh_var() [1/2]

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
template<typename T = double, Entity_kind ONWHAT, template< int, Entity_kind, class, class, class, class, class, template< class, int, class, class > class, class, class, class > class Interpolate>
void Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >::interpolate_mesh_var ( std::string  srcvarname,
std::string  trgvarname,
Portage::vector< std::vector< Weights_t >> const &  sources_and_weights,
Portage::vector< Vector< D >> *  gradients = nullptr 
)
inline

Interpolate a mesh variable of type T residing on entity kind ONWHAT using previously computed intersection weights

Template Parameters
Ttype of variable
ONWHATEntity_kind that field resides on
InterpolateFunctor for doing the interpolate from mesh to mesh
Parameters
[in]srcvarnameVariable name on source mesh
[in]trgvarnameVariable name on target mesh
[in]gradientsGradients of variable on source mesh (can be nullptr for 1st order remap)

◆ interpolate_mesh_var() [2/2]

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
template<typename T = double, Entity_kind ONWHAT, template< int, Entity_kind, class, class, class, class, class, template< class, int, class, class > class, class, class, class > class Interpolate>
std::enable_if<ONWHAT == CELL, void> Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >::interpolate_mesh_var ( std::string  srcvarname,
std::string  trgvarname,
Portage::vector< std::vector< Weights_t >> const &  sources_and_weights,
const PartPair< D, SourceMesh, SourceState, TargetMesh, TargetState > *  parts_pair,
Portage::vector< Vector< D >> *  gradients = nullptr 
)
inline

(Part-by-part) Interpolate a cell variable of type T residing on entity kind ONWHAT using previously computed intersection weights

Template Parameters
Ttype of variable
ONWHATEntity_kind that field resides on (enabled only for CELLs)
InterpolateFunctor for doing the interpolate from mesh to mesh
Parameters
[in]srcvarnameVariable name on source mesh
[in]trgvarnameVariable name on target mesh
[in]parts_pairPair of parts between which we have to remap
[in]gradientsGradients of variable on source mesh (can be nullptr for 1st order remap)

Enabled only for cells using SFINAE

◆ intersect_materials()

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
template<template< Entity_kind, class, class, class, template< class, int, class, class > class, class, class > class Intersect>
std::vector<Portage::vector<std::vector<Portage::Weights_t> > > Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >::intersect_materials ( Portage::vector< std::vector< int >> const &  intersection_candidates)
inline

intersect target cells with source material polygons

Template Parameters
Intersectintersect functor
Parameters
[in]intersection_candidatesvector of intersection candidates for each target entity
Returns
material-wise vector of intersection moments for each target cell

◆ intersect_meshes()

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
template<Entity_kind ONWHAT, template< Entity_kind, class, class, class, template< class, int, class, class > class, class, class > class Intersect>
Portage::vector<std::vector<Portage::Weights_t> > Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >::intersect_meshes ( Portage::vector< std::vector< int >> const &  intersection_candidates)
inline

intersect target entities with candidate source entities

Template Parameters
Entity_kindwhat kind of entity are we searching on/for
Intersectintersect functor
Parameters
[in]intersection_candidatesvector of intersection candidates for each target entity
Returns
vector of intersection moments for each target entity

◆ onwhat()

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
virtual Entity_kind Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >::onwhat ( )
pure virtual

◆ search()

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
template<Entity_kind ONWHAT, template< int, Entity_kind, class, class > class Search>
Portage::vector<std::vector<int> > Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >::search ( )
inline

search for candidate source entities whose control volumes (cells, dual cells) overlap the control volumes of target cells

Template Parameters
Entity_kindwhat kind of entity are we searching on/for
Searchsearch functor
Returns
vector of intersection candidates for each target entity

◆ set_num_tols() [1/2]

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
template<Entity_kind ONWHAT>
void Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >::set_num_tols ( const double  min_absolute_distance,
const double  min_absolute_volume 
)
inline

Set numerical tolerances for small distances and volumes.

Parameters
Entity_kindwhat kind of entity are we setting for
min_absolute_distanceselected minimal distance
min_absolute_volumeselected minimal volume

◆ set_num_tols() [2/2]

template<int D, class SourceMesh, class SourceState, class TargetMesh, class TargetState, template< class, int, class, class > class InterfaceReconstructorType, class Matpoly_Splitter, class Matpoly_Clipper, class CoordSys>
template<Entity_kind ONWHAT>
void Portage::CoreDriverBase< D, SourceMesh, SourceState, TargetMesh, TargetState, InterfaceReconstructorType, Matpoly_Splitter, Matpoly_Clipper, CoordSys >::set_num_tols ( const NumericTolerances_t num_tols)
inline

Set numerical tolerances for small volumes, distances, etc.

Template Parameters
Entity_kindwhat kind of entity are we setting for
num_tolsstruct of selected numerical tolerances

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