diff --git a/ProcessLib/ThermoMechanics/LocalAssemblerInterface.h b/ProcessLib/ThermoMechanics/LocalAssemblerInterface.h
index 36d82fe137186abd2d2a99cbe72012a7462c5968..f3af5f732a494d378b75576ba909c19598386145 100644
--- a/ProcessLib/ThermoMechanics/LocalAssemblerInterface.h
+++ b/ProcessLib/ThermoMechanics/LocalAssemblerInterface.h
@@ -22,77 +22,24 @@ struct ThermoMechanicsLocalAssemblerInterface
: public ProcessLib::LocalAssemblerInterface,
public NumLib::ExtrapolatableElement
{
- virtual std::vector<double> const& getIntPtSigmaXX(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const = 0;
-
- virtual std::vector<double> const& getIntPtSigmaYY(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const = 0;
-
- virtual std::vector<double> const& getIntPtSigmaZZ(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const = 0;
-
- virtual std::vector<double> const& getIntPtSigmaXY(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const = 0;
+ virtual std::size_t setIPDataInitialConditions(
+ std::string const& name, double const* values,
+ int const integration_order) = 0;
- virtual std::vector<double> const& getIntPtSigmaXZ(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const = 0;
+ virtual std::vector<double> getSigma() const = 0;
- virtual std::vector<double> const& getIntPtSigmaYZ(
+ virtual std::vector<double> const& getIntPtSigma(
const double /*t*/,
GlobalVector const& /*current_solution*/,
NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
std::vector<double>& cache) const = 0;
- virtual std::vector<double> const& getIntPtEpsilonXX(
+ virtual std::vector<double> const& getIntPtEpsilon(
const double /*t*/,
GlobalVector const& /*current_solution*/,
NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
std::vector<double>& cache) const = 0;
- virtual std::vector<double> const& getIntPtEpsilonYY(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const = 0;
-
- virtual std::vector<double> const& getIntPtEpsilonZZ(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const = 0;
-
- virtual std::vector<double> const& getIntPtEpsilonXY(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const = 0;
-
- virtual std::vector<double> const& getIntPtEpsilonXZ(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const = 0;
-
- virtual std::vector<double> const& getIntPtEpsilonYZ(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const = 0;
};
} // namespace ThermoMechanics
diff --git a/ProcessLib/ThermoMechanics/Tests.cmake b/ProcessLib/ThermoMechanics/Tests.cmake
index 8adab6a0ee881d53a901d2b7e31dab732f94af27..f09859b087d653b1d3e75c3e055be1397962bea0 100644
--- a/ProcessLib/ThermoMechanics/Tests.cmake
+++ b/ProcessLib/ThermoMechanics/Tests.cmake
@@ -1,5 +1,5 @@
AddTest(
- NAME 3D_ThermoElastic_Stress_Analysis
+ NAME ThermoMechanics_3D_ThermoElastic_Stress_Analysis
PATH ThermoMechanics
EXECUTABLE ogs
EXECUTABLE_ARGS cube_1e3.prj
@@ -7,15 +7,15 @@ AddTest(
TESTER vtkdiff
REQUIREMENTS NOT (OGS_USE_LIS OR OGS_USE_MPI)
DIFF_DATA
- stress_analytical.vtu cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu sigma_xx sigma_xx 1e-10 1e-11
- stress_analytical.vtu cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu sigma_yy sigma_yy 1e-10 1e-11
- stress_analytical.vtu cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu sigma_zz sigma_zz 1e-10 1e-11
+ stress_analytical.vtu cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu sigma sigma 1e-5 1e-12
expected_cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu displacement displacement 1e-10 1e-12
expected_cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu temperature temperature 1e-10 1e-12
+ expected_cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu sigma sigma 1e-6 1e-12
+ expected_cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu epsilon epsilon 1e-16 0
)
AddTest(
- NAME 2D_ThermoElastic_IGLU_Plane_Strain
+ NAME ThermoMechanics_2D_ThermoElastic_IGLU_Plane_Strain
PATH ThermoMechanics
EXECUTABLE ogs
EXECUTABLE_ARGS iglu_quarter_plane_strain.prj
@@ -25,13 +25,12 @@ AddTest(
DIFF_DATA
expected_tm_q_pcs_0_ts_20_t_20000.000000.vtu tm_q_pcs_0_ts_20_t_20000.000000.vtu displacement displacement 1e-9 1e-15
expected_tm_q_pcs_0_ts_20_t_20000.000000.vtu tm_q_pcs_0_ts_20_t_20000.000000.vtu temperature temperature 2e-6 1e-15
- expected_tm_q_pcs_0_ts_20_t_20000.000000.vtu tm_q_pcs_0_ts_20_t_20000.000000.vtu sigma_xx sigma_xx 5e-6 1e-15
- expected_tm_q_pcs_0_ts_20_t_20000.000000.vtu tm_q_pcs_0_ts_20_t_20000.000000.vtu sigma_yy sigma_yy 5e-6 1e-15
- expected_tm_q_pcs_0_ts_20_t_20000.000000.vtu tm_q_pcs_0_ts_20_t_20000.000000.vtu sigma_zz sigma_zz 5e-6 1e-15
+ expected_tm_q_pcs_0_ts_20_t_20000.000000.vtu tm_q_pcs_0_ts_20_t_20000.000000.vtu sigma sigma 5e-6 1e-15
+ expected_tm_q_pcs_0_ts_20_t_20000.000000.vtu tm_q_pcs_0_ts_20_t_20000.000000.vtu epsilon epsilon 5e-6 1e-15
)
AddTest(
- NAME 2D_ThermoElastic_IGLU_Axisymmetric_Plane_Strain
+ NAME ThermoMechanics_2D_ThermoElastic_IGLU_Axisymmetric_Plane_Strain
PATH ThermoMechanics
EXECUTABLE ogs
EXECUTABLE_ARGS iglu_axisymmetric_plane_strain.prj
@@ -41,13 +40,12 @@ AddTest(
DIFF_DATA
expected_tm_a_pcs_0_ts_20_t_20000.000000.vtu tm_a_pcs_0_ts_20_t_20000.000000.vtu displacement displacement 1e-9 1e-15
expected_tm_a_pcs_0_ts_20_t_20000.000000.vtu tm_a_pcs_0_ts_20_t_20000.000000.vtu temperature temperature 1e-10 1e-8
- expected_tm_a_pcs_0_ts_20_t_20000.000000.vtu tm_a_pcs_0_ts_20_t_20000.000000.vtu sigma_xx sigma_xx 1e-10 1e-6
- expected_tm_a_pcs_0_ts_20_t_20000.000000.vtu tm_a_pcs_0_ts_20_t_20000.000000.vtu sigma_yy sigma_yy 1e-10 1e-6
- expected_tm_a_pcs_0_ts_20_t_20000.000000.vtu tm_a_pcs_0_ts_20_t_20000.000000.vtu sigma_zz sigma_zz 1e-10 1e-6
+ expected_tm_a_pcs_0_ts_20_t_20000.000000.vtu tm_a_pcs_0_ts_20_t_20000.000000.vtu sigma sigma 1e-6 0
+ expected_tm_a_pcs_0_ts_20_t_20000.000000.vtu tm_a_pcs_0_ts_20_t_20000.000000.vtu epsilon epsilon 1e-10 0
)
AddTest(
- NAME LARGE_2D_ThermoElastic_IGLU_Plane_Strain_Quadratic_Mesh
+ NAME ThermoMechanics_LARGE_2D_ThermoElastic_IGLU_Plane_Strain_Quadratic_Mesh
PATH ThermoMechanics
EXECUTABLE ogs
EXECUTABLE_ARGS iglu_quarter_plane_strain_quad.prj
@@ -57,13 +55,12 @@ AddTest(
DIFF_DATA
expected_tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu displacement displacement 5e-10 1e-15
expected_tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu temperature temperature 2e-6 1e-15
- expected_tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu sigma_xx sigma_xx 5e-6 0
- expected_tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu sigma_yy sigma_yy 6e-6 0
- expected_tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu sigma_zz sigma_zz 5e-6 0
+ expected_tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu sigma sigma 5e-6 0
+ expected_tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu epsilon epsilon 6e-6 0
)
AddTest(
- NAME 2D_ThermoElastic_IGLU_Axisymmetric_Plane_Strain_Quadratic_Mesh
+ NAME ThermoMechanics_2D_ThermoElastic_IGLU_Axisymmetric_Plane_Strain_Quadratic_Mesh
PATH ThermoMechanics
EXECUTABLE ogs
EXECUTABLE_ARGS iglu_axisymmetric_plane_strain_quad.prj
@@ -73,9 +70,8 @@ AddTest(
DIFF_DATA
expected_tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu displacement displacement 5e-10 1e-15
expected_tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu temperature temperature 5e-10 1e-7
- expected_tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu sigma_xx sigma_xx 5e-10 1e-6
- expected_tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu sigma_yy sigma_yy 5e-10 1e-8
- expected_tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu sigma_zz sigma_zz 5e-10 1e-7
+ expected_tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu sigma sigma 1e-6 0
+ expected_tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu epsilon epsilon 5e-10 1e-8
)
AddTest(
@@ -87,8 +83,9 @@ AddTest(
TESTER vtkdiff
REQUIREMENTS NOT (OGS_USE_LIS OR OGS_USE_MPI)
DIFF_DATA
- expected_SimpleAxisymmetricCreep_pcs_0_ts_370_t_360.000000.vtu SimpleAxisymmetricCreep_pcs_0_ts_370_t_360.000000.vtu displacement displacement 1e-14 1e-10
- expected_SimpleAxisymmetricCreep_pcs_0_ts_370_t_360.000000.vtu SimpleAxisymmetricCreep_pcs_0_ts_370_t_360.000000.vtu sigma_yy sigma_yy 1e-16 1e-10
+ expected_SimpleAxisymmetricCreep_pcs_0_ts_370_t_360.000000.vtu SimpleAxisymmetricCreep_pcs_0_ts_370_t_360.000000.vtu displacement displacement 1e-14 1e-10
+ expected_SimpleAxisymmetricCreep_pcs_0_ts_370_t_360.000000.vtu SimpleAxisymmetricCreep_pcs_0_ts_370_t_360.000000.vtu sigma sigma 1e-7 0
+ expected_SimpleAxisymmetricCreep_pcs_0_ts_370_t_360.000000.vtu SimpleAxisymmetricCreep_pcs_0_ts_370_t_360.000000.vtu epsilon epsilon 1e-12 0
)
AddTest(
@@ -100,7 +97,7 @@ AddTest(
TESTER vtkdiff
REQUIREMENTS NOT (OGS_USE_LIS OR OGS_USE_MPI)
DIFF_DATA
- SimpleAxisymmetricCreepWithAnalyticSolution.vtu SimpleAxisymmetricCreepWithAnalyticalSolution_pcs_0_ts_1000_t_100.000000.vtu analytic_strain_yy epsilon_yy 1e-10 1e-4
+ SimpleAxisymmetricCreepWithAnalyticSolution.vtu SimpleAxisymmetricCreepWithAnalyticalSolution_pcs_0_ts_1000_t_100.000000.vtu analytic_strain epsilon 1e-7 0
)
AddTest(
@@ -112,7 +109,7 @@ AddTest(
TESTER vtkdiff
REQUIREMENTS NOT (OGS_USE_LIS OR OGS_USE_MPI)
DIFF_DATA
- ExpectedCreepAfterExcavation_pcs_0_ts_61_t_4320000.000000.vtu CreepAfterExcavation_pcs_0_ts_61_t_4320000.000000.vtu sigma_xx sigma_xx 1e-16 1e-11
- ExpectedCreepAfterExcavation_pcs_0_ts_61_t_4320000.000000.vtu CreepAfterExcavation_pcs_0_ts_61_t_4320000.000000.vtu sigma_yy sigma_yy 1e-16 1e-11
+ ExpectedCreepAfterExcavation_pcs_0_ts_61_t_4320000.000000.vtu CreepAfterExcavation_pcs_0_ts_61_t_4320000.000000.vtu sigma sigma 5e-6 0
+ ExpectedCreepAfterExcavation_pcs_0_ts_61_t_4320000.000000.vtu CreepAfterExcavation_pcs_0_ts_61_t_4320000.000000.vtu epsilon epsilon 1e-15 0
ExpectedCreepAfterExcavation_pcs_0_ts_61_t_4320000.000000.vtu CreepAfterExcavation_pcs_0_ts_61_t_4320000.000000.vtu displacement displacement 1e-16 1e-9
)
diff --git a/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h b/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h
index bf4abac34f6b5565de65462e3795401c73f21bc4..0cae1175f97cb6abdd68cc1bc9128f9941e0b963 100644
--- a/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h
+++ b/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h
@@ -154,6 +154,29 @@ public:
}
}
+ /// Returns number of read integration points.
+ std::size_t setIPDataInitialConditions(std::string const& name,
+ double const* values,
+ int const integration_order) override
+ {
+ if (integration_order !=
+ static_cast<int>(_integration_method.getIntegrationOrder()))
+ {
+ OGS_FATAL(
+ "Setting integration point initial conditions; The integration "
+ "order of the local assembler for element %d is different from "
+ "the integration order in the initial condition.",
+ _element.getID());
+ }
+
+ if (name == "sigma_ip")
+ {
+ return setSigma(values);
+ }
+
+ return 0;
+ }
+
void assemble(double const /*t*/, std::vector<double> const& /*local_x*/,
std::vector<double>& /*local_M_data*/,
std::vector<double>& /*local_K_data*/,
@@ -365,148 +388,106 @@ public:
return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
}
- std::vector<double> const& getIntPtSigmaXX(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const override
- {
- return getIntPtSigma(cache, 0);
- }
-
- std::vector<double> const& getIntPtSigmaYY(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const override
- {
- return getIntPtSigma(cache, 1);
- }
-
- std::vector<double> const& getIntPtSigmaZZ(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const override
- {
- return getIntPtSigma(cache, 2);
- }
-
- std::vector<double> const& getIntPtSigmaXY(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const override
+private:
+ std::size_t setSigma(double const* values)
{
- return getIntPtSigma(cache, 3);
- }
+ auto const kelvin_vector_size =
+ MathLib::KelvinVector::KelvinVectorDimensions<
+ DisplacementDim>::value;
+ unsigned const n_integration_points =
+ _integration_method.getNumberOfPoints();
- std::vector<double> const& getIntPtSigmaYZ(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const override
- {
- assert(DisplacementDim == 3);
- return getIntPtSigma(cache, 4);
- }
+ std::vector<double> ip_sigma_values;
+ auto sigma_values =
+ Eigen::Map<Eigen::Matrix<double, kelvin_vector_size, Eigen::Dynamic,
+ Eigen::ColMajor> const>(
+ values, kelvin_vector_size, n_integration_points);
- std::vector<double> const& getIntPtSigmaXZ(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const override
- {
- assert(DisplacementDim == 3);
- return getIntPtSigma(cache, 5);
- }
+ for (unsigned ip = 0; ip < n_integration_points; ++ip)
+ {
+ _ip_data[ip].sigma =
+ MathLib::KelvinVector::symmetricTensorToKelvinVector(
+ sigma_values.col(ip));
+ }
- std::vector<double> const& getIntPtEpsilonXX(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const override
- {
- return getIntPtEpsilon(cache, 0);
+ return n_integration_points;
}
- std::vector<double> const& getIntPtEpsilonYY(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const override
+ // TODO (naumov) This method is same as getIntPtSigma but for arguments and
+ // the ordering of the cache_mat.
+ // There should be only one.
+ std::vector<double> getSigma() const override
{
- return getIntPtEpsilon(cache, 1);
- }
+ auto const kelvin_vector_size =
+ MathLib::KelvinVector::KelvinVectorDimensions<
+ DisplacementDim>::value;
+ unsigned const n_integration_points =
+ _integration_method.getNumberOfPoints();
- std::vector<double> const& getIntPtEpsilonZZ(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const override
- {
- return getIntPtEpsilon(cache, 2);
- }
+ std::vector<double> ip_sigma_values;
+ auto cache_mat = MathLib::createZeroedMatrix<Eigen::Matrix<
+ double, Eigen::Dynamic, kelvin_vector_size, Eigen::RowMajor>>(
+ ip_sigma_values, n_integration_points, kelvin_vector_size);
- std::vector<double> const& getIntPtEpsilonXY(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const override
- {
- return getIntPtEpsilon(cache, 3);
- }
+ for (unsigned ip = 0; ip < n_integration_points; ++ip)
+ {
+ auto const& sigma = _ip_data[ip].sigma;
+ cache_mat.row(ip) =
+ MathLib::KelvinVector::kelvinVectorToSymmetricTensor(sigma);
+ }
- std::vector<double> const& getIntPtEpsilonYZ(
- const double /*t*/,
- GlobalVector const& /*current_solution*/,
- NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
- std::vector<double>& cache) const override
- {
- assert(DisplacementDim == 3);
- return getIntPtEpsilon(cache, 4);
+ return ip_sigma_values;
}
- std::vector<double> const& getIntPtEpsilonXZ(
+ std::vector<double> const& getIntPtSigma(
const double /*t*/,
GlobalVector const& /*current_solution*/,
NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
std::vector<double>& cache) const override
{
- assert(DisplacementDim == 3);
- return getIntPtEpsilon(cache, 5);
- }
+ static const int kelvin_vector_size =
+ MathLib::KelvinVector::KelvinVectorDimensions<
+ DisplacementDim>::value;
+ unsigned const n_integration_points =
+ _integration_method.getNumberOfPoints();
-private:
- std::vector<double> const& getIntPtSigma(std::vector<double>& cache,
- std::size_t const component) const
- {
cache.clear();
- cache.reserve(_ip_data.size());
+ auto cache_mat = MathLib::createZeroedMatrix<Eigen::Matrix<
+ double, kelvin_vector_size, Eigen::Dynamic, Eigen::RowMajor>>(
+ cache, kelvin_vector_size, n_integration_points);
- for (auto const& ip_data : _ip_data)
+ for (unsigned ip = 0; ip < n_integration_points; ++ip)
{
- if (component < 3) // xx, yy, zz components
- cache.push_back(ip_data.sigma[component]);
- else // mixed xy, yz, xz components
- cache.push_back(ip_data.sigma[component] / std::sqrt(2));
+ auto const& sigma = _ip_data[ip].sigma;
+ cache_mat.col(ip) =
+ MathLib::KelvinVector::kelvinVectorToSymmetricTensor(sigma);
}
return cache;
}
- std::vector<double> const& getIntPtEpsilon(
- std::vector<double>& cache, std::size_t const component) const
+ virtual std::vector<double> const& getIntPtEpsilon(
+ const double /*t*/,
+ GlobalVector const& /*current_solution*/,
+ NumLib::LocalToGlobalIndexMap const& /*dof_table*/,
+ std::vector<double>& cache) const override
{
+ auto const kelvin_vector_size =
+ MathLib::KelvinVector::KelvinVectorDimensions<
+ DisplacementDim>::value;
+ unsigned const n_integration_points =
+ _integration_method.getNumberOfPoints();
+
cache.clear();
- cache.reserve(_ip_data.size());
+ auto cache_mat = MathLib::createZeroedMatrix<Eigen::Matrix<
+ double, kelvin_vector_size, Eigen::Dynamic, Eigen::RowMajor>>(
+ cache, kelvin_vector_size, n_integration_points);
- for (auto const& ip_data : _ip_data)
+ for (unsigned ip = 0; ip < n_integration_points; ++ip)
{
- if (component < 3) // xx, yy, zz components
- cache.push_back(ip_data.eps[component]);
- else // mixed xy, yz, xz components
- cache.push_back(ip_data.eps[component] / std::sqrt(2));
+ auto const& eps = _ip_data[ip].eps;
+ cache_mat.col(ip) =
+ MathLib::KelvinVector::kelvinVectorToSymmetricTensor(eps);
}
return cache;
diff --git a/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp b/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp
index d0fa7d2226509d5c1a34b696d37586096b55a700..cf73f8ca353c6154e81e45b3abe21e46b2395526 100644
--- a/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp
+++ b/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp
@@ -38,6 +38,24 @@ ThermoMechanicsProcess<DisplacementDim>::ThermoMechanicsProcess(
use_monolithic_scheme),
_process_data(std::move(process_data))
{
+ _integration_point_writer.emplace_back(
+ std::make_unique<SigmaIntegrationPointWriter>(
+ static_cast<int>(mesh.getDimension() == 2 ? 4 : 6) /*n components*/,
+ 2 /*integration order*/, [this]() {
+ // Result containing integration point data for each local
+ // assembler.
+ std::vector<std::vector<double>> result;
+ result.resize(_local_assemblers.size());
+
+ for (std::size_t i = 0; i < _local_assemblers.size(); ++i)
+ {
+ auto const& local_asm = *_local_assemblers[i];
+
+ result[i] = local_asm.getSigma();
+ }
+
+ return result;
+ }));
}
template <int DisplacementDim>
@@ -68,79 +86,69 @@ void ThermoMechanicsProcess<DisplacementDim>::initializeConcreteProcess(
NumLib::ComponentOrder::BY_LOCATION));
_secondary_variables.addSecondaryVariable(
- "sigma_xx",
+ "sigma",
makeExtrapolator(
- 1, getExtrapolator(), _local_assemblers,
- &ThermoMechanicsLocalAssemblerInterface::getIntPtSigmaXX));
+ MathLib::KelvinVector::KelvinVectorType<
+ DisplacementDim>::RowsAtCompileTime,
+ getExtrapolator(), _local_assemblers,
+ &ThermoMechanicsLocalAssemblerInterface::getIntPtSigma));
_secondary_variables.addSecondaryVariable(
- "sigma_yy",
+ "epsilon",
makeExtrapolator(
- 1, getExtrapolator(), _local_assemblers,
- &ThermoMechanicsLocalAssemblerInterface::getIntPtSigmaYY));
+ MathLib::KelvinVector::KelvinVectorType<
+ DisplacementDim>::RowsAtCompileTime,
+ getExtrapolator(), _local_assemblers,
+ &ThermoMechanicsLocalAssemblerInterface::getIntPtEpsilon));
- _secondary_variables.addSecondaryVariable(
- "sigma_zz",
- makeExtrapolator(
- 1, getExtrapolator(), _local_assemblers,
- &ThermoMechanicsLocalAssemblerInterface::getIntPtSigmaZZ));
-
- _secondary_variables.addSecondaryVariable(
- "sigma_xy",
- makeExtrapolator(
- 1, getExtrapolator(), _local_assemblers,
- &ThermoMechanicsLocalAssemblerInterface::getIntPtSigmaXY));
-
- if (DisplacementDim == 3)
- {
- _secondary_variables.addSecondaryVariable(
- "sigma_xz",
- makeExtrapolator(
- 1, getExtrapolator(), _local_assemblers,
- &ThermoMechanicsLocalAssemblerInterface::getIntPtSigmaXZ));
-
- _secondary_variables.addSecondaryVariable(
- "sigma_yz",
- makeExtrapolator(
- 1, getExtrapolator(), _local_assemblers,
- &ThermoMechanicsLocalAssemblerInterface::getIntPtSigmaYZ));
- }
- _secondary_variables.addSecondaryVariable(
- "epsilon_xx",
- makeExtrapolator(
- 1, getExtrapolator(), _local_assemblers,
- &ThermoMechanicsLocalAssemblerInterface::getIntPtEpsilonXX));
-
- _secondary_variables.addSecondaryVariable(
- "epsilon_yy",
- makeExtrapolator(
- 1, getExtrapolator(), _local_assemblers,
- &ThermoMechanicsLocalAssemblerInterface::getIntPtEpsilonYY));
-
- _secondary_variables.addSecondaryVariable(
- "epsilon_zz",
- makeExtrapolator(
- 1, getExtrapolator(), _local_assemblers,
- &ThermoMechanicsLocalAssemblerInterface::getIntPtEpsilonZZ));
-
- _secondary_variables.addSecondaryVariable(
- "epsilon_xy",
- makeExtrapolator(
- 1, getExtrapolator(), _local_assemblers,
- &ThermoMechanicsLocalAssemblerInterface::getIntPtEpsilonXY));
- if (DisplacementDim == 3)
+ // Set initial conditions for integration point data.
+ for (auto const& ip_writer : _integration_point_writer)
{
- _secondary_variables.addSecondaryVariable(
- "epsilon_yz",
- makeExtrapolator(
- 1, getExtrapolator(), _local_assemblers,
- &ThermoMechanicsLocalAssemblerInterface::getIntPtEpsilonYZ));
-
- _secondary_variables.addSecondaryVariable(
- "epsilon_xz",
- makeExtrapolator(
- 1, getExtrapolator(), _local_assemblers,
- &ThermoMechanicsLocalAssemblerInterface::getIntPtEpsilonXZ));
+ // Find the mesh property with integration point writer's name.
+ auto const& name = ip_writer->name();
+ if (!mesh.getProperties().existsPropertyVector<double>(name))
+ {
+ continue;
+ }
+ auto const& mesh_property =
+ *mesh.getProperties().template getPropertyVector<double>(name);
+
+ // The mesh property must be defined on integration points.
+ if (mesh_property.getMeshItemType() !=
+ MeshLib::MeshItemType::IntegrationPoint)
+ {
+ continue;
+ }
+
+ auto const ip_meta_data = getIntegrationPointMetaData(mesh, name);
+
+ // Check the number of components.
+ if (ip_meta_data.n_components != mesh_property.getNumberOfComponents())
+ {
+ OGS_FATAL(
+ "Different number of components in meta data (%d) than in "
+ "the integration point field data for \"%s\": %d.",
+ ip_meta_data.n_components, name.c_str(),
+ mesh_property.getNumberOfComponents());
+ }
+
+ // Now we have a properly named vtk's field data array and the
+ // corresponding meta data.
+ std::size_t position = 0;
+ for (auto& local_asm : _local_assemblers)
+ {
+ std::size_t const integration_points_read =
+ local_asm->setIPDataInitialConditions(
+ name, &mesh_property[position],
+ ip_meta_data.integration_order);
+ if (integration_points_read == 0)
+ {
+ OGS_FATAL(
+ "No integration points read in the integration point "
+ "initial conditions set function.");
+ }
+ position += integration_points_read * ip_meta_data.n_components;
+ }
}
}
diff --git a/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.h b/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.h
index 25b6f37fe049eb8ca80d94bfec84399a82467358..77a55ba1be47ae5daee210eeae4eac550427d7c9 100644
--- a/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.h
+++ b/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.h
@@ -20,6 +20,41 @@ namespace ThermoMechanics
{
struct ThermoMechanicsLocalAssemblerInterface;
+struct SigmaIntegrationPointWriter final : public IntegrationPointWriter
+{
+ explicit SigmaIntegrationPointWriter(
+ int const n_components,
+ int const integration_order,
+ std::function<std::vector<std::vector<double>>()>
+ callback)
+ : _n_components(n_components),
+ _integration_order(integration_order),
+ _callback(callback)
+ {
+ }
+
+ int numberOfComponents() const override { return _n_components; }
+ int integrationOrder() const override { return _integration_order; }
+
+ std::string name() const override
+ {
+ // 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";
+ }
+
+ std::vector<std::vector<double>> values() const override
+ {
+ return _callback();
+ }
+
+private:
+ int const _n_components;
+ int const _integration_order;
+ std::function<std::vector<std::vector<double>>()> _callback;
+};
+
template <int DisplacementDim>
class ThermoMechanicsProcess final : public Process
{
diff --git a/Tests/Data/ThermoMechanics/CreepBGRa/CreepAfterExcavation/CreepAfterExcavation.prj b/Tests/Data/ThermoMechanics/CreepBGRa/CreepAfterExcavation/CreepAfterExcavation.prj
index f25851f5eaa10d3f2b7258e608f7f16421a48a7b..1e00e2049f10e33f1844e0f77b470c43f53c486f 100644
--- a/Tests/Data/ThermoMechanics/CreepBGRa/CreepAfterExcavation/CreepAfterExcavation.prj
+++ b/Tests/Data/ThermoMechanics/CreepBGRa/CreepAfterExcavation/CreepAfterExcavation.prj
@@ -33,14 +33,8 @@
<temperature>temperature</temperature>
</process_variables>
<secondary_variables>
- <secondary_variable type="static" internal_name="sigma_xx" output_name="sigma_xx"/>
- <secondary_variable type="static" internal_name="sigma_yy" output_name="sigma_yy"/>
- <secondary_variable type="static" internal_name="sigma_zz" output_name="sigma_zz"/>
- <secondary_variable type="static" internal_name="sigma_xy" output_name="sigma_xy"/>
- <secondary_variable type="static" internal_name="epsilon_xx" output_name="epsilon_xx"/>
- <secondary_variable type="static" internal_name="epsilon_yy" output_name="epsilon_yy"/>
- <secondary_variable type="static" internal_name="epsilon_zz" output_name="epsilon_zz"/>
- <secondary_variable type="static" internal_name="epsilon_xy" output_name="epsilon_xy"/>
+ <secondary_variable type="static" internal_name="sigma" output_name="sigma"/>
+ <secondary_variable type="static" internal_name="epsilon" output_name="epsilon"/>
</secondary_variables>
<specific_body_force>0 -9.81</specific_body_force>
</process>
@@ -61,14 +55,8 @@
<variables>
<variable>displacement</variable>
<variable>temperature</variable>
- <variable>sigma_xx</variable>
- <variable>sigma_yy</variable>
- <variable>sigma_zz</variable>
- <variable>sigma_xy</variable>
- <variable>epsilon_xx</variable>
- <variable>epsilon_yy</variable>
- <variable>epsilon_zz</variable>
- <variable>epsilon_xy</variable>
+ <variable>sigma</variable>
+ <variable>epsilon</variable>
</variables>
</output>
<time_stepping>
diff --git a/Tests/Data/ThermoMechanics/CreepBGRa/CreepAfterExcavation/ExpectedCreepAfterExcavation_pcs_0_ts_61_t_4320000.000000.vtu b/Tests/Data/ThermoMechanics/CreepBGRa/CreepAfterExcavation/ExpectedCreepAfterExcavation_pcs_0_ts_61_t_4320000.000000.vtu
index 4e6b50bc16400b58edc36580b1bb45088dcede11..4838dd2b62316cdf07842adea13127e7e220d0c2 100644
--- a/Tests/Data/ThermoMechanics/CreepBGRa/CreepAfterExcavation/ExpectedCreepAfterExcavation_pcs_0_ts_61_t_4320000.000000.vtu
+++ b/Tests/Data/ThermoMechanics/CreepBGRa/CreepAfterExcavation/ExpectedCreepAfterExcavation_pcs_0_ts_61_t_4320000.000000.vtu
@@ -1,3 +1,3 @@
version https://git-lfs.github.com/spec/v1
-oid sha256:f7ff2e11894119a96e8b7110d7a51b7d117378cc9ddbac4da8ad93cafe9ac8ee
-size 185317
+oid sha256:a32fbd57e28fc24616d41156c573b70ed13e231829336f70f7bd2b894f5ad618
+size 186595
diff --git a/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/SimpleAxisymmetricCreep.prj b/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/SimpleAxisymmetricCreep.prj
index e7752f3824a86fa74d7f87c013eab4ede5e2d5a2..4baf5a6f8e8f3da5b2e2bf7bbcecef5e74b6072f 100644
--- a/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/SimpleAxisymmetricCreep.prj
+++ b/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/SimpleAxisymmetricCreep.prj
@@ -17,7 +17,7 @@
<q>Q</q>
<nonlinear_solver>
<maximum_iterations>1000</maximum_iterations>
- <error_tolerance>1e-8</error_tolerance>
+ <error_tolerance>1e-15</error_tolerance>
</nonlinear_solver>
</constitutive_relation>
<reference_solid_density>rho_sr</reference_solid_density>
@@ -29,14 +29,8 @@
<temperature>temperature</temperature>
</process_variables>
<secondary_variables>
- <secondary_variable type="static" internal_name="sigma_xx" output_name="sigma_xx"/>
- <secondary_variable type="static" internal_name="sigma_yy" output_name="sigma_yy"/>
- <secondary_variable type="static" internal_name="sigma_zz" output_name="sigma_zz"/>
- <secondary_variable type="static" internal_name="sigma_xy" output_name="sigma_xy"/>
- <secondary_variable type="static" internal_name="epsilon_xx" output_name="epsilon_xx"/>
- <secondary_variable type="static" internal_name="epsilon_yy" output_name="epsilon_yy"/>
- <secondary_variable type="static" internal_name="epsilon_zz" output_name="epsilon_zz"/>
- <secondary_variable type="static" internal_name="epsilon_xy" output_name="epsilon_xy"/>
+ <secondary_variable type="static" internal_name="sigma" output_name="sigma"/>
+ <secondary_variable type="static" internal_name="epsilon" output_name="epsilon"/>
</secondary_variables>
<specific_body_force>0 0</specific_body_force>
</process>
@@ -46,9 +40,9 @@
<process ref="ThermoMechanics">
<nonlinear_solver>basic_newton</nonlinear_solver>
<convergence_criterion>
- <type>DeltaX</type>
+ <type>PerComponentDeltaX</type>
<norm_type>NORM2</norm_type>
- <abstol>1e-10</abstol>
+ <abstols>5e-12 1e-16 1e-16</abstols>
</convergence_criterion>
<time_discretization>
<type>BackwardEuler</type>
@@ -57,14 +51,8 @@
<variables>
<variable>displacement</variable>
<variable>temperature</variable>
- <variable>sigma_xx</variable>
- <variable>sigma_yy</variable>
- <variable>sigma_zz</variable>
- <variable>sigma_xy</variable>
- <variable>epsilon_xx</variable>
- <variable>epsilon_yy</variable>
- <variable>epsilon_zz</variable>
- <variable>epsilon_xy</variable>
+ <variable>sigma</variable>
+ <variable>epsilon</variable>
</variables>
</output>
<time_stepping>
@@ -237,7 +225,7 @@
<solver_type>BiCGSTAB</solver_type>
<precon_type>DIAGONAL</precon_type>
<max_iteration_step>10000</max_iteration_step>
- <error_tolerance>1e-11</error_tolerance>
+ <error_tolerance>1e-17</error_tolerance>
</eigen>
</linear_solver>
</linear_solvers>
diff --git a/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/SimpleAxisymmetricCreepWithAnalyticSolution.prj b/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/SimpleAxisymmetricCreepWithAnalyticSolution.prj
index 7c66a11b8bb2af76efcf1496eb201b197d4f7a9b..3c4c9aaab1932885b5f377395a4376990b390b29 100644
--- a/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/SimpleAxisymmetricCreepWithAnalyticSolution.prj
+++ b/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/SimpleAxisymmetricCreepWithAnalyticSolution.prj
@@ -17,7 +17,7 @@
<q>Q</q>
<nonlinear_solver>
<maximum_iterations>1000</maximum_iterations>
- <error_tolerance>1e-8</error_tolerance>
+ <error_tolerance>1e-15</error_tolerance>
</nonlinear_solver>
</constitutive_relation>
<reference_solid_density>rho_sr</reference_solid_density>
@@ -29,14 +29,8 @@
<temperature>temperature</temperature>
</process_variables>
<secondary_variables>
- <secondary_variable type="static" internal_name="sigma_xx" output_name="sigma_xx"/>
- <secondary_variable type="static" internal_name="sigma_yy" output_name="sigma_yy"/>
- <secondary_variable type="static" internal_name="sigma_zz" output_name="sigma_zz"/>
- <secondary_variable type="static" internal_name="sigma_xy" output_name="sigma_xy"/>
- <secondary_variable type="static" internal_name="epsilon_xx" output_name="epsilon_xx"/>
- <secondary_variable type="static" internal_name="epsilon_yy" output_name="epsilon_yy"/>
- <secondary_variable type="static" internal_name="epsilon_zz" output_name="epsilon_zz"/>
- <secondary_variable type="static" internal_name="epsilon_xy" output_name="epsilon_xy"/>
+ <secondary_variable type="static" internal_name="sigma" output_name="sigma"/>
+ <secondary_variable type="static" internal_name="epsilon" output_name="epsilon"/>
</secondary_variables>
<specific_body_force>0 0</specific_body_force>
</process>
@@ -46,9 +40,9 @@
<process ref="ThermoMechanics">
<nonlinear_solver>basic_newton</nonlinear_solver>
<convergence_criterion>
- <type>DeltaX</type>
+ <type>PerComponentDeltaX</type>
<norm_type>NORM2</norm_type>
- <abstol>1e-9</abstol>
+ <abstols>5e-12 1e-16 1e-16</abstols>
</convergence_criterion>
<time_discretization>
<type>BackwardEuler</type>
@@ -57,14 +51,8 @@
<variables>
<variable>displacement</variable>
<variable>temperature</variable>
- <variable>sigma_xx</variable>
- <variable>sigma_yy</variable>
- <variable>sigma_zz</variable>
- <variable>sigma_xy</variable>
- <variable>epsilon_xx</variable>
- <variable>epsilon_yy</variable>
- <variable>epsilon_zz</variable>
- <variable>epsilon_xy</variable>
+ <variable>sigma</variable>
+ <variable>epsilon</variable>
</variables>
</output>
<time_stepping>
@@ -230,7 +218,7 @@
<solver_type>BiCGSTAB</solver_type>
<precon_type>DIAGONAL</precon_type>
<max_iteration_step>10000</max_iteration_step>
- <error_tolerance>1e-11</error_tolerance>
+ <error_tolerance>1e-17</error_tolerance>
</eigen>
</linear_solver>
</linear_solvers>
diff --git a/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/SimpleAxisymmetricCreepWithAnalyticSolution.vtu b/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/SimpleAxisymmetricCreepWithAnalyticSolution.vtu
index 4219698e0f48a19875417cacce64ceacb2076473..add6828ac1a7ff139e916cfa4af62df9ca494dfb 100644
--- a/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/SimpleAxisymmetricCreepWithAnalyticSolution.vtu
+++ b/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/SimpleAxisymmetricCreepWithAnalyticSolution.vtu
@@ -1,3 +1,3 @@
version https://git-lfs.github.com/spec/v1
-oid sha256:bc97d747694a6674b536b756548d1438b404ab4ee6d12cce4639fd3abaa1e036
-size 17588
+oid sha256:cc232efee4cda8c4485ed20bc32684de7ccdeb7b712b171cb3547eb7ecf7e67d
+size 6194
diff --git a/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/expected_SimpleAxisymmetricCreep_pcs_0_ts_370_t_360.000000.vtu b/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/expected_SimpleAxisymmetricCreep_pcs_0_ts_370_t_360.000000.vtu
index 84345b33ea1feca66cdcb75346c042b6d908efd3..9f451d777848baf9d5bd8a30e62a0efac9a884aa 100644
--- a/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/expected_SimpleAxisymmetricCreep_pcs_0_ts_370_t_360.000000.vtu
+++ b/Tests/Data/ThermoMechanics/CreepBGRa/SimpleAxisymmetricCreep/expected_SimpleAxisymmetricCreep_pcs_0_ts_370_t_360.000000.vtu
@@ -1,3 +1,3 @@
version https://git-lfs.github.com/spec/v1
-oid sha256:1ba0c28933aba289c38d4558b0f03d1213bf4eab27e24e571b16442f70459fa2
-size 22088
+oid sha256:b8af569b4f1b150cc06cb2ef85c7eddfe77ae6514bc8ec77458c602a67d7cc83
+size 21369
diff --git a/Tests/Data/ThermoMechanics/cube_1e3.prj b/Tests/Data/ThermoMechanics/cube_1e3.prj
index 3373f1f04d1e352a6716d5db046a8a325ddaf4b7..9fbe7b31c4f40c8207813752b2374a31f298ae59 100644
--- a/Tests/Data/ThermoMechanics/cube_1e3.prj
+++ b/Tests/Data/ThermoMechanics/cube_1e3.prj
@@ -21,14 +21,8 @@
<temperature>temperature</temperature>
</process_variables>
<secondary_variables>
- <secondary_variable type="static" internal_name="sigma_xx" output_name="sigma_xx"/>
- <secondary_variable type="static" internal_name="sigma_yy" output_name="sigma_yy"/>
- <secondary_variable type="static" internal_name="sigma_zz" output_name="sigma_zz"/>
- <secondary_variable type="static" internal_name="sigma_xy" output_name="sigma_xy"/>
- <secondary_variable type="static" internal_name="epsilon_xx" output_name="epsilon_xx"/>
- <secondary_variable type="static" internal_name="epsilon_yy" output_name="epsilon_yy"/>
- <secondary_variable type="static" internal_name="epsilon_zz" output_name="epsilon_zz"/>
- <secondary_variable type="static" internal_name="epsilon_xy" output_name="epsilon_xy"/>
+ <secondary_variable type="static" internal_name="sigma" output_name="sigma"/>
+ <secondary_variable type="static" internal_name="epsilon" output_name="epsilon"/>
</secondary_variables>
<specific_body_force>0 0 0</specific_body_force>
</process>
@@ -40,7 +34,7 @@
<convergence_criterion>
<type>DeltaX</type>
<norm_type>NORM2</norm_type>
- <abstol>1e-10</abstol>
+ <reltol>1e-15</reltol>
</convergence_criterion>
<time_discretization>
<type>BackwardEuler</type>
@@ -49,14 +43,8 @@
<variables>
<variable>displacement</variable>
<variable>temperature</variable>
- <variable>sigma_xx</variable>
- <variable>sigma_yy</variable>
- <variable>sigma_zz</variable>
- <variable>sigma_xy</variable>
- <variable>epsilon_xx</variable>
- <variable>epsilon_yy</variable>
- <variable>epsilon_zz</variable>
- <variable>epsilon_xy</variable>
+ <variable>sigma</variable>
+ <variable>epsilon</variable>
</variables>
</output>
<time_stepping>
@@ -82,7 +70,7 @@
<timesteps>
<pair>
<repeat>1</repeat>
- <each_steps>1</each_steps>
+ <each_steps>100000</each_steps>
</pair>
</timesteps>
</output>
@@ -228,7 +216,7 @@
<solver_type>BiCGSTAB</solver_type>
<precon_type>DIAGONAL</precon_type>
<max_iteration_step>10000</max_iteration_step>
- <error_tolerance>1e-11</error_tolerance>
+ <error_tolerance>1e-16</error_tolerance>
</eigen>
</linear_solver>
</linear_solvers>
diff --git a/Tests/Data/ThermoMechanics/cube_1e3_tm_pcs_0_ts_80_t_72000.000000.vtu b/Tests/Data/ThermoMechanics/cube_1e3_tm_pcs_0_ts_80_t_72000.000000.vtu
index 81c300cb5b2449da6c12baf8d75334bf0a14afbf..a21150b7e6ca4ef92020e163a99416eadd263f4b 100644
--- a/Tests/Data/ThermoMechanics/cube_1e3_tm_pcs_0_ts_80_t_72000.000000.vtu
+++ b/Tests/Data/ThermoMechanics/cube_1e3_tm_pcs_0_ts_80_t_72000.000000.vtu
@@ -1,3 +1,3 @@
version https://git-lfs.github.com/spec/v1
-oid sha256:2b3081955f0728044e72d5b64e7bb012667f3f8438a593bb254e5dcb7725e105
-size 133572
+oid sha256:0f435d15001708190703c067dcf3b7031634c4609dd28b4b7e5f138a94374898
+size 138620
diff --git a/Tests/Data/ThermoMechanics/expected_cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu b/Tests/Data/ThermoMechanics/expected_cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu
index 81c300cb5b2449da6c12baf8d75334bf0a14afbf..a21150b7e6ca4ef92020e163a99416eadd263f4b 100644
--- a/Tests/Data/ThermoMechanics/expected_cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu
+++ b/Tests/Data/ThermoMechanics/expected_cube_1e3_tm_pcs_0_ts_17_t_72000.000000.vtu
@@ -1,3 +1,3 @@
version https://git-lfs.github.com/spec/v1
-oid sha256:2b3081955f0728044e72d5b64e7bb012667f3f8438a593bb254e5dcb7725e105
-size 133572
+oid sha256:0f435d15001708190703c067dcf3b7031634c4609dd28b4b7e5f138a94374898
+size 138620
diff --git a/Tests/Data/ThermoMechanics/expected_tm_a_pcs_0_ts_20_t_20000.000000.vtu b/Tests/Data/ThermoMechanics/expected_tm_a_pcs_0_ts_20_t_20000.000000.vtu
index 043713412f8d6b115eb7b591d56413580c604ca9..0721cc18f6a32509479a3dc57271d41305a67279 100644
--- a/Tests/Data/ThermoMechanics/expected_tm_a_pcs_0_ts_20_t_20000.000000.vtu
+++ b/Tests/Data/ThermoMechanics/expected_tm_a_pcs_0_ts_20_t_20000.000000.vtu
@@ -1,3 +1,3 @@
version https://git-lfs.github.com/spec/v1
-oid sha256:75103bda0b248b181c3809d8738197bed765905d26161674d71db2a80ba8a321
-size 9522
+oid sha256:1136a16d828b5a18dd49ee3c373e56751a4f5392cfe7d4e3080069279fe3a7d1
+size 8795
diff --git a/Tests/Data/ThermoMechanics/expected_tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu b/Tests/Data/ThermoMechanics/expected_tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu
index 6ada2faf1548a530cf71af69a664c1a47475e24b..3e7218cf42200ea8af8f9cfe557d7ab1ebefc51b 100644
--- a/Tests/Data/ThermoMechanics/expected_tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu
+++ b/Tests/Data/ThermoMechanics/expected_tm_a_quad_pcs_0_ts_20_t_20000.000000.vtu
@@ -1,3 +1,3 @@
version https://git-lfs.github.com/spec/v1
-oid sha256:b7fed6966fb40ac53db83252c83b4062891de560201259fddfdadd291aac57a2
-size 22558
+oid sha256:d6227e63ed991c0b294e0fed6b9c2df395c023136d3258490f31c09c0c3d41dc
+size 21939
diff --git a/Tests/Data/ThermoMechanics/expected_tm_q_pcs_0_ts_20_t_20000.000000.vtu b/Tests/Data/ThermoMechanics/expected_tm_q_pcs_0_ts_20_t_20000.000000.vtu
index 989af0921228c9e4f39fe4320ec8f0a1f9d74c57..e48289533ad67c8bac4f25dec96aac23445f2ca4 100644
--- a/Tests/Data/ThermoMechanics/expected_tm_q_pcs_0_ts_20_t_20000.000000.vtu
+++ b/Tests/Data/ThermoMechanics/expected_tm_q_pcs_0_ts_20_t_20000.000000.vtu
@@ -1,3 +1,3 @@
version https://git-lfs.github.com/spec/v1
-oid sha256:5c5ab528497cb5654f7a646920abb3121b93647aa04cf520b0ca900f1c2fa192
-size 181403
+oid sha256:c8d3b7a7695280afee8630652048f2ec58df8c6129bd0f42c7c7b1f63fe4d262
+size 180820
diff --git a/Tests/Data/ThermoMechanics/expected_tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu b/Tests/Data/ThermoMechanics/expected_tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu
index 9418d398329db78113443094365a0b00b4de570c..dc9ef35addc5699c3ae4a914255e6cd47f6aed58 100644
--- a/Tests/Data/ThermoMechanics/expected_tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu
+++ b/Tests/Data/ThermoMechanics/expected_tm_q_quad_pcs_0_ts_20_t_20000.000000.vtu
@@ -1,3 +1,3 @@
version https://git-lfs.github.com/spec/v1
-oid sha256:c4c410e23ce8b496a90e88add0e6106318781290be6e4fe6172226a98e37cf03
-size 538964
+oid sha256:02cd7b7f1612d950e44f0acd11235876a228d8ef65ec54c37be0dab023ae8a45
+size 538804
diff --git a/Tests/Data/ThermoMechanics/iglu_axisymmetric_plane_strain.prj b/Tests/Data/ThermoMechanics/iglu_axisymmetric_plane_strain.prj
index 0c46219fe51c4d9aee611c34316547e7e0b30eeb..3c47888ff72cfea47d7c56171b0843415bdae639 100644
--- a/Tests/Data/ThermoMechanics/iglu_axisymmetric_plane_strain.prj
+++ b/Tests/Data/ThermoMechanics/iglu_axisymmetric_plane_strain.prj
@@ -21,14 +21,8 @@
<temperature>temperature</temperature>
</process_variables>
<secondary_variables>
- <secondary_variable type="static" internal_name="sigma_xx" output_name="sigma_xx"/>
- <secondary_variable type="static" internal_name="sigma_yy" output_name="sigma_yy"/>
- <secondary_variable type="static" internal_name="sigma_zz" output_name="sigma_zz"/>
- <secondary_variable type="static" internal_name="sigma_xy" output_name="sigma_xy"/>
- <secondary_variable type="static" internal_name="epsilon_xx" output_name="epsilon_xx"/>
- <secondary_variable type="static" internal_name="epsilon_yy" output_name="epsilon_yy"/>
- <secondary_variable type="static" internal_name="epsilon_zz" output_name="epsilon_zz"/>
- <secondary_variable type="static" internal_name="epsilon_xy" output_name="epsilon_xy"/>
+ <secondary_variable type="static" internal_name="sigma" output_name="sigma"/>
+ <secondary_variable type="static" internal_name="epsilon" output_name="epsilon"/>
</secondary_variables>
<specific_body_force>0 0</specific_body_force>
</process>
@@ -49,14 +43,8 @@
<variables>
<variable>displacement</variable>
<variable>temperature</variable>
- <variable>sigma_xx</variable>
- <variable>sigma_yy</variable>
- <variable>sigma_zz</variable>
- <variable>sigma_xy</variable>
- <variable>epsilon_xx</variable>
- <variable>epsilon_yy</variable>
- <variable>epsilon_zz</variable>
- <variable>epsilon_xy</variable>
+ <variable>sigma</variable>
+ <variable>epsilon</variable>
</variables>
</output>
<time_stepping>
diff --git a/Tests/Data/ThermoMechanics/iglu_axisymmetric_plane_strain_quad.prj b/Tests/Data/ThermoMechanics/iglu_axisymmetric_plane_strain_quad.prj
index 1b3497aa72c126544435859a7149d950c1bba3ad..92e216d176dcf67abfa379401ae6ae14ff5871d8 100644
--- a/Tests/Data/ThermoMechanics/iglu_axisymmetric_plane_strain_quad.prj
+++ b/Tests/Data/ThermoMechanics/iglu_axisymmetric_plane_strain_quad.prj
@@ -21,14 +21,8 @@
<temperature>temperature</temperature>
</process_variables>
<secondary_variables>
- <secondary_variable type="static" internal_name="sigma_xx" output_name="sigma_xx"/>
- <secondary_variable type="static" internal_name="sigma_yy" output_name="sigma_yy"/>
- <secondary_variable type="static" internal_name="sigma_zz" output_name="sigma_zz"/>
- <secondary_variable type="static" internal_name="sigma_xy" output_name="sigma_xy"/>
- <secondary_variable type="static" internal_name="epsilon_xx" output_name="epsilon_xx"/>
- <secondary_variable type="static" internal_name="epsilon_yy" output_name="epsilon_yy"/>
- <secondary_variable type="static" internal_name="epsilon_zz" output_name="epsilon_zz"/>
- <secondary_variable type="static" internal_name="epsilon_xy" output_name="epsilon_xy"/>
+ <secondary_variable type="static" internal_name="sigma" output_name="sigma"/>
+ <secondary_variable type="static" internal_name="epsilon" output_name="epsilon"/>
</secondary_variables>
<specific_body_force>0 0</specific_body_force>
</process>
@@ -49,14 +43,8 @@
<variables>
<variable>displacement</variable>
<variable>temperature</variable>
- <variable>sigma_xx</variable>
- <variable>sigma_yy</variable>
- <variable>sigma_zz</variable>
- <variable>sigma_xy</variable>
- <variable>epsilon_xx</variable>
- <variable>epsilon_yy</variable>
- <variable>epsilon_zz</variable>
- <variable>epsilon_xy</variable>
+ <variable>sigma</variable>
+ <variable>epsilon</variable>
</variables>
</output>
<time_stepping>
diff --git a/Tests/Data/ThermoMechanics/iglu_quarter_plane_strain.prj b/Tests/Data/ThermoMechanics/iglu_quarter_plane_strain.prj
index e48dc31792d354359b27f2948cd99c717f636122..52e17dcfe411ecf3df0e4aa02f58449ebc974f84 100644
--- a/Tests/Data/ThermoMechanics/iglu_quarter_plane_strain.prj
+++ b/Tests/Data/ThermoMechanics/iglu_quarter_plane_strain.prj
@@ -25,14 +25,8 @@
<temperature>temperature</temperature>
</process_variables>
<secondary_variables>
- <secondary_variable type="static" internal_name="sigma_xx" output_name="sigma_xx"/>
- <secondary_variable type="static" internal_name="sigma_yy" output_name="sigma_yy"/>
- <secondary_variable type="static" internal_name="sigma_zz" output_name="sigma_zz"/>
- <secondary_variable type="static" internal_name="sigma_xy" output_name="sigma_xy"/>
- <secondary_variable type="static" internal_name="epsilon_xx" output_name="epsilon_xx"/>
- <secondary_variable type="static" internal_name="epsilon_yy" output_name="epsilon_yy"/>
- <secondary_variable type="static" internal_name="epsilon_zz" output_name="epsilon_zz"/>
- <secondary_variable type="static" internal_name="epsilon_xy" output_name="epsilon_xy"/>
+ <secondary_variable type="static" internal_name="sigma" output_name="sigma"/>
+ <secondary_variable type="static" internal_name="epsilon" output_name="epsilon"/>
</secondary_variables>
<specific_body_force>0 0</specific_body_force>
</process>
@@ -53,14 +47,8 @@
<variables>
<variable>displacement</variable>
<variable>temperature</variable>
- <variable>sigma_xx</variable>
- <variable>sigma_yy</variable>
- <variable>sigma_zz</variable>
- <variable>sigma_xy</variable>
- <variable>epsilon_xx</variable>
- <variable>epsilon_yy</variable>
- <variable>epsilon_zz</variable>
- <variable>epsilon_xy</variable>
+ <variable>sigma</variable>
+ <variable>epsilon</variable>
</variables>
</output>
<time_stepping>
diff --git a/Tests/Data/ThermoMechanics/iglu_quarter_plane_strain_quad.prj b/Tests/Data/ThermoMechanics/iglu_quarter_plane_strain_quad.prj
index b9033f35c964fca01214b4a5386431d3a629d642..e6bf0cbde154c9702c62af5bda8ef6b67c1afd4a 100644
--- a/Tests/Data/ThermoMechanics/iglu_quarter_plane_strain_quad.prj
+++ b/Tests/Data/ThermoMechanics/iglu_quarter_plane_strain_quad.prj
@@ -25,14 +25,8 @@
<temperature>temperature</temperature>
</process_variables>
<secondary_variables>
- <secondary_variable type="static" internal_name="sigma_xx" output_name="sigma_xx"/>
- <secondary_variable type="static" internal_name="sigma_yy" output_name="sigma_yy"/>
- <secondary_variable type="static" internal_name="sigma_zz" output_name="sigma_zz"/>
- <secondary_variable type="static" internal_name="sigma_xy" output_name="sigma_xy"/>
- <secondary_variable type="static" internal_name="epsilon_xx" output_name="epsilon_xx"/>
- <secondary_variable type="static" internal_name="epsilon_yy" output_name="epsilon_yy"/>
- <secondary_variable type="static" internal_name="epsilon_zz" output_name="epsilon_zz"/>
- <secondary_variable type="static" internal_name="epsilon_xy" output_name="epsilon_xy"/>
+ <secondary_variable type="static" internal_name="sigma" output_name="sigma"/>
+ <secondary_variable type="static" internal_name="epsilon" output_name="epsilon"/>
</secondary_variables>
<specific_body_force>0 0</specific_body_force>
</process>
@@ -53,14 +47,8 @@
<variables>
<variable>displacement</variable>
<variable>temperature</variable>
- <variable>sigma_xx</variable>
- <variable>sigma_yy</variable>
- <variable>sigma_zz</variable>
- <variable>sigma_xy</variable>
- <variable>epsilon_xx</variable>
- <variable>epsilon_yy</variable>
- <variable>epsilon_zz</variable>
- <variable>epsilon_xy</variable>
+ <variable>sigma</variable>
+ <variable>epsilon</variable>
</variables>
</output>
<time_stepping>
diff --git a/Tests/Data/ThermoMechanics/stress_analytical.vtu b/Tests/Data/ThermoMechanics/stress_analytical.vtu
index 70017b976b85097d62a2f4da5e140250d3eebe40..f3e85ef2c13a985877a9264551e09e68810d9a67 100644
--- a/Tests/Data/ThermoMechanics/stress_analytical.vtu
+++ b/Tests/Data/ThermoMechanics/stress_analytical.vtu
@@ -1,3 +1,3 @@
version https://git-lfs.github.com/spec/v1
-oid sha256:242c7bacb711c10cd0b058b40d7436a6a887b588f834f40d1dad92e22e73258c
-size 318250
+oid sha256:2f8801e5ccd2048d4cfa45d736bc3e46b04000c3d0062230e97f4218acdca1e6
+size 22086
diff --git a/web/content/docs/benchmarks/creepbgra/CreepBRGa.pandoc b/web/content/docs/benchmarks/creepbgra/CreepBRGa.pandoc
index 76bd21ac23c81bec8570649874c7470c4be02512..83ca08206ea5177a04c9dbe1fd314935963dc5dc 100644
--- a/web/content/docs/benchmarks/creepbgra/CreepBRGa.pandoc
+++ b/web/content/docs/benchmarks/creepbgra/CreepBRGa.pandoc
@@ -33,7 +33,7 @@ To establish equivalence between the three-dimensional and the uniaxial formulat
we use the effective stress defined by
${ \bar\sigma}={\sqrt{{\frac{3}{2}}}}{\left\Vert{\mathbf s}\right\Vert}$
with ${\mathbf s}= { \mathbf \sigma}-\frac{1}{3}{ \mathrm{tr}(\mathbf\sigma}){\mathbf I}$,
-the deviatoric stress.
+the deviatoric stress.
The creep strain rate is then expressed as
$$\begin{gathered}
\dot { \mathbf \epsilon}^c ({ \sigma})= {\dfrac{\partial g^c}{\partial {\bar\sigma}}}
@@ -48,7 +48,7 @@ the creep rate equation to a uniaxial stress state ${\mathbf \sigma} = \mathrm{d
$$\begin{gathered}
\dot { \epsilon_1}^c = {\dfrac{\partial g^c}{\partial { \bar\sigma}}}=Ae^{-Q/R_uT}\left(\dfrac{{ \sigma_1}}{{ \sigma}_f}\right)^m
\end{gathered}$$
-which says
+which says
$$\begin{gathered}
{\dfrac{\partial g^c}{\partial { \bar\sigma}}}=Ae^{-Q/R_u T}\left(\dfrac{{ \sigma_1}}{{ \sigma}_f}\right)^m
\end{gathered}$$
@@ -72,21 +72,21 @@ $$\begin{gathered}
\dot { \mathbf \epsilon}= \dot { \mathbf \epsilon}^e + \dot { \mathbf \epsilon}^T
+ \dot { \mathbf \epsilon}^c
\end{gathered}$$ Due to Hook’s law, the stress rate is
-given by
+given by
$$\begin{gathered}
\dot { \mathbf \sigma}= \mathbf{C} \dot { \mathbf \epsilon}^e = \mathbf{C}
(\dot { \mathbf \epsilon}- \dot { \mathbf \epsilon}^T- \dot { \mathbf \epsilon}^c)
-\end{gathered}$$
+\end{gathered}$$
where
-$\mathbf{C}:= \lambda \mathcal{J} + 2G \mathbf I \otimes \mathbf I $
+$\mathbf{C}:= \lambda \mathcal{J} + 2G \mathbf I \otimes \mathbf I $
with $\mathcal{J}$ the forth order identity, $\mathbf I$ the second order identity,
$\lambda$ the Lamé constant, $G$ the shear modulus, and $\otimes$ the tensor
product notation.
is a fourth order tensor. Substituting equation and the expression of $C$
-into the stress rate expression, equation , we have
+into the stress rate expression, equation, we have
$$\begin{gathered}
-\dot { \mathbf \sigma}= \mathbf{C} (\dot { \mathbf \epsilon}- \dot { \mathbf \epsilon}^T)- 2bG {\left\Vert{\mathbf s}\right\Vert}^{m-1}
+\dot { \mathbf \sigma}= \mathbf{C} (\dot { \mathbf \epsilon}- \dot { \mathbf \epsilon}^T)- 2bG {\left\Vert{\mathbf s}\right\Vert}^{m-1}
{\mathbf s}
\end{gathered}$$
@@ -105,12 +105,12 @@ Newton-Raphson is applied to .
{ \mathbf \sigma}^{n} - \mathbf{C} (\Delta { \mathbf \epsilon} - \alpha_T \Delta T \mathbf I)
+ 2bG \Delta t {\left\Vert{\mathbf s}^{n+1}\right\Vert}^{m-1}
{\mathbf s}^{n+1}
-\end{gathered}$$
+\end{gathered}$$
be the residual. The Jacobian
of is derived as
$$\begin{gathered}
- \mathbf{J}_{{ \mathbf \sigma}}= {\dfrac{\partial \mathbf{r}}{\partial { \mathbf \sigma}^{n+1}}} = \mathcal{I} + 2bG\Delta t {\left\Vert{\mathbf s}^{n+1}\right\Vert}^{m-1}
+ \mathbf{J}_{{ \mathbf \sigma}}= {\dfrac{\partial \mathbf{r}}{\partial { \mathbf \sigma}^{n+1}}} = \mathcal{I} + 2bG\Delta t {\left\Vert{\mathbf s}^{n+1}\right\Vert}^{m-1}
\left(\mathcal{I} - \frac{1}{3} \mathbf{I} \otimes \mathbf{I} + (m-1){\left\Vert{\mathbf s}^{n+1}\right\Vert}^{-2} {\mathbf s}^{n+1}\otimes {\mathbf s}^{n+1}\right)
\end{gathered}$$
@@ -120,10 +120,10 @@ Consistent tangent
Once the converged stress integration is obtained, the tangential of
stress with respect to strain can be obtained straightforwardly by
applying the partial derivative to equation with respect to strain
-increment as
+increment as
$$\begin{aligned}
{\dfrac{\partial { \mathbf \sigma}^{n+1}}{\partial \Delta { \mathbf \epsilon}}} = \mathbf{C} - {\dfrac{\partial \left(2bG \Delta t {\left\Vert{\mathbf s}^{n+1}\right\Vert}^{m-1} {\mathbf s}^{n+1}\right)}{\partial \Delta { \mathbf \epsilon}}}
- \end{aligned}$$
+ \end{aligned}$$
We see that $$\begin{aligned}
{\dfrac{\partial { \mathbf \sigma}^{n+1}}{\partial { \mathbf \epsilon}^{n+1}}} = \mathbf{J}_{{ \mathbf \sigma}}^{-1} \mathbf{C}
\end{aligned}$$
@@ -136,7 +136,7 @@ of the global Jacobian can be derived as
\end{aligned}$$
*Note*: The above rate form of stress integration is implemented in ogs6.
- Alternatively, one can use a absolute stress integration form, which can be found in the attached
+ Alternatively, one can use a absolute stress integration form, which can be found in the attached
[PDF](../doku_BGRa.pdf).
Example
@@ -147,7 +147,7 @@ pressure of ${ \sigma}_0=5 ~\mbox{MPa}$. The values of the parameters are
given as $A=0.18 \mbox{d}^{-1}$, $m=5$, $Q=54 ~\mbox{kJ/mol}$,
$E=25000 ~\mbox{MPa}$, and the temperature is constant everywhere of
100 $^{\circ}\mbox{C}$. The analytical solution of the strain is given
-straightforward as
+straightforward as
$$\begin{gathered}
{ \epsilon}=-\dfrac{{ \sigma}_0}{E}-Ae^{-Q/RT}{ \sigma}_0^m t
\end{gathered}$$
@@ -183,3 +183,49 @@ $$\begin{gathered}
obtained by the present multidimensional scheme with the analytical
solution is shown in the following figure:
{{< img src="../bgra0.png" >}}
+
+
+Python code
+-----------
+A short python snippet, to compute the values.
+<details>
+<summary>
+Insert this into Paraview's ProgrammableFilter:
+</summary>
+```python
+A = self.GetInputDataObject(0, 0)
+numPoints = A.GetNumberOfPoints()
+outSxx = vtk.vtkDoubleArray()
+outSxx.SetName('analytic_strain_xx')
+outSyy = vtk.vtkDoubleArray()
+outSyy.SetName('analytic_strain_yy')
+outSzz = vtk.vtkDoubleArray()
+outSzz.SetName('analytic_strain_zz')
+outSxy = vtk.vtkDoubleArray()
+outSxy.SetName('analytic_strain_xy')
+sigma0=5.0
+A=0.18
+m=5
+Q=54000
+E=25000
+R=8.314472
+sf=1.
+t=100
+e0=-sigma0/E
+c=-A*exp(-Q/(R*373.15))*pow(sigma0,m)
+strain_zz=e0+c*t
+nv=0.27
+G=0.5*E/(1+nv)
+strain_rr=strain_zz + 0.5*sigma0/G+1.5*A*exp(-Q/(R*373.15))*pow(sigma0, m)*t
+for i in range(0, numPoints):
+ outSxx.InsertNextValue(strain_rr)
+ outSyy.InsertNextValue(strain_zz)
+ outSzz.InsertNextValue(strain_rr)
+ outSxy.InsertNextValue(0)
+output = self.GetOutput()
+output.GetPointData().AddArray(outSxx)
+output.GetPointData().AddArray(outSyy)
+output.GetPointData().AddArray(outSzz)
+output.GetPointData().AddArray(outSxy)
+```
+</details>
diff --git a/web/content/docs/benchmarks/thermo-mechanics/thermomechanics.pandoc b/web/content/docs/benchmarks/thermo-mechanics/thermomechanics.pandoc
index 07643749c5f0cb0fd4be9d2dd9f4ca5075de8ce9..db1ec79db78b84efba715d2e2704c37007d7202c 100644
--- a/web/content/docs/benchmarks/thermo-mechanics/thermomechanics.pandoc
+++ b/web/content/docs/benchmarks/thermo-mechanics/thermomechanics.pandoc
@@ -15,23 +15,27 @@ weight = 156
## Problem description
-We solve a thermo-mechanical homogeneous model in cube domain. The dimensions of this cube model are 1\,m in all directions. The boundary conditions and temperature loadings, as well as the material can refer Chapter 14 in Kolditz et al. for detailed problem description.
+We solve a thermo-mechanical homogeneous model in cube domain. The dimensions of
+this cube model are 1 m in all directions. The boundary conditions and
+temperature loadings, as well as the material can refer Chapter 14 in Kolditz et
+al. \cite Kolditz2012 for detailed problem description.
## Results and evaluation
-Result showing temperature and stresses development with time in the centre node of the model:
+Result showing temperature and stresses development with time in the centre node
+of the model:
{{< img src="../temperature.png" >}}
{{< img src="../stress.png" >}}
The analytical solution of stresses after heating is:
-$$
-\begin{equation}
-\sigma_{xx} = \sigma_{yy} = \sigma_{zz} = - \frac{\alpha \Delta T E}{1 - 2 \nu} = - 3.260869\, \mathrm{MPa}
-\end{equation}
-$$
+$$\begin{equation}
+\sigma_{xx} = \sigma_{yy} = \sigma_{zz} = - \frac{\alpha \Delta T E}{1 - 2 \nu}
+= - 3.260869\, \textrm{MPa}
+\end{equation}$$
-The relative error between the numerical simulation and the analytical solution is $9.2 \cdot 10^{-13}$.
+The relative error between the numerical simulation and the analytical solution
+is 9.2<span class="math inline">â‹…10<sup>-13</sup></span>.
## References