ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim > Class Template Reference

Detailed Description

template<typename ShapeFunction, typename IntegrationMethod, int GlobalDim>
class ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >

Definition at line 76 of file RichardsComponentTransportFEM.h.

#include <RichardsComponentTransportFEM.h>

Inheritance diagram for ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >:

Collaboration diagram for ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >:

Public Member Functions
	LocalAssemblerData (MeshLib::Element const &element, std::size_t const local_matrix_size, bool is_axially_symmetric, unsigned const integration_order, RichardsComponentTransportProcessData const &process_data, ProcessVariable const &transport_process_variable)
void	assemble (double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_xdot, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data) override
std::vector< double > const &	getIntPtDarcyVelocity (const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
Eigen::Map< const Eigen::RowVectorXd >	getShapeMatrix (const unsigned integration_point) const override
	Provides the shape matrix at the given integration point. More...
std::vector< double > const &	getIntPtSaturation (const double t, std::vector< GlobalVector * > const &x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_table, std::vector< double > &cache) const override
Public Member Functions inherited from ProcessLib::LocalAssemblerInterface
virtual	~LocalAssemblerInterface ()=default
void	setInitialConditions (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, double const t, bool const use_monolithic_scheme, int const process_id)
virtual void	initialize (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table)
virtual void	preAssemble (double const, double const, std::vector< double > const &)
virtual void	assembleForStaggeredScheme (double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_xdot, int const process_id, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data)
virtual void	assembleWithJacobian (double const t, double const dt, std::vector< double > const &local_x, std::vector< double > const &local_xdot, const double dxdot_dx, const double dx_dx, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
virtual void	assembleWithJacobianForStaggeredScheme (double const t, double const dt, Eigen::VectorXd const &local_x, Eigen::VectorXd const &local_xdot, const double dxdot_dx, const double dx_dx, int const process_id, std::vector< double > &local_M_data, std::vector< double > &local_K_data, std::vector< double > &local_b_data, std::vector< double > &local_Jac_data)
virtual void	computeSecondaryVariable (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, double const t, double const dt, std::vector< GlobalVector * > const &x, GlobalVector const &x_dot, int const process_id)
virtual void	preTimestep (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, double const t, double const delta_t)
virtual void	postTimestep (std::size_t const mesh_item_id, std::vector< NumLib::LocalToGlobalIndexMap const * > const &dof_tables, std::vector< GlobalVector * > const &x, double const t, double const dt)
void	postNonLinearSolver (std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const &dof_table, GlobalVector const &x, GlobalVector const &xdot, double const t, double const dt, bool const use_monolithic_scheme, int const process_id)
virtual std::vector< double >	interpolateNodalValuesToIntegrationPoints (std::vector< double > const &)
virtual Eigen::Vector3d	getFlux (MathLib::Point3d const &, double const, std::vector< double > const &) const
virtual Eigen::Vector3d	getFlux (MathLib::Point3d const &, double const, std::vector< std::vector< double >> const &) const
	Fits to staggered scheme. More...
Public Member Functions inherited from NumLib::ExtrapolatableElement
virtual	~ExtrapolatableElement ()=default

Private Types
using	ShapeMatricesType = ShapeMatrixPolicyType< ShapeFunction, GlobalDim >
using	ShapeMatrices = typename ShapeMatricesType::ShapeMatrices
using	LocalMatrixType = typename ShapeMatricesType::template MatrixType< NUM_NODAL_DOF *ShapeFunction::NPOINTS, NUM_NODAL_DOF *ShapeFunction::NPOINTS >
using	LocalVectorType = typename ShapeMatricesType::template VectorType< NUM_NODAL_DOF *ShapeFunction::NPOINTS >
using	NodalVectorType = typename ShapeMatricesType::NodalVectorType
using	NodalRowVectorType = typename ShapeMatricesType::NodalRowVectorType
using	GlobalDimVectorType = typename ShapeMatricesType::GlobalDimVectorType
using	GlobalDimNodalMatrixType = typename ShapeMatricesType::GlobalDimNodalMatrixType
using	NodalMatrixType = typename ShapeMatricesType::NodalMatrixType
using	GlobalDimMatrixType = typename ShapeMatricesType::GlobalDimMatrixType

Private Attributes
unsigned const	_element_id
RichardsComponentTransportProcessData const &	_process_data
IntegrationMethod const	_integration_method
ProcessVariable const &	_transport_process_variable
std::vector< IntegrationPointData< NodalRowVectorType, GlobalDimNodalMatrixType, NodalMatrixType >, Eigen::aligned_allocator< IntegrationPointData< NodalRowVectorType, GlobalDimNodalMatrixType, NodalMatrixType > > >	_ip_data

Static Private Attributes
static const int	concentration_index = 0
static const int	concentration_size = ShapeFunction::NPOINTS
static const int	pressure_index = ShapeFunction::NPOINTS
static const int	pressure_size = ShapeFunction::NPOINTS

Member Typedef Documentation

GlobalDimMatrixType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

using ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::GlobalDimMatrixType = typename ShapeMatricesType::GlobalDimMatrixType

private

Definition at line 96 of file RichardsComponentTransportFEM.h.

GlobalDimNodalMatrixType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

using ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::GlobalDimNodalMatrixType = typename ShapeMatricesType::GlobalDimNodalMatrixType

private

Definition at line 93 of file RichardsComponentTransportFEM.h.

GlobalDimVectorType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

using ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::GlobalDimVectorType = typename ShapeMatricesType::GlobalDimVectorType

private

Definition at line 92 of file RichardsComponentTransportFEM.h.

LocalMatrixType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

using ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::LocalMatrixType = typename ShapeMatricesType::template MatrixType< NUM_NODAL_DOF * ShapeFunction::NPOINTS, NUM_NODAL_DOF * ShapeFunction::NPOINTS>

private

Definition at line 82 of file RichardsComponentTransportFEM.h.

LocalVectorType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

using ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::LocalVectorType = typename ShapeMatricesType::template VectorType<NUM_NODAL_DOF * ShapeFunction::NPOINTS>

private

Definition at line 85 of file RichardsComponentTransportFEM.h.

NodalMatrixType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

using ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::NodalMatrixType = typename ShapeMatricesType::NodalMatrixType

private

Definition at line 95 of file RichardsComponentTransportFEM.h.

NodalRowVectorType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

using ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::NodalRowVectorType = typename ShapeMatricesType::NodalRowVectorType

private

Definition at line 90 of file RichardsComponentTransportFEM.h.

NodalVectorType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

using ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::NodalVectorType = typename ShapeMatricesType::NodalVectorType

private

Definition at line 89 of file RichardsComponentTransportFEM.h.

ShapeMatrices

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

using ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::ShapeMatrices = typename ShapeMatricesType::ShapeMatrices

private

Definition at line 80 of file RichardsComponentTransportFEM.h.

ShapeMatricesType

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

using ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::ShapeMatricesType = ShapeMatrixPolicyType<ShapeFunction, GlobalDim>

private

Definition at line 79 of file RichardsComponentTransportFEM.h.

Constructor & Destructor Documentation

LocalAssemblerData()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::LocalAssemblerData	(	MeshLib::Element const &	element,
		std::size_t const	local_matrix_size,
		bool	is_axially_symmetric,
		unsigned const	integration_order,
		RichardsComponentTransportProcessData const &	process_data,
		ProcessVariable const &	transport_process_variable
	)

Definition at line 23 of file RichardsComponentTransportFEM-impl.h.

    : _element_id(element.getID()),
      _process_data(process_data),
      _integration_method(integration_order),
      _transport_process_variable(transport_process_variable)
{
    // This assertion is valid only if all nodal d.o.f. use the same shape
    // matrices.
    assert(local_matrix_size == ShapeFunction::NPOINTS * NUM_NODAL_DOF);
    (void)local_matrix_size;
 
    unsigned const n_integration_points =
        _integration_method.getNumberOfPoints();
    _ip_data.reserve(n_integration_points);
 
    auto const shape_matrices =
        NumLib::initShapeMatrices<ShapeFunction, ShapeMatricesType, GlobalDim>(
            element, is_axially_symmetric, _integration_method);
 
    for (unsigned ip = 0; ip < n_integration_points; ip++)
    {
        auto const& sm = shape_matrices[ip];
        double const integration_factor = sm.integralMeasure * sm.detJ;
        _ip_data.emplace_back(
            sm.N, sm.dNdx,
            _integration_method.getWeightedPoint(ip).getWeight() *
                integration_factor,
            sm.N.transpose() * sm.N * integration_factor *
                _integration_method.getWeightedPoint(ip).getWeight());
    }
}

References ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_integration_method, ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_ip_data, and ProcessLib::RichardsComponentTransport::NUM_NODAL_DOF.

Member Function Documentation

assemble()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

void ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::assemble ( double const  t, double const  dt, std::vector< double > const &  local_x, std::vector< double > const &  local_xdot, std::vector< double > &  local_M_data, std::vector< double > &  local_K_data, std::vector< double > &  local_b_data  )

overridevirtual

Reimplemented from ProcessLib::LocalAssemblerInterface.

Definition at line 63 of file RichardsComponentTransportFEM-impl.h.

{
    auto const local_matrix_size = local_x.size();
    // This assertion is valid only if all nodal d.o.f. use the same shape
    // matrices.
    assert(local_matrix_size == ShapeFunction::NPOINTS * NUM_NODAL_DOF);
 
    auto local_M = MathLib::createZeroedMatrix<LocalMatrixType>(
        local_M_data, local_matrix_size, local_matrix_size);
    auto local_K = MathLib::createZeroedMatrix<LocalMatrixType>(
        local_K_data, local_matrix_size, local_matrix_size);
    auto local_b = MathLib::createZeroedVector<LocalVectorType>(
        local_b_data, local_matrix_size);
 
    unsigned const n_integration_points =
        _integration_method.getNumberOfPoints();
 
    ParameterLib::SpatialPosition pos;
    pos.setElementID(_element_id);
 
    auto p_nodal_values = Eigen::Map<const NodalVectorType>(
        &local_x[pressure_index], pressure_size);
 
    auto const& b = _process_data.specific_body_force;
 
    MaterialPropertyLib::VariableArray vars;
 
    GlobalDimMatrixType const& I(
        GlobalDimMatrixType::Identity(GlobalDim, GlobalDim));
 
    // Get material properties
    auto const& medium = *_process_data.media_map->getMedium(_element_id);
    auto const& phase = medium.phase("AqueousLiquid");
    auto const& component =
        phase.component(_transport_process_variable.getName());
 
    auto KCC = local_K.template block<concentration_size, concentration_size>(
        concentration_index, concentration_index);
    auto MCC = local_M.template block<concentration_size, concentration_size>(
        concentration_index, concentration_index);
    auto Kpp = local_K.template block<pressure_size, pressure_size>(
        pressure_index, pressure_index);
    auto Mpp = local_M.template block<pressure_size, pressure_size>(
        pressure_index, pressure_index);
    auto Bp = local_b.template block<pressure_size, 1>(pressure_index, 0);
 
    for (std::size_t ip(0); ip < n_integration_points; ++ip)
    {
        pos.setIntegrationPoint(ip);
 
        auto const& ip_data = _ip_data[ip];
        auto const& N = ip_data.N;
        auto const& dNdx = ip_data.dNdx;
        auto const& w = ip_data.integration_weight;
 
        double C_int_pt = 0.0;
        double p_int_pt = 0.0;
        // Order matters: First C, then p!
        NumLib::shapeFunctionInterpolate(local_x, N, C_int_pt, p_int_pt);
 
        vars[static_cast<int>(
            MaterialPropertyLib::Variable::capillary_pressure)] = -p_int_pt;
        auto const Sw = medium[MaterialPropertyLib::PropertyType::saturation]
                            .template value<double>(vars, pos, t, dt);
 
        double const dSw_dpc =
            medium[MaterialPropertyLib::PropertyType::saturation]
                .template dValue<double>(
                    vars, MaterialPropertyLib::Variable::capillary_pressure,
                    pos, t, dt);
 
        vars[static_cast<int>(MaterialPropertyLib::Variable::concentration)] =
            C_int_pt;
        vars[static_cast<int>(MaterialPropertyLib::Variable::phase_pressure)] =
            p_int_pt;
 
        // \todo the argument to getValue() has to be changed for non
        // constant storage model
        auto const specific_storage =
            medium[MaterialPropertyLib::PropertyType::storage]
                .template value<double>(vars, pos, t, dt);
        // \todo the first argument has to be changed for non constant
        // porosity model
        auto const porosity =
            medium[MaterialPropertyLib::PropertyType::porosity]
                .template value<double>(vars, pos, t, dt);
 
        auto const retardation_factor =
            component[MaterialPropertyLib::PropertyType::retardation_factor]
                .template value<double>(vars, pos, t, dt);
 
        auto const solute_dispersivity_transverse =
            medium[MaterialPropertyLib::PropertyType::transversal_dispersivity]
                .template value<double>(vars, pos, t, dt);
        auto const solute_dispersivity_longitudinal =
            medium[MaterialPropertyLib::PropertyType::longitudinal_dispersivity]
                .template value<double>(vars, pos, t, dt);
 
        // Use the fluid density model to compute the density
        auto const density = phase[MaterialPropertyLib::PropertyType::density]
                                 .template value<double>(vars, pos, t, dt);
        auto const decay_rate =
            component[MaterialPropertyLib::PropertyType::decay_rate]
                .template value<double>(vars, pos, t, dt);
        auto const pore_diffusion_coefficient =
            MaterialPropertyLib::formEigenTensor<GlobalDim>(
                component[MaterialPropertyLib::PropertyType::pore_diffusion]
                    .value(vars, pos, t, dt));
 
        auto const K = MaterialPropertyLib::formEigenTensor<GlobalDim>(
            medium[MaterialPropertyLib::PropertyType::permeability].value(
                vars, pos, t, dt));
        vars[static_cast<int>(
            MaterialPropertyLib::Variable::liquid_saturation)] = Sw;
        auto const k_rel =
            medium[MaterialPropertyLib::PropertyType::relative_permeability]
                .template value<double>(vars, pos, t, dt);
        // Use the viscosity model to compute the viscosity
        auto const mu = phase[MaterialPropertyLib::PropertyType::viscosity]
                            .template value<double>(vars, pos, t, dt);
        auto const K_times_k_rel_over_mu = K * (k_rel / mu);
 
        GlobalDimVectorType const velocity =
            _process_data.has_gravity
                ? GlobalDimVectorType(-K_times_k_rel_over_mu *
                                      (dNdx * p_nodal_values - density * b))
                : GlobalDimVectorType(-K_times_k_rel_over_mu * dNdx *
                                      p_nodal_values);
 
        double const velocity_magnitude = velocity.norm();
        GlobalDimMatrixType const hydrodynamic_dispersion =
            velocity_magnitude != 0.0
                ? GlobalDimMatrixType(
                      porosity * pore_diffusion_coefficient +
                      solute_dispersivity_transverse * velocity_magnitude * I +
                      (solute_dispersivity_longitudinal -
                       solute_dispersivity_transverse) /
                          velocity_magnitude * velocity * velocity.transpose())
                : GlobalDimMatrixType(porosity * pore_diffusion_coefficient +
                                      solute_dispersivity_transverse *
                                          velocity_magnitude * I);
 
        // matrix assembly
        KCC.noalias() +=
            (dNdx.transpose() * hydrodynamic_dispersion * dNdx +
             N.transpose() * velocity.transpose() * dNdx +
             N.transpose() * decay_rate * porosity * retardation_factor * N) *
            w;
        MCC.noalias() += w * N.transpose() * porosity * retardation_factor * N;
        Kpp.noalias() += w * dNdx.transpose() * K_times_k_rel_over_mu * dNdx;
        // \TODO Extend to pressure dependent density.
        double const drhow_dp(0.0);
        Mpp.noalias() += (specific_storage * Sw + porosity * Sw * drhow_dp -
                          porosity * dSw_dpc) *
                         ip_data.mass_operator;
 
        if (_process_data.has_gravity)
        {
            Bp += w * density * dNdx.transpose() * K_times_k_rel_over_mu * b;
        }
        /* with Oberbeck-Boussing assumption density difference only exists
         * in buoyancy effects */
    }
}

References MaterialPropertyLib::capillary_pressure, MaterialPropertyLib::concentration, MaterialPropertyLib::decay_rate, MaterialPropertyLib::density, MaterialPropertyLib::liquid_saturation, MaterialPropertyLib::longitudinal_dispersivity, ProcessLib::RichardsComponentTransport::NUM_NODAL_DOF, MaterialPropertyLib::permeability, MaterialPropertyLib::phase_pressure, MaterialPropertyLib::pore_diffusion, MaterialPropertyLib::porosity, MaterialPropertyLib::relative_permeability, MaterialPropertyLib::retardation_factor, MaterialPropertyLib::saturation, ParameterLib::SpatialPosition::setElementID(), ParameterLib::SpatialPosition::setIntegrationPoint(), NumLib::shapeFunctionInterpolate(), MaterialPropertyLib::storage, MaterialPropertyLib::transversal_dispersivity, and MaterialPropertyLib::viscosity.

getIntPtDarcyVelocity()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

std::vector< double > const & ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::getIntPtDarcyVelocity ( const double  t, std::vector< GlobalVector * > const &  x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &  dof_table, std::vector< double > &  cache  ) const

overridevirtual

Implements ProcessLib::RichardsComponentTransport::RichardsComponentTransportLocalAssemblerInterface.

Definition at line 234 of file RichardsComponentTransportFEM-impl.h.

{
    unsigned const n_integration_points =
        _integration_method.getNumberOfPoints();
 
    constexpr int process_id = 0;  // monolithic scheme
    auto const indices =
        NumLib::getIndices(_element_id, *dof_table[process_id]);
    assert(!indices.empty());
    auto const local_x = x[process_id]->get(indices);
 
    cache.clear();
    auto cache_mat = MathLib::createZeroedMatrix<
        Eigen::Matrix<double, GlobalDim, Eigen::Dynamic, Eigen::RowMajor>>(
        cache, GlobalDim, n_integration_points);
 
    ParameterLib::SpatialPosition pos;
    pos.setElementID(_element_id);
 
    MaterialPropertyLib::VariableArray vars;
 
    // Get material properties
    auto const& medium = *_process_data.media_map->getMedium(_element_id);
    auto const& phase = medium.phase("AqueousLiquid");
 
    auto const p_nodal_values = Eigen::Map<const NodalVectorType>(
        &local_x[ShapeFunction::NPOINTS], ShapeFunction::NPOINTS);
 
    for (unsigned ip = 0; ip < n_integration_points; ++ip)
    {
        auto const& ip_data = _ip_data[ip];
        auto const& N = ip_data.N;
        auto const& dNdx = ip_data.dNdx;
 
        pos.setIntegrationPoint(ip);
 
        auto const dt = std::numeric_limits<double>::quiet_NaN();
        auto const K = MaterialPropertyLib::formEigenTensor<GlobalDim>(
            medium[MaterialPropertyLib::PropertyType::permeability].value(
                vars, pos, t, dt));
        auto const mu = phase[MaterialPropertyLib::PropertyType::viscosity]
                            .template value<double>(vars, pos, t, dt);
 
        double C_int_pt = 0.0;
        double p_int_pt = 0.0;
        NumLib::shapeFunctionInterpolate(local_x, N, C_int_pt, p_int_pt);
 
        // saturation
        vars[static_cast<int>(
            MaterialPropertyLib::Variable::capillary_pressure)] = -p_int_pt;
        auto const Sw = medium[MaterialPropertyLib::PropertyType::saturation]
                            .template value<double>(vars, pos, t, dt);
 
        vars[static_cast<int>(
            MaterialPropertyLib::Variable::liquid_saturation)] = Sw;
        auto const k_rel =
            medium[MaterialPropertyLib::PropertyType::relative_permeability]
                .template value<double>(vars, pos, t, dt);
 
        cache_mat.col(ip).noalias() = -dNdx * p_nodal_values;
        if (_process_data.has_gravity)
        {
            vars[static_cast<int>(
                MaterialPropertyLib::Variable::concentration)] = C_int_pt;
            vars[static_cast<int>(
                MaterialPropertyLib::Variable::phase_pressure)] = p_int_pt;
            auto const rho_w = phase[MaterialPropertyLib::PropertyType::density]
                                   .template value<double>(vars, pos, t, dt);
            auto const b = _process_data.specific_body_force;
            // here it is assumed that the vector b is directed 'downwards'
            cache_mat.col(ip).noalias() += rho_w * b;
        }
        cache_mat.col(ip).noalias() = k_rel / mu * (K * cache_mat.col(ip));
    }
    return cache;
}

References MaterialPropertyLib::capillary_pressure, MaterialPropertyLib::concentration, MathLib::createZeroedMatrix(), MaterialPropertyLib::density, NumLib::getIndices(), MaterialPropertyLib::liquid_saturation, MaterialPropertyLib::permeability, MaterialPropertyLib::phase_pressure, MaterialPropertyLib::relative_permeability, MaterialPropertyLib::saturation, ParameterLib::SpatialPosition::setElementID(), NumLib::shapeFunctionInterpolate(), and MaterialPropertyLib::viscosity.

getIntPtSaturation()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

std::vector< double > const & ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::getIntPtSaturation ( const double  t, std::vector< GlobalVector * > const &  x, std::vector< NumLib::LocalToGlobalIndexMap const * > const &  dof_table, std::vector< double > &  cache  ) const

overridevirtual

Implements ProcessLib::RichardsComponentTransport::RichardsComponentTransportLocalAssemblerInterface.

Definition at line 329 of file RichardsComponentTransportFEM-impl.h.

{
    ParameterLib::SpatialPosition pos;
    pos.setElementID(_element_id);
 
    MaterialPropertyLib::VariableArray vars;
 
    auto const& medium = *_process_data.media_map->getMedium(_element_id);
 
    unsigned const n_integration_points =
        _integration_method.getNumberOfPoints();
 
    constexpr int process_id = 0;  // monolithic scheme
    auto const indices =
        NumLib::getIndices(_element_id, *dof_table[process_id]);
    assert(!indices.empty());
    auto const local_x = x[process_id]->get(indices);
 
    cache.clear();
    auto cache_vec = MathLib::createZeroedVector<LocalVectorType>(
        cache, n_integration_points);
 
    for (unsigned ip = 0; ip < n_integration_points; ++ip)
    {
        auto const& ip_data = _ip_data[ip];
        auto const& N = ip_data.N;
 
        double C_int_pt = 0.0;
        double p_int_pt = 0.0;
        NumLib::shapeFunctionInterpolate(local_x, N, C_int_pt, p_int_pt);
 
        // saturation
        vars[static_cast<int>(
            MaterialPropertyLib::Variable::capillary_pressure)] = -p_int_pt;
        auto const dt = std::numeric_limits<double>::quiet_NaN();
        auto const Sw = medium[MaterialPropertyLib::PropertyType::saturation]
                            .template value<double>(vars, pos, t, dt);
        cache_vec[ip] = Sw;
    }
 
    return cache;
}

References MaterialPropertyLib::capillary_pressure, NumLib::getIndices(), MaterialPropertyLib::saturation, ParameterLib::SpatialPosition::setElementID(), and NumLib::shapeFunctionInterpolate().

getShapeMatrix()

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

Eigen::Map< const Eigen::RowVectorXd > ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::getShapeMatrix ( const unsigned  integration_point ) const

overridevirtual

Provides the shape matrix at the given integration point.

Implements NumLib::ExtrapolatableElement.

Definition at line 318 of file RichardsComponentTransportFEM-impl.h.

{
    auto const& N = _ip_data[integration_point].N;
 
    // assumes N is stored contiguously in memory
    return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
}

Member Data Documentation

_element_id

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

unsigned const ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_element_id

private

Definition at line 130 of file RichardsComponentTransportFEM.h.

_integration_method

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

IntegrationMethod const ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_integration_method

private

Definition at line 133 of file RichardsComponentTransportFEM.h.

Referenced by ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::LocalAssemblerData().

_ip_data

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

std::vector< IntegrationPointData<NodalRowVectorType, GlobalDimNodalMatrixType, NodalMatrixType>, Eigen::aligned_allocator<IntegrationPointData< NodalRowVectorType, GlobalDimNodalMatrixType, NodalMatrixType> > > ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_ip_data

private

Definition at line 141 of file RichardsComponentTransportFEM.h.

Referenced by ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::LocalAssemblerData().

_process_data

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

RichardsComponentTransportProcessData const& ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_process_data

private

Definition at line 131 of file RichardsComponentTransportFEM.h.

_transport_process_variable

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

ProcessVariable const& ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::_transport_process_variable

private

Definition at line 134 of file RichardsComponentTransportFEM.h.

concentration_index

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

const int ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::concentration_index = 0

staticprivate

Definition at line 143 of file RichardsComponentTransportFEM.h.

concentration_size

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

const int ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::concentration_size = ShapeFunction::NPOINTS

staticprivate

Definition at line 144 of file RichardsComponentTransportFEM.h.

pressure_index

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

const int ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::pressure_index = ShapeFunction::NPOINTS

staticprivate

Definition at line 145 of file RichardsComponentTransportFEM.h.

pressure_size

template<typename ShapeFunction , typename IntegrationMethod , int GlobalDim>

const int ProcessLib::RichardsComponentTransport::LocalAssemblerData< ShapeFunction, IntegrationMethod, GlobalDim >::pressure_size = ShapeFunction::NPOINTS

staticprivate

Definition at line 146 of file RichardsComponentTransportFEM.h.

---

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

• ProcessLib/RichardsComponentTransport/RichardsComponentTransportFEM.h
• ProcessLib/RichardsComponentTransport/RichardsComponentTransportFEM-impl.h