ug4
|
Global assembler based on the ghost-fluid method with a level-set function. More...
#include <dom_disc_embb.h>
Public Types | |
typedef TAlgebra | algebra_type |
Algebra type. More... | |
typedef ApproximationSpace< domain_type > | approx_space_type |
type of approximation space More... | |
typedef TDomain | domain_type |
Domain type. More... | |
typedef TExtrapolation | extrapolation_type |
Extrapolation type. More... | |
typedef extrapolation_type::ls_grid_func_type | ls_grid_func_type |
Grid function type for the LSF. More... | |
typedef algebra_type::matrix_type | matrix_type |
Matrix type in the algebra. More... | |
typedef algebra_type::vector_type | vector_type |
Vector type in the algebra. More... | |
Public Member Functions | |
template<typename TElem , typename TIterator > | |
void | AssembleDefect (const std::vector< IElemDisc< domain_type > * > &vElemDisc, ConstSmartPtr< domain_type > spDomain, ConstSmartPtr< DoFDistribution > dd, TIterator iterBegin, TIterator iterEnd, int si, bool bNonRegularGrid, vector_type &d, const vector_type &u, ConstSmartPtr< AssemblingTuner< TAlgebra > > spAssTuner) |
template<typename TElem , typename TIterator > | |
void | AssembleDefect (const std::vector< IElemDisc< domain_type > * > &vElemDisc, ConstSmartPtr< domain_type > spDomain, ConstSmartPtr< DoFDistribution > dd, TIterator iterBegin, TIterator iterEnd, int si, bool bNonRegularGrid, vector_type &d, ConstSmartPtr< VectorTimeSeries< vector_type > > vSol, const std::vector< number > &vScaleMass, const std::vector< number > &vScaleStiff, ConstSmartPtr< AssemblingTuner< TAlgebra > > spAssTuner) |
template<typename TElem , typename TIterator > | |
void | AssembleJacobian (const std::vector< IElemDisc< domain_type > * > &vElemDisc, ConstSmartPtr< domain_type > spDomain, ConstSmartPtr< DoFDistribution > dd, TIterator iterBegin, TIterator iterEnd, int si, bool bNonRegularGrid, matrix_type &J, const vector_type &u, ConstSmartPtr< AssemblingTuner< TAlgebra > > spAssTuner) |
template<typename TElem , typename TIterator > | |
void | AssembleJacobian (const std::vector< IElemDisc< domain_type > * > &vElemDisc, ConstSmartPtr< domain_type > spDomain, ConstSmartPtr< DoFDistribution > dd, TIterator iterBegin, TIterator iterEnd, int si, bool bNonRegularGrid, matrix_type &J, ConstSmartPtr< VectorTimeSeries< vector_type > > vSol, number s_a0, ConstSmartPtr< AssemblingTuner< TAlgebra > > spAssTuner) |
template<typename TElem , typename TIterator > | |
void | AssembleLinear (const std::vector< IElemDisc< domain_type > * > &vElemDisc, ConstSmartPtr< domain_type > spDomain, ConstSmartPtr< DoFDistribution > dd, TIterator iterBegin, TIterator iterEnd, int si, bool bNonRegularGrid, matrix_type &A, vector_type &rhs, ConstSmartPtr< AssemblingTuner< TAlgebra > > spAssTuner) |
template<typename TElem , typename TIterator > | |
void | AssembleLinear (const std::vector< IElemDisc< domain_type > * > &vElemDisc, ConstSmartPtr< domain_type > spDomain, ConstSmartPtr< DoFDistribution > dd, TIterator iterBegin, TIterator iterEnd, int si, bool bNonRegularGrid, matrix_type &A, vector_type &rhs, ConstSmartPtr< VectorTimeSeries< vector_type > > vSol, const std::vector< number > &vScaleMass, const std::vector< number > &vScaleStiff, ConstSmartPtr< AssemblingTuner< TAlgebra > > spAssTuner) |
template<typename TElem , typename TIterator > | |
void | AssembleMassMatrix (const std::vector< IElemDisc< domain_type > * > &vElemDisc, ConstSmartPtr< domain_type > spDomain, ConstSmartPtr< DoFDistribution > dd, TIterator iterBegin, TIterator iterEnd, int si, bool bNonRegularGrid, matrix_type &M, const vector_type &u, ConstSmartPtr< AssemblingTuner< TAlgebra > > spAssTuner) |
template<typename TElem , typename TIterator > | |
void | AssembleStiffnessMatrix (const std::vector< IElemDisc< domain_type > * > &vElemDisc, ConstSmartPtr< domain_type > spDomain, ConstSmartPtr< DoFDistribution > dd, TIterator iterBegin, TIterator iterEnd, int si, bool bNonRegularGrid, matrix_type &A, const vector_type &u, ConstSmartPtr< AssemblingTuner< TAlgebra > > spAssTuner) |
int | check_elem_lsf (size_t n_co, GridObject *pElem, int si, int g_level, bool use_hanging, const MathVector< dim > vCornerCoords[], number time) |
checks whether the element is intersected by the interface, or what, and prepares the data More... | |
virtual void | clear_outer_values (vector_type &d, const DoFDistribution *dd) const |
sets the values at the outer vertices to 0 More... | |
bool | corner_inside (size_t co) const |
returns true if the corner is "inside" (use after check_elem_lsf) More... | |
void | exclude_subsets (SmartPtr< approx_space_type > spApproxSpace, const char *subset_names) |
excludes a (boundary) subsets from the extrapolation More... | |
void | extrapolate_by_lsf (size_t num_co, number *u, size_t fct) const |
extrapolates a component the solution to the vertices behind the interface (by averaging) More... | |
void | extrapolate_by_lsf (size_t num_co, size_t base_co, number *u, size_t fct) const |
extrapolates a component the solution to the vertices behind the interface (w.r.t. a base corner) More... | |
extrapolation_type & | extrapolation () |
returns the extrapolation More... | |
void | FinishTimestep (const std::vector< IElemDisc< domain_type > * > &vElemDisc, ConstSmartPtr< DoFDistribution > dd, bool bNonRegularGrid, ConstSmartPtr< VectorTimeSeries< vector_type > > vSol, ConstSmartPtr< AssemblingTuner< TAlgebra > > spAssTuner) |
template<typename TElem , typename TIterator > | |
void | FinishTimestepElem (const std::vector< IElemDisc< domain_type > * > &vElemDisc, ConstSmartPtr< domain_type > spDomain, ConstSmartPtr< DoFDistribution > dd, TIterator iterBegin, TIterator iterEnd, int si, bool bNonRegularGrid, ConstSmartPtr< VectorTimeSeries< vector_type > > vSol, ConstSmartPtr< AssemblingTuner< TAlgebra > > spAssTuner) |
template<typename TElem , typename TIterator > | |
void | InitAllExports (const std::vector< IElemDisc< domain_type > * > &vElemDisc, ConstSmartPtr< DoFDistribution > dd, TIterator iterBegin, TIterator iterEnd, int si, bool bNonRegularGrid, bool bAsTimeDependent) |
number | lsf_at (size_t co) const |
returns the effective value of the LSF at a corner (use after check_elem_lsf) More... | |
LSGFGlobAssembler () | |
class constructor (may not have any arguments!) More... | |
void | prepare_interface_bc (SmartPtr< approx_space_type > spApproxSpace) |
prepares the boundary conditions at the interface: sets all them to Dirichlet-0 More... | |
void | PrepareTimestep (const std::vector< IElemDisc< domain_type > * > &vElemDisc, ConstSmartPtr< DoFDistribution > dd, bool bNonRegularGrid, ConstSmartPtr< VectorTimeSeries< vector_type > > vSol, number future_time, ConstSmartPtr< AssemblingTuner< TAlgebra > > spAssTuner) |
template<typename TElem , typename TIterator > | |
void | PrepareTimestepElem (const std::vector< IElemDisc< domain_type > * > &vElemDisc, ConstSmartPtr< domain_type > spDomain, ConstSmartPtr< DoFDistribution > dd, TIterator iterBegin, TIterator iterEnd, int si, bool bNonRegularGrid, ConstSmartPtr< VectorTimeSeries< vector_type > > vSol, ConstSmartPtr< AssemblingTuner< TAlgebra > > spAssTuner) |
void | project_LSF () |
project the level-set function to the coarse levels More... | |
void | set_assemble_only_cut (bool b) |
set the "assemble only in cut elements" flag More... | |
void | set_Dirichlet_on_if_for (SmartPtr< approx_space_type > spApproxSpace, const char *fct_name, number value) |
adds a Dirichlet BC with a given value on the interface More... | |
void | set_Dirichlet_on_if_for (SmartPtr< approx_space_type > spApproxSpace, const char *fct_name, SmartPtr< CplUserData< number, dim > > func) |
adds a Dirichlet BC with a given value on the interface More... | |
void | set_LSF (SmartPtr< ls_grid_func_type > spLSF) |
set the level-set function and check it More... | |
void | set_Neumann0_on_if_for (SmartPtr< approx_space_type > spApproxSpace, const char *fct_name) |
adds a Neumann-0 with on the interface More... | |
virtual void | set_outer_matrices (matrix_type &A, const DoFDistribution *dd) const |
sets the matrices at outer vertices to identity More... | |
virtual void | set_outer_values (vector_type &u, const DoFDistribution *dd, number time) |
sets the values at the outer vertices to given values More... | |
void | set_plain_Dirichlet_on_if_for (SmartPtr< approx_space_type > spApproxSpace, const char *fct_name, number value) |
adds a "plain" Dirichlet BC with a given value on the interface More... | |
void | set_plain_Dirichlet_on_if_for (SmartPtr< approx_space_type > spApproxSpace, const char *fct_name, SmartPtr< CplUserData< number, dim > > func) |
adds a "plain" Dirichlet BC with a given value on the interface More... | |
virtual | ~LSGFGlobAssembler () |
virtual destructor More... | |
Static Public Member Functions | |
template<typename TElem , typename TIterator > | |
static void | AssembleErrorEstimator (const std::vector< IElemError< domain_type > * > &vElemDisc, ConstSmartPtr< domain_type > spDomain, ConstSmartPtr< DoFDistribution > dd, TIterator iterBegin, TIterator iterEnd, int si, bool bNonRegularGrid, const vector_type &u) |
template<typename TElem , typename TIterator > | |
static void | AssembleErrorEstimator (const std::vector< IElemError< domain_type > * > &vElemDisc, ConstSmartPtr< domain_type > spDomain, ConstSmartPtr< DoFDistribution > dd, TIterator iterBegin, TIterator iterEnd, int si, bool bNonRegularGrid, std::vector< number > vScaleMass, std::vector< number > vScaleStiff, ConstSmartPtr< VectorTimeSeries< vector_type > > vSol) |
template<typename TElem , typename TIterator > | |
static void | AssembleRhs (const std::vector< IElemDisc< domain_type > * > &vElemDisc, ConstSmartPtr< domain_type > spDomain, ConstSmartPtr< DoFDistribution > dd, TIterator iterBegin, TIterator iterEnd, int si, bool bNonRegularGrid, vector_type &rhs, const vector_type &u, ConstSmartPtr< AssemblingTuner< TAlgebra > > spAssTuner) |
template<typename TElem , typename TIterator > | |
static void | AssembleRhs (const std::vector< IElemDisc< domain_type > * > &vElemDisc, ConstSmartPtr< domain_type > spDomain, ConstSmartPtr< DoFDistribution > dd, TIterator iterBegin, TIterator iterEnd, int si, bool bNonRegularGrid, vector_type &rhs, ConstSmartPtr< VectorTimeSeries< vector_type > > vSol, const std::vector< number > &vScaleMass, const std::vector< number > &vScaleStiff, ConstSmartPtr< AssemblingTuner< TAlgebra > > spAssTuner) |
Static Public Attributes | |
static const int | dim = TDomain::dim |
world dimension More... | |
Private Attributes | |
bool | m_bAssembleOnlyCut |
extrapolation_type | m_extrapol |
the extrapolation at the interface More... | |
Global assembler based on the ghost-fluid method with a level-set function.
Template class of the global assembler based on the ghost-fluid method with a piecewise linear level-set function.
TDomain | domain type |
TAlgebra | algebra type |
TExtrapolation | extrapolation class |
typedef TAlgebra ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::algebra_type |
Algebra type.
typedef ApproximationSpace<domain_type> ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::approx_space_type |
type of approximation space
typedef TDomain ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::domain_type |
Domain type.
typedef TExtrapolation ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::extrapolation_type |
Extrapolation type.
typedef extrapolation_type::ls_grid_func_type ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::ls_grid_func_type |
Grid function type for the LSF.
typedef algebra_type::matrix_type ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::matrix_type |
Matrix type in the algebra.
typedef algebra_type::vector_type ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::vector_type |
Vector type in the algebra.
|
inline |
class constructor (may not have any arguments!)
|
inlinevirtual |
virtual destructor
void ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::AssembleDefect | ( | const std::vector< IElemDisc< domain_type > * > & | vElemDisc, |
ConstSmartPtr< domain_type > | spDomain, | ||
ConstSmartPtr< DoFDistribution > | dd, | ||
TIterator | iterBegin, | ||
TIterator | iterEnd, | ||
int | si, | ||
bool | bNonRegularGrid, | ||
vector_type & | d, | ||
const vector_type & | u, | ||
ConstSmartPtr< AssemblingTuner< TAlgebra > > | spAssTuner | ||
) |
This function adds the contributions of all passed element discretizations on one given subset to the global Defect in the stationary case. (This version processes elements in a given interval.)
[in] | vElemDisc | element discretizations |
[in] | spDomain | domain |
[in] | dd | DoF Distribution |
[in] | iterBegin | element iterator |
[in] | iterEnd | element iterator |
[in] | si | subset index |
[in] | bNonRegularGrid | flag to indicate if non regular grid is used |
[in,out] | d | defect |
[in] | u | solution |
[in] | spAssTuner | assemble adapter |
References ug::bridge::AlgebraCommon::Eval(), ug::FillCornerCoordinates(), ug::GetLocalVector(), ug::hex_rules::NUM_VERTICES, ug::LocalVector::resize(), ug::RHS, ug::LocalVector::scale_append(), ug::STIFF, ug::GridLevel::TOP, and UG_CATCH_THROW.
void ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::AssembleDefect | ( | const std::vector< IElemDisc< domain_type > * > & | vElemDisc, |
ConstSmartPtr< domain_type > | spDomain, | ||
ConstSmartPtr< DoFDistribution > | dd, | ||
TIterator | iterBegin, | ||
TIterator | iterEnd, | ||
int | si, | ||
bool | bNonRegularGrid, | ||
vector_type & | d, | ||
ConstSmartPtr< VectorTimeSeries< vector_type > > | vSol, | ||
const std::vector< number > & | vScaleMass, | ||
const std::vector< number > & | vScaleStiff, | ||
ConstSmartPtr< AssemblingTuner< TAlgebra > > | spAssTuner | ||
) |
This function adds the contributions of all passed element discretizations on one given subset to the global Defect in the instationary case. (This version processes elements in a given interval.)
[in] | vElemDisc | element discretizations |
[in] | spDomain | domain |
[in] | dd | DoF Distribution |
[in] | iterBegin | element iterator |
[in] | iterEnd | element iterator |
[in] | si | subset index |
[in] | bNonRegularGrid | flag to indicate if non regular grid is used |
[in,out] | d | defect |
[in] | vSol | current and previous solutions |
[in] | vScaleMass | scaling factors for mass part |
[in] | vScaleStiff | scaling factors for stiffness part |
[in] | spAssTuner | assemble adapter |
References ug::bridge::AlgebraCommon::Eval(), ug::EXPL, ug::FillCornerCoordinates(), ug::MASS, ug::hex_rules::NUM_VERTICES, ug::PT_INSTATIONARY, ug::PT_STATIONARY, ug::LocalVectorTimeSeries::read_times(), ug::LocalVectorTimeSeries::read_values(), ug::LocalVector::resize(), ug::RHS, ug::LocalVector::scale_append(), ug::LocalVectorTimeSeries::solution(), ug::STIFF, ug::GridLevel::TOP, UG_CATCH_THROW, and UG_THROW.
|
inlinestatic |
References UG_THROW.
|
inlinestatic |
References UG_THROW.
void ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::AssembleJacobian | ( | const std::vector< IElemDisc< domain_type > * > & | vElemDisc, |
ConstSmartPtr< domain_type > | spDomain, | ||
ConstSmartPtr< DoFDistribution > | dd, | ||
TIterator | iterBegin, | ||
TIterator | iterEnd, | ||
int | si, | ||
bool | bNonRegularGrid, | ||
matrix_type & | J, | ||
const vector_type & | u, | ||
ConstSmartPtr< AssemblingTuner< TAlgebra > > | spAssTuner | ||
) |
This function adds the contributions of all passed element discretizations on one given subset to the global Jacobian in the stationary case. (This version processes elements in a given interval.)
[in] | vElemDisc | element discretizations |
[in] | spDomain | domain |
[in] | iterBegin | element iterator |
[in] | iterEnd | element iterator |
[in] | si | subset index |
[in] | bNonRegularGrid | flag to indicate if non regular grid is used |
[in,out] | J | jacobian |
[in] | u | solution |
[in] | spAssTuner | assemble adapter |
References ug::bridge::AlgebraCommon::Eval(), ug::FillCornerCoordinates(), ug::GetLocalVector(), ug::hex_rules::NUM_VERTICES, ug::LocalVector::resize(), ug::LocalMatrix::resize(), ug::RHS, ug::STIFF, ug::GridLevel::TOP, and UG_CATCH_THROW.
void ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::AssembleJacobian | ( | const std::vector< IElemDisc< domain_type > * > & | vElemDisc, |
ConstSmartPtr< domain_type > | spDomain, | ||
ConstSmartPtr< DoFDistribution > | dd, | ||
TIterator | iterBegin, | ||
TIterator | iterEnd, | ||
int | si, | ||
bool | bNonRegularGrid, | ||
matrix_type & | J, | ||
ConstSmartPtr< VectorTimeSeries< vector_type > > | vSol, | ||
number | s_a0, | ||
ConstSmartPtr< AssemblingTuner< TAlgebra > > | spAssTuner | ||
) |
This function adds the contributions of all passed element discretizations on one given subset to the global Jacobian in the time-dependent case. Note, that s_m0 == 1 (This version processes elements in a given interval.)
[in] | vElemDisc | element discretizations |
[in] | spDomain | domain |
[in] | dd | DoF Distribution |
[in] | iterBegin | element iterator |
[in] | iterEnd | element iterator |
[in] | si | subset index |
[in] | bNonRegularGrid | flag to indicate if non regular grid is used |
[in,out] | J | jacobian |
[in] | vSol | current and previous solutions |
[in] | s_a0 | scaling factor for stiffness part |
[in] | spAssTuner | assemble adapter |
References ug::bridge::AlgebraCommon::Eval(), ug::FillCornerCoordinates(), ug::GetLocalVector(), ug::MASS, ug::hex_rules::NUM_VERTICES, ug::PT_INSTATIONARY, ug::PT_STATIONARY, ug::LocalVectorTimeSeries::read_times(), ug::LocalVectorTimeSeries::read_values(), ug::LocalVector::resize(), ug::LocalMatrix::resize(), ug::RHS, ug::LocalVectorTimeSeries::size(), ug::LocalVectorTimeSeries::solution(), ug::STIFF, ug::GridLevel::TOP, and UG_CATCH_THROW.
void ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::AssembleLinear | ( | const std::vector< IElemDisc< domain_type > * > & | vElemDisc, |
ConstSmartPtr< domain_type > | spDomain, | ||
ConstSmartPtr< DoFDistribution > | dd, | ||
TIterator | iterBegin, | ||
TIterator | iterEnd, | ||
int | si, | ||
bool | bNonRegularGrid, | ||
matrix_type & | A, | ||
vector_type & | rhs, | ||
ConstSmartPtr< AssemblingTuner< TAlgebra > > | spAssTuner | ||
) |
This function adds the contributions of all passed element discretizations on one given subset to the global Matrix and the global Right-Hand Side of the Linear problem in the stationary case. (This version processes elements in a given interval.)
[in] | vElemDisc | element discretizations |
[in] | spDomain | domain |
[in] | dd | DoF Distribution |
[in] | iterBegin | element iterator |
[in] | iterEnd | element iterator |
[in] | si | subset index |
[in] | bNonRegularGrid | flag to indicate if non regular grid is used |
[in,out] | A | Matrix |
[in,out] | rhs | Right-hand side |
[in] | spAssTuner | assemble adapter |
References ug::bridge::AlgebraCommon::Eval(), ug::FillCornerCoordinates(), ug::hex_rules::NUM_VERTICES, ug::LocalVector::resize(), ug::LocalMatrix::resize(), ug::RHS, ug::STIFF, ug::GridLevel::TOP, and UG_CATCH_THROW.
void ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::AssembleLinear | ( | const std::vector< IElemDisc< domain_type > * > & | vElemDisc, |
ConstSmartPtr< domain_type > | spDomain, | ||
ConstSmartPtr< DoFDistribution > | dd, | ||
TIterator | iterBegin, | ||
TIterator | iterEnd, | ||
int | si, | ||
bool | bNonRegularGrid, | ||
matrix_type & | A, | ||
vector_type & | rhs, | ||
ConstSmartPtr< VectorTimeSeries< vector_type > > | vSol, | ||
const std::vector< number > & | vScaleMass, | ||
const std::vector< number > & | vScaleStiff, | ||
ConstSmartPtr< AssemblingTuner< TAlgebra > > | spAssTuner | ||
) |
This function adds the contributions of all passed element discretizations on one given subset to the global Matrix and the global Right-Hand Side of the Linear problem in the time-dependent case. (This version processes elements in a given interval.)
[in] | vElemDisc | element discretizations |
[in] | spDomain | domain |
[in] | dd | DoF Distribution |
[in] | iterBegin | element iterator |
[in] | iterEnd | element iterator |
[in] | si | subset index |
[in] | bNonRegularGrid | flag to indicate if non regular grid is used |
[in,out] | A | Matrix |
[in,out] | rhs | Right-hand side |
[in] | vSol | current and previous solutions |
[in] | vScaleMass | scaling factors for mass part |
[in] | vScaleStiff | scaling factors for stiffness part |
[in] | spAssTuner | assemble adapter |
References ug::bridge::AlgebraCommon::Eval(), ug::FillCornerCoordinates(), ug::MASS, ug::hex_rules::NUM_VERTICES, ug::PT_INSTATIONARY, ug::PT_STATIONARY, ug::LocalVectorTimeSeries::read_times(), ug::LocalVectorTimeSeries::read_values(), ug::LocalVector::resize(), ug::LocalMatrix::resize(), ug::RHS, ug::LocalVector::scale_append(), ug::LocalMatrix::scale_append(), ug::LocalVectorTimeSeries::solution(), ug::STIFF, ug::GridLevel::TOP, UG_CATCH_THROW, and UG_THROW.
void ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::AssembleMassMatrix | ( | const std::vector< IElemDisc< domain_type > * > & | vElemDisc, |
ConstSmartPtr< domain_type > | spDomain, | ||
ConstSmartPtr< DoFDistribution > | dd, | ||
TIterator | iterBegin, | ||
TIterator | iterEnd, | ||
int | si, | ||
bool | bNonRegularGrid, | ||
matrix_type & | M, | ||
const vector_type & | u, | ||
ConstSmartPtr< AssemblingTuner< TAlgebra > > | spAssTuner | ||
) |
This function adds the contributions of all passed element discretizations on one given subset to the global Mass matrix. (This version processes elements in a given interval.)
[in] | vElemDisc | element discretizations |
[in] | spDomain | domain |
[in] | dd | DoF Distribution |
[in] | iterBegin | element iterator |
[in] | iterEnd | element iterator |
[in] | si | subset index |
[in] | bNonRegularGrid | flag to indicate if non regular grid is used |
[in,out] | M | Mass matrix |
[in] | u | solution |
[in] | spAssTuner | assemble adapter |
References ug::bridge::AlgebraCommon::Eval(), ug::FillCornerCoordinates(), ug::GetLocalVector(), ug::MASS, ug::hex_rules::NUM_VERTICES, ug::LocalVector::resize(), ug::LocalMatrix::resize(), ug::GridLevel::TOP, and UG_CATCH_THROW.
|
inlinestatic |
References UG_THROW.
|
inlinestatic |
References UG_THROW.
void ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::AssembleStiffnessMatrix | ( | const std::vector< IElemDisc< domain_type > * > & | vElemDisc, |
ConstSmartPtr< domain_type > | spDomain, | ||
ConstSmartPtr< DoFDistribution > | dd, | ||
TIterator | iterBegin, | ||
TIterator | iterEnd, | ||
int | si, | ||
bool | bNonRegularGrid, | ||
matrix_type & | A, | ||
const vector_type & | u, | ||
ConstSmartPtr< AssemblingTuner< TAlgebra > > | spAssTuner | ||
) |
This function adds the contributions of all passed element discretizations on one given subset to the global Stiffness matrix. (This version processes elements in a given interval.)
[in] | vElemDisc | element discretizations |
[in] | spDomain | domain |
[in] | dd | DoF Distribution |
[in] | iterBegin | element iterator |
[in] | iterEnd | element iterator |
[in] | si | subset index |
[in] | bNonRegularGrid | flag to indicate if non regular grid is used |
[in,out] | A | Stiffness matrix |
[in] | u | solution |
[in] | spAssTuner | assemble adapter |
References ug::bridge::AlgebraCommon::Eval(), ug::FillCornerCoordinates(), ug::GetLocalVector(), ug::hex_rules::NUM_VERTICES, ug::LocalVector::resize(), ug::LocalMatrix::resize(), ug::STIFF, ug::GridLevel::TOP, and UG_CATCH_THROW.
|
inline |
checks whether the element is intersected by the interface, or what, and prepares the data
n_co | number of the corners of the element |
pElem | the element to process |
si | subset of the element |
g_level | grid level of the element |
use_hanging | if there can be hanging nodes |
vCornerCoords | coordinates of the corners of the element |
time | the phisical time |
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
Referenced by ug::LSGFDomainDiscretization< TDomain, TAlgebra, TExtrapolation >::check_elem_lsf().
|
inlinevirtual |
sets the values at the outer vertices to 0
d | the vector where to set |
dd | dof distribution of the grid function to reset |
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
|
inline |
returns true if the corner is "inside" (use after check_elem_lsf)
co | the corner |
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
Referenced by ug::LSGFDomainDiscretization< TDomain, TAlgebra, TExtrapolation >::corner_inside().
|
inline |
excludes a (boundary) subsets from the extrapolation
spApproxSpace | the approximation space of the domain discretization |
subset_names | names of the subsets to exclude |
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
Referenced by ug::LSGFDomainDiscretization< TDomain, TAlgebra, TExtrapolation >::exclude_subsets().
|
inline |
extrapolates a component the solution to the vertices behind the interface (by averaging)
num_co | number of the corners |
u | nodal values to extrapolate |
fct | index of the function (to identify to type of the extrapolation) |
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
|
inline |
extrapolates a component the solution to the vertices behind the interface (w.r.t. a base corner)
num_co | number of the corners |
base_co | the base corner |
u | nodal values to extrapolate |
fct | index of the function (to identify to type of the extrapolation) |
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
Referenced by ug::LSGFDomainDiscretization< TDomain, TAlgebra, TExtrapolation >::extrapolate_by_lsf().
|
inline |
returns the extrapolation
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
void ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::FinishTimestep | ( | const std::vector< IElemDisc< domain_type > * > & | vElemDisc, |
ConstSmartPtr< DoFDistribution > | dd, | ||
bool | bNonRegularGrid, | ||
ConstSmartPtr< VectorTimeSeries< vector_type > > | vSol, | ||
ConstSmartPtr< AssemblingTuner< TAlgebra > > | spAssTuner | ||
) |
This function finishes the global discretization for a time-stepping scheme by calling the "finish_timestep" methods of all passed element discretizations.
[in] | vElemDisc | element discretizations |
[in] | dd | DoF Distribution |
[in] | bNonRegularGrid | flag to indicate if non regular grid is used |
[in] | vSol | current and previous solutions |
[in] | spAssTuner | assemble adapter |
References ug::bridge::AlgebraCommon::Eval(), ug::bridge::AlgebraTypeIDProvider::id(), ug::bridge::AlgebraTypeIDProvider::instance(), ug::NONE, ug::LocalVectorTimeSeries::read_times(), and UG_CATCH_THROW.
void ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::FinishTimestepElem | ( | const std::vector< IElemDisc< domain_type > * > & | vElemDisc, |
ConstSmartPtr< domain_type > | spDomain, | ||
ConstSmartPtr< DoFDistribution > | dd, | ||
TIterator | iterBegin, | ||
TIterator | iterEnd, | ||
int | si, | ||
bool | bNonRegularGrid, | ||
ConstSmartPtr< VectorTimeSeries< vector_type > > | vSol, | ||
ConstSmartPtr< AssemblingTuner< TAlgebra > > | spAssTuner | ||
) |
This function finalizes the global discretization in a time-stepping scheme by calling the "finish_timestep_elem" methods of all passed element discretizations on one given subset. (This version processes elements in a given interval.)
REMARK: This version does not extrapolate anything. It merely skips outer elements.
[in] | vElemDisc | element discretizations |
[in] | dd | DoF Distribution |
[in] | iterBegin | element iterator |
[in] | iterEnd | element iterator |
[in] | si | subset index |
[in] | bNonRegularGrid | flag to indicate if non regular grid is used |
[in] | vSol | current and previous solutions |
[in] | spAssTuner | assemble adapter |
References ug::bridge::AlgebraCommon::Eval(), ug::FillCornerCoordinates(), ug::GetLocalVector(), ug::NONE, ug::hex_rules::NUM_VERTICES, ug::LocalVectorTimeSeries::read_times(), ug::LocalVectorTimeSeries::read_values(), ug::LocalVector::resize(), ug::LocalVectorTimeSeries::time(), ug::GridLevel::TOP, and UG_CATCH_THROW.
void ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::InitAllExports | ( | const std::vector< IElemDisc< domain_type > * > & | vElemDisc, |
ConstSmartPtr< DoFDistribution > | dd, | ||
TIterator | iterBegin, | ||
TIterator | iterEnd, | ||
int | si, | ||
bool | bNonRegularGrid, | ||
bool | bAsTimeDependent | ||
) |
This function finalizes the global discretization in a time-stepping scheme by calling the "finish_timestep_elem" methods of all passed element discretizations on one given subset. (This version processes elements in a given interval.)
[in] | vElemDisc | element discretizations |
[in] | dd | DoF Distribution |
[in] | iterBegin | element iterator |
[in] | iterEnd | element iterator |
[in] | si | subset index |
[in] | bNonRegularGrid | flag to indicate if non regular grid is used |
[in] | bAsTimeDependent | flag to simulate the instationary case for the initialization |
References ug::bridge::AlgebraCommon::Eval(), ug::MASS, ug::RHS, ug::STIFF, and UG_CATCH_THROW.
|
inline |
returns the effective value of the LSF at a corner (use after check_elem_lsf)
co | the corner |
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
Referenced by ug::LSGFDomainDiscretization< TDomain, TAlgebra, TExtrapolation >::lsf_at().
|
inline |
prepares the boundary conditions at the interface: sets all them to Dirichlet-0
spApproxSpace | the approximation space of the domain discretization |
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
Referenced by ug::LSGFDomainDiscretization< TDomain, TAlgebra, TExtrapolation >::LSGFDomainDiscretization().
void ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::PrepareTimestep | ( | const std::vector< IElemDisc< domain_type > * > & | vElemDisc, |
ConstSmartPtr< DoFDistribution > | dd, | ||
bool | bNonRegularGrid, | ||
ConstSmartPtr< VectorTimeSeries< vector_type > > | vSol, | ||
number | future_time, | ||
ConstSmartPtr< AssemblingTuner< TAlgebra > > | spAssTuner | ||
) |
This function prepares the global discretization for a time-stepping scheme by calling the "prepare_timestep" methods of all passed element discretizations.
[in] | vElemDisc | element discretizations |
[in] | dd | DoF Distribution |
[in] | bNonRegularGrid | flag to indicate if non regular grid is used |
[in] | vSol | current and previous solutions |
[in] | spAssTuner | assemble adapter |
References ug::bridge::AlgebraCommon::Eval(), ug::bridge::AlgebraTypeIDProvider::id(), ug::bridge::AlgebraTypeIDProvider::instance(), ug::NONE, ug::LocalVectorTimeSeries::read_times(), and UG_CATCH_THROW.
void ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::PrepareTimestepElem | ( | const std::vector< IElemDisc< domain_type > * > & | vElemDisc, |
ConstSmartPtr< domain_type > | spDomain, | ||
ConstSmartPtr< DoFDistribution > | dd, | ||
TIterator | iterBegin, | ||
TIterator | iterEnd, | ||
int | si, | ||
bool | bNonRegularGrid, | ||
ConstSmartPtr< VectorTimeSeries< vector_type > > | vSol, | ||
ConstSmartPtr< AssemblingTuner< TAlgebra > > | spAssTuner | ||
) |
This function prepares the global discretization for a time-stepping scheme by calling the "prepare_timestep_elem" methods of all passed element discretizations on one given subset. (This version processes elements in a given interval.)
REMARK: This version does not extrapolate anything. It merely skips outer elements.
[in] | vElemDisc | element discretizations |
[in] | spDomain | domain |
[in] | dd | DoF Distribution |
[in] | iterBegin | element iterator |
[in] | iterEnd | element iterator |
[in] | si | subset index |
[in] | bNonRegularGrid | flag to indicate if non regular grid is used |
[in] | vSol | current and previous solutions |
[in] | spAssTuner | assemble adapter |
References ug::bridge::AlgebraCommon::Eval(), ug::FillCornerCoordinates(), ug::GetLocalVector(), ug::NONE, ug::hex_rules::NUM_VERTICES, ug::LocalVectorTimeSeries::read_times(), ug::LocalVectorTimeSeries::read_values(), ug::LocalVector::resize(), ug::GridLevel::TOP, and UG_CATCH_THROW.
|
inline |
project the level-set function to the coarse levels
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
Referenced by ug::LSGFDomainDiscretization< TDomain, TAlgebra, TExtrapolation >::project_LSF().
|
inline |
set the "assemble only in cut elements" flag
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_bAssembleOnlyCut.
Referenced by ug::LSGFDomainDiscretization< TDomain, TAlgebra, TExtrapolation >::set_assemble_only_cut().
|
inline |
adds a Dirichlet BC with a given value on the interface
spApproxSpace | the approximation space of the domain discretization |
fct_name | function to impose the condition for |
value | the Dirichlet value |
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
Referenced by ug::LSGFDomainDiscretization< TDomain, TAlgebra, TExtrapolation >::set_Dirichlet_on_if_for().
|
inline |
adds a Dirichlet BC with a given value on the interface
spApproxSpace | the approximation space of the domain discretization |
fct_name | function to impose the condition for |
func | the Dirichlet function |
References func(), and ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
|
inline |
set the level-set function and check it
spLSF | the function to set |
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
Referenced by ug::LSGFDomainDiscretization< TDomain, TAlgebra, TExtrapolation >::set_LSF().
|
inline |
adds a Neumann-0 with on the interface
spApproxSpace | the approximation space of the domain discretization |
fct_name | function to impose the condition for |
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
Referenced by ug::LSGFDomainDiscretization< TDomain, TAlgebra, TExtrapolation >::set_Neumann0_on_if_for().
|
inlinevirtual |
sets the matrices at outer vertices to identity
A | the vector where to set |
dd | dof distribution of the grid function to reset |
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
|
inlinevirtual |
sets the values at the outer vertices to given values
u | the vector where to set |
dd | dof distribution of the grid function to reset |
time | the physical time |
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
|
inline |
adds a "plain" Dirichlet BC with a given value on the interface
spApproxSpace | the approximation space of the domain discretization |
fct_name | function to impose the condition for |
value | the Dirichlet value |
References ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
|
inline |
adds a "plain" Dirichlet BC with a given value on the interface
spApproxSpace | the approximation space of the domain discretization |
fct_name | function to impose the condition for |
func | the Dirichlet function |
References func(), and ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::m_extrapol.
Referenced by ug::LSGFDomainDiscretization< TDomain, TAlgebra, TExtrapolation >::set_plain_Dirichlet_on_if_for().
|
static |
world dimension
|
private |
The following flag is used only for research purposes (measuring the volume sources/sinks at the embedded interface etc.). It switches off assembing in all inner elements (i.e. makes the procedures to assemble only in the cut elements:
Referenced by ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::set_assemble_only_cut().
|
private |
the extrapolation at the interface
Referenced by ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::check_elem_lsf(), ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::clear_outer_values(), ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::corner_inside(), ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::exclude_subsets(), ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::extrapolate_by_lsf(), ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::extrapolation(), ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::lsf_at(), ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::prepare_interface_bc(), ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::project_LSF(), ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::set_Dirichlet_on_if_for(), ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::set_LSF(), ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::set_Neumann0_on_if_for(), ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::set_outer_matrices(), ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::set_outer_values(), and ug::LSGFGlobAssembler< TDomain, TAlgebra, TExtrapolation >::set_plain_Dirichlet_on_if_for().