diff --git a/ProcessLib/ThermoRichardsMechanics/LocalAssemblerInterface.h b/ProcessLib/ThermoRichardsMechanics/LocalAssemblerInterface.h index 7e87c3ab59923d33fc66e3f187a273c6e42a8cd3..965c03a15ea422fb0c30a040228e4ec8607c78da 100644 --- a/ProcessLib/ThermoRichardsMechanics/LocalAssemblerInterface.h +++ b/ProcessLib/ThermoRichardsMechanics/LocalAssemblerInterface.h @@ -18,7 +18,6 @@ #include "ProcessLib/ThermoRichardsMechanics/ConstitutiveSetting.h" #include "ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcessData.h" #include "ProcessLib/Utils/SetOrGetIntegrationPointData.h" -#include "ProcessLib/Utils/TransposeInPlace.h" namespace ProcessLib::ThermoRichardsMechanics { @@ -119,212 +118,6 @@ struct LocalAssemblerInterface : public ProcessLib::LocalAssemblerInterface, return 0; } -private: - std::vector<double> getSigma() const - { - constexpr int kelvin_vector_size = - MathLib::KelvinVector::kelvin_vector_dimensions(DisplacementDim); - - return transposeInPlace<kelvin_vector_size>( - [this](std::vector<double>& values) - { return getIntPtSigma(0, {}, {}, values); }); - } - - std::vector<double> const& getIntPtSigma( - const double /*t*/, - std::vector<GlobalVector*> const& /*x*/, - std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/, - std::vector<double>& cache) const - { - return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>( - current_states_, [](auto const& cs) -> auto const& { - return cs.s_mech_data.sigma_eff; - }, - cache); - } - - std::vector<double> getSwellingStress() const - { - constexpr int kelvin_vector_size = - MathLib::KelvinVector::kelvin_vector_dimensions(DisplacementDim); - - return transposeInPlace<kelvin_vector_size>( - [this](std::vector<double>& values) - { return getIntPtSwellingStress(0, {}, {}, values); }); - } - - std::vector<double> const& getIntPtSwellingStress( - const double /*t*/, - std::vector<GlobalVector*> const& /*x*/, - std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/, - std::vector<double>& cache) const - { - constexpr int kelvin_vector_size = - MathLib::KelvinVector::kelvin_vector_dimensions(DisplacementDim); - auto const n_integration_points = current_states_.size(); - - cache.clear(); - auto cache_mat = MathLib::createZeroedMatrix<Eigen::Matrix< - double, kelvin_vector_size, Eigen::Dynamic, Eigen::RowMajor>>( - cache, kelvin_vector_size, n_integration_points); - - for (unsigned ip = 0; ip < n_integration_points; ++ip) - { - auto const& sigma_sw = current_states_[ip].swelling_data.sigma_sw; - cache_mat.col(ip) = - MathLib::KelvinVector::kelvinVectorToSymmetricTensor(sigma_sw); - } - - return cache; - } - - std::vector<double> getEpsilon() const - { - constexpr int kelvin_vector_size = - MathLib::KelvinVector::kelvin_vector_dimensions(DisplacementDim); - - return transposeInPlace<kelvin_vector_size>( - [this](std::vector<double>& values) - { return getIntPtEpsilon(0, {}, {}, values); }); - } - - std::vector<double> const& getIntPtEpsilon( - const double /*t*/, - std::vector<GlobalVector*> const& /*x*/, - std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/, - std::vector<double>& cache) const - { - return ProcessLib::getIntegrationPointKelvinVectorData<DisplacementDim>( - current_states_, - [](auto const& cs) -> auto const& { return cs.eps_data.eps; }, - cache); - } - - 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 - { - unsigned const n_integration_points = - integration_method_.getNumberOfPoints(); - - cache.clear(); - auto cache_matrix = MathLib::createZeroedMatrix<Eigen::Matrix< - double, DisplacementDim, Eigen::Dynamic, Eigen::RowMajor>>( - cache, DisplacementDim, n_integration_points); - - for (unsigned ip = 0; ip < n_integration_points; ip++) - { - cache_matrix.col(ip).noalias() = - output_data_[ip].darcy_data.v_darcy; - } - - return cache; - } - - std::vector<double> getSaturation() const - { - std::vector<double> result; - getIntPtSaturation(0, {}, {}, result); - return result; - } - - 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 - { - return ProcessLib::getIntegrationPointScalarData( - current_states_, - [](auto const& state) -> auto const& { return state.S_L_data.S_L; }, - cache); - } - - std::vector<double> getPorosity() const - { - std::vector<double> result; - getIntPtPorosity(0, {}, {}, result); - return result; - } - - std::vector<double> const& getIntPtPorosity( - const double /*t*/, - std::vector<GlobalVector*> const& /*x*/, - std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/, - std::vector<double>& cache) const - { - return ProcessLib::getIntegrationPointScalarData( - current_states_, [](auto const& state) -> auto const& { - return state.poro_data.phi; - }, - cache); - } - - std::vector<double> getTransportPorosity() const - { - std::vector<double> result; - getIntPtTransportPorosity(0, {}, {}, result); - return result; - } - - std::vector<double> const& getIntPtTransportPorosity( - const double /*t*/, - std::vector<GlobalVector*> const& /*x*/, - std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/, - std::vector<double>& cache) const - { - return ProcessLib::getIntegrationPointScalarData( - current_states_, - [](auto const& state) -> auto const& { - return state.transport_poro_data.phi; - }, - cache); - } - - std::vector<double> const& getIntPtDryDensitySolid( - const double /*t*/, - std::vector<GlobalVector*> const& /*x*/, - std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/, - std::vector<double>& cache) const - { - return ProcessLib::getIntegrationPointScalarData( - output_data_, - [](auto const& out) -> auto const& { - return out.rho_S_data.dry_density_solid; - }, - cache); - } - - std::vector<double> const& getIntPtLiquidDensity( - const double /*t*/, - std::vector<GlobalVector*> const& /*x*/, - std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/, - std::vector<double>& cache) const - { - return ProcessLib::getIntegrationPointScalarData( - output_data_, - [](auto const& out) -> auto const& { - return out.rho_L_data.rho_LR; - }, - cache); - } - - std::vector<double> const& getIntPtViscosity( - const double /*t*/, - std::vector<GlobalVector*> const& /*x*/, - std::vector<NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/, - std::vector<double>& cache) const - { - return ProcessLib::getIntegrationPointScalarData( - output_data_, [](auto const& out) -> auto const& { - return out.mu_L_data.viscosity; - }, - cache); - } - -public: // TODO move to NumLib::ExtrapolatableElement unsigned getNumberOfIntegrationPoints() const { @@ -354,65 +147,12 @@ public: prev_states_ = current_states_; } - struct IPDataAccessorForExtrapolation - { - std::string name; - unsigned num_comp; - std::function<std::vector<double> const&( - LocalAssemblerInterface<DisplacementDim> const&, - const double /*t*/, - std::vector<GlobalVector*> const& /*x*/, - std::vector< - NumLib::LocalToGlobalIndexMap const*> const& /*dof_table*/, - std::vector<double>& /*cache*/)> - ip_data_accessor; - }; - - static std::vector<IPDataAccessorForExtrapolation> - getIPDataAccessorsForExtrapolation() - { - using Self = LocalAssemblerInterface<DisplacementDim>; - constexpr auto kv_size = - MathLib::KelvinVector::kelvin_vector_dimensions(DisplacementDim); - - return {{"sigma", kv_size, &Self::getIntPtSigma}, - {"swelling_stress", kv_size, &Self::getIntPtSwellingStress}, - {"epsilon", kv_size, &Self::getIntPtEpsilon}, - {"velocity", DisplacementDim, &Self::getIntPtDarcyVelocity}, - {"saturation", 1, &Self::getIntPtSaturation}, - {"porosity", 1, &Self::getIntPtPorosity}, - {"transport_porosity", 1, &Self::getIntPtTransportPorosity}, - {"dry_density_solid", 1, &Self::getIntPtDryDensitySolid}, - {"liquid_density", 1, &Self::getIntPtLiquidDensity}, - {"viscosity", 1, &Self::getIntPtViscosity}}; - } - - struct IPDataAccessorForIPWriter - { - std::string name; - unsigned num_comp; - std::vector<double> (LocalAssemblerInterface<DisplacementDim>::* - ip_data_accessor)() const; - }; - - static std::vector<IPDataAccessorForIPWriter> - getIPDataAccessorsForIPWriter() + static auto getReflectionDataForOutput() { using Self = LocalAssemblerInterface<DisplacementDim>; - constexpr auto kv_size = - MathLib::KelvinVector::kelvin_vector_dimensions(DisplacementDim); - // TODO (naumov) remove ip suffix. Probably needs modification of the - // mesh properties, s.t. there is no "overlapping" with cell/point data. - // See getOrCreateMeshProperty. - return { - {"sigma_ip", kv_size, &Self::getSigma}, - {"saturation_ip", 1, &Self::getSaturation}, - {"porosity_ip", 1, &Self::getPorosity}, - {"transport_porosity_ip", 1, &Self::getTransportPorosity}, - {"swelling_stress_ip", kv_size, &Self::getSwellingStress}, - {"epsilon_ip", kv_size, &Self::getEpsilon}, - }; + return ProcessLib::Reflection::reflectWithoutName( + &Self::current_states_, &Self::output_data_); } protected: diff --git a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.cpp b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.cpp index fc1d1b0dbd1c8567249042589f7f8cef084cdaca..439a664a56c5cc73a74e6cb3fdc9cab4ef9bc986 100644 --- a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.cpp +++ b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.cpp @@ -17,6 +17,8 @@ #include "NumLib/DOF/DOFTableUtil.h" #include "ProcessLib/Deformation/SolidMaterialInternalToSecondaryVariables.h" #include "ProcessLib/Process.h" +#include "ProcessLib/Reflection/ReflectionForExtrapolation.h" +#include "ProcessLib/Reflection/ReflectionForIPWriters.h" #include "ProcessLib/Utils/CreateLocalAssemblersTaylorHood.h" #include "ProcessLib/Utils/SetIPDataInitialConditions.h" #include "ThermoRichardsMechanicsFEM.h" @@ -51,14 +53,10 @@ ThermoRichardsMechanicsProcess<DisplacementDim>::ThermoRichardsMechanicsProcess( heat_flux_ = MeshLib::getOrCreateMeshProperty<double>( mesh, "HeatFlowRate", MeshLib::MeshItemType::Node, 1); - for (auto& [name, num_comp, ip_data_accessor] : - LocalAssemblerIF::getIPDataAccessorsForIPWriter()) - { - _integration_point_writer.emplace_back( - std::make_unique<MeshLib::IntegrationPointWriter>( - name, num_comp, integration_order, local_assemblers_, - ip_data_accessor)); - } + ProcessLib::Reflection::addReflectedIntegrationPointWriters< + DisplacementDim>(LocalAssemblerIF::getReflectionDataForOutput(), + _integration_point_writer, integration_order, + local_assemblers_); } template <int DisplacementDim> @@ -143,11 +141,9 @@ void ThermoRichardsMechanicsProcess<DisplacementDim>::initializeConcreteProcess( std::move(get_ip_values_function))); }; - for (auto& [name, num_comp, fct] : LocalAssemblerInterface< - DisplacementDim>::getIPDataAccessorsForExtrapolation()) - { - add_secondary_variable(name, num_comp, fct); - }; + ProcessLib::Reflection::addReflectedSecondaryVariables<DisplacementDim>( + LocalAssemblerIF::getReflectionDataForOutput(), _secondary_variables, + getExtrapolator(), local_assemblers_); // // enable output of internal variables defined by material models