From 5cd47d34c4883396f305673f6a2ebba39d54380e Mon Sep 17 00:00:00 2001
From: labourasse <labourassee@gmail.com>
Date: Thu, 16 Jan 2025 16:19:32 +0100
Subject: [PATCH] fix a bug in limiter and br calculations
---
src/scheme/Order2HyperelasticSolver.cpp | 650 ++++++++++++------------
1 file changed, 339 insertions(+), 311 deletions(-)
diff --git a/src/scheme/Order2HyperelasticSolver.cpp b/src/scheme/Order2HyperelasticSolver.cpp
index 83b802f51..e771f81b8 100644
--- a/src/scheme/Order2HyperelasticSolver.cpp
+++ b/src/scheme/Order2HyperelasticSolver.cpp
@@ -29,7 +29,8 @@
#include <vector>
template <MeshConcept MeshType>
-class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public Order2HyperelasticSolverHandler::IOrder2HyperelasticSolver
+class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final
+ : public Order2HyperelasticSolverHandler::IOrder2HyperelasticSolver
{
private:
static constexpr size_t Dimension = MeshType::Dimension;
@@ -190,7 +191,7 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(mesh);
const NodeValuePerCell<const Rd>& Cjr = mesh_data.Cjr();
- const auto& xr = mesh.xr();
+ const auto& xr = mesh.xr();
const auto& node_to_cell_matrix = mesh.connectivity().nodeToCellMatrix();
const auto& node_local_numbers_in_their_cells = mesh.connectivity().nodeLocalNumbersInTheirCells();
@@ -206,7 +207,8 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
for (size_t j = 0; j < node_to_cell.size(); ++j) {
const CellId J = node_to_cell[j];
const unsigned int R = node_local_number_in_its_cells[j];
- br += Ajr(J, R) * DPk_u[J](xr[r]) - DPk_sigma[j](xr[r]) * Cjr(J, R);
+ // br += Ajr(J, R) * DPk_u[J](xr[r]) - DPk_sigma[j](xr[r]) * Cjr(J, R);
+ br += Ajr(J, R) * DPk_u[J](xr[r]) - DPk_sigma[J](xr[r]) * Cjr(J, R);
}
b[r] = br;
@@ -337,13 +339,13 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
_computeFjr(const MeshType& mesh,
const NodeValuePerCell<const Rdxd>& Ajr,
const NodeValue<const Rd>& ur,
- DiscreteFunctionDPk<Dimension,const Rd> DPk_uh,
- DiscreteFunctionDPk<Dimension,const Rdxd> DPk_sigma) const
+ DiscreteFunctionDPk<Dimension, const Rd> DPk_uh,
+ DiscreteFunctionDPk<Dimension, const Rdxd> DPk_sigma) const
{
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(mesh);
const NodeValuePerCell<const Rd> Cjr = mesh_data.Cjr();
- const NodeValue<const Rd>& xr = mesh.xr();
+ const NodeValue<const Rd>& xr = mesh.xr();
const auto& node_to_cell_matrix = mesh.connectivity().nodeToCellMatrix();
const auto& node_local_numbers_in_their_cells = mesh.connectivity().nodeLocalNumbersInTheirCells();
@@ -354,17 +356,17 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
const auto& node_to_cell = node_to_cell_matrix[r];
const auto& node_local_number_in_its_cells = node_local_numbers_in_their_cells.itemArray(r);
- for(size_t j = 0; j < node_to_cell.size(); ++j){
+ for (size_t j = 0; j < node_to_cell.size(); ++j) {
const CellId J = node_to_cell[j];
const unsigned int R = node_local_number_in_its_cells[j];
- F(J,R) = -Ajr(J,R) * (DPk_uh[J](xr[r]) - ur[r]) + DPk_sigma[J](xr[r]) * Cjr(J,R);
+ F(J, R) = -Ajr(J, R) * (DPk_uh[J](xr[r]) - ur[r]) + DPk_sigma[J](xr[r]) * Cjr(J, R);
}
});
return F;
}
- void
+ void
_tensor_limiter(const MeshType& mesh,
const DiscreteFunctionP0<const Rdxd>& f,
DiscreteFunctionDPk<Dimension, Rdxd>& DPk_fh) const
@@ -372,83 +374,85 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
MeshData<MeshType>& mesh_data = MeshDataManager::instance().getMeshData(mesh);
StencilManager::BoundaryDescriptorList symmetry_boundary_descriptor_list;
StencilDescriptor stencil_descriptor{1, StencilDescriptor::ConnectionType::by_nodes};
- auto stencil = StencilManager::instance().getCellToCellStencilArray(mesh.connectivity(), stencil_descriptor, symmetry_boundary_descriptor_list);
+ auto stencil = StencilManager::instance().getCellToCellStencilArray(mesh.connectivity(), stencil_descriptor,
+ symmetry_boundary_descriptor_list);
const auto xj = mesh_data.xj();
- parallel_for(mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id){
- const Rdxd fj = f[cell_id];
-
- Rdxd f_min = fj;
- Rdxd f_max = fj;
-
- const auto cell_stencil = stencil[cell_id];
- for(size_t i_cell = 0; i_cell < cell_stencil.size(); ++i_cell){
- for(size_t row = 0; row < Dimension; ++row){
- for(size_t col = 0; col < Dimension; ++col){
- f_min(row,col) = std::min(f_min(row,col), f[cell_stencil[i_cell]](row,col));
- f_max(row,col) = std::max(f_max(row,col), f[cell_stencil[i_cell]](row,col));
+ parallel_for(
+ mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id) {
+ const Rdxd fj = f[cell_id];
+
+ Rdxd f_min = fj;
+ Rdxd f_max = fj;
+
+ const auto cell_stencil = stencil[cell_id];
+ for (size_t i_cell = 0; i_cell < cell_stencil.size(); ++i_cell) {
+ for (size_t row = 0; row < Dimension; ++row) {
+ for (size_t col = 0; col < Dimension; ++col) {
+ f_min(row, col) = std::min(f_min(row, col), f[cell_stencil[i_cell]](row, col));
+ f_max(row, col) = std::max(f_max(row, col), f[cell_stencil[i_cell]](row, col));
+ }
}
}
- }
- Rdxd f_bar_min = fj;
- Rdxd f_bar_max = fj;
+ Rdxd f_bar_min = fj;
+ Rdxd f_bar_max = fj;
- for(size_t i_cell = 0; i_cell < cell_stencil.size(); ++i_cell){
- const CellId cell_k_id = cell_stencil[i_cell];
- const Rdxd f_xk = DPk_fh[cell_id](xj[cell_k_id]);
+ for (size_t i_cell = 0; i_cell < cell_stencil.size(); ++i_cell) {
+ const CellId cell_k_id = cell_stencil[i_cell];
+ const Rdxd f_xk = DPk_fh[cell_id](xj[cell_k_id]);
- for(size_t row = 0; row < Dimension; ++row){
- for(size_t col = 0; col < Dimension; ++col){
- f_bar_min(row,col) = std::min(f_bar_min(row,col), f_xk(row,col));
- f_bar_max(row,col) = std::max(f_bar_max(row,col), f_xk(row,col));
+ for (size_t row = 0; row < Dimension; ++row) {
+ for (size_t col = 0; col < Dimension; ++col) {
+ f_bar_min(row, col) = std::min(f_bar_min(row, col), f_xk(row, col));
+ f_bar_max(row, col) = std::max(f_bar_max(row, col), f_xk(row, col));
+ }
}
}
- }
- const double eps = 1E-14;
- Rdxd coef1;
- for(size_t row = 0; row < Dimension; ++row){
- for(size_t col = 0; col < Dimension; ++col){
- coef1(row,col) = 1;
- if(std::abs(f_bar_min(row,col) - fj(row,col)) > eps){
- coef1(row,col) = (f_max(row,col) - fj(row,col)) / (f_bar_max(row,col) - fj(row,col));
+ const double eps = 1E-14;
+ Rdxd coef1;
+ for (size_t row = 0; row < Dimension; ++row) {
+ for (size_t col = 0; col < Dimension; ++col) {
+ coef1(row, col) = 1;
+ if (std::abs(f_bar_max(row, col) - fj(row, col)) > eps) {
+ coef1(row, col) = (f_max(row, col) - fj(row, col)) / (f_bar_max(row, col) - fj(row, col));
+ }
}
}
- }
- Rdxd coef2;
- for(size_t row = 0; row < Dimension; ++row){
- for(size_t col = 0; col < Dimension; ++col){
- coef2(row,col) = 1;
- if (std::abs(f_bar_min(row,col) - fj(row,col)) > eps) {
- coef2(row,col) = (f_min(row,col) - fj(row,col)) / ((f_bar_min(row,col) - fj(row,col)));
+ Rdxd coef2;
+ for (size_t row = 0; row < Dimension; ++row) {
+ for (size_t col = 0; col < Dimension; ++col) {
+ coef2(row, col) = 1;
+ if (std::abs(f_bar_min(row, col) - fj(row, col)) > eps) {
+ coef2(row, col) = (f_min(row, col) - fj(row, col)) / ((f_bar_min(row, col) - fj(row, col)));
+ }
}
}
- }
- double min_coef1 = coef1(0,0);
- double min_coef2 = coef2(0,0);
- for(size_t row = 0; row < Dimension; ++row){
- for(size_t col = 0; col < Dimension; ++col){
- min_coef1 = std::min(min_coef1,coef1(row,col));
- min_coef2 = std::min(min_coef2,coef2(row,col));
+ double min_coef1 = coef1(0, 0);
+ double min_coef2 = coef2(0, 0);
+ for (size_t row = 0; row < Dimension; ++row) {
+ for (size_t col = 0; col < Dimension; ++col) {
+ min_coef1 = std::min(min_coef1, coef1(row, col));
+ min_coef2 = std::min(min_coef2, coef2(row, col));
+ }
}
- }
- const double lambda = std::max(0., std::min(1.,std::min(min_coef1, min_coef2)));
+ const double lambda = std::max(0., std::min(1., std::min(min_coef1, min_coef2)));
- auto coefficients = DPk_fh.coefficients(cell_id);
+ auto coefficients = DPk_fh.coefficients(cell_id);
- coefficients[0] = (1 - lambda) * f[cell_id] + lambda * coefficients[0];
- for (size_t i = 1; i < coefficients.size(); ++i) {
- coefficients[i] *= lambda;
- }
- });
+ coefficients[0] = (1 - lambda) * f[cell_id] + lambda * coefficients[0];
+ for (size_t i = 1; i < coefficients.size(); ++i) {
+ coefficients[i] *= lambda;
+ }
+ });
}
- void
+ void
_vector_limiter(const MeshType& mesh,
const DiscreteFunctionP0<const Rd>& f,
DiscreteFunctionDPk<Dimension, Rd>& DPk_fh) const
@@ -456,73 +460,75 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
MeshData<MeshType>& mesh_data = MeshDataManager::instance().getMeshData(mesh);
StencilManager::BoundaryDescriptorList symmetry_boundary_descriptor_list;
StencilDescriptor stencil_descriptor{1, StencilDescriptor::ConnectionType::by_nodes};
- auto stencil = StencilManager::instance().getCellToCellStencilArray(mesh.connectivity(), stencil_descriptor, symmetry_boundary_descriptor_list);
+ auto stencil = StencilManager::instance().getCellToCellStencilArray(mesh.connectivity(), stencil_descriptor,
+ symmetry_boundary_descriptor_list);
const auto xj = mesh_data.xj();
- parallel_for(mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id){
- const Rd fj = f[cell_id];
+ parallel_for(
+ mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id) {
+ const Rd fj = f[cell_id];
- Rd f_min = fj;
- Rd f_max = fj;
+ Rd f_min = fj;
+ Rd f_max = fj;
- const auto cell_stencil = stencil[cell_id];
- for(size_t i_cell = 0; i_cell < cell_stencil.size(); ++i_cell){
- for(size_t dim = 0; dim < Dimension; ++dim){
- f_min[dim] = std::min(f_min[dim], f[cell_stencil[i_cell]][dim]);
- f_max[dim] = std::max(f_max[dim], f[cell_stencil[i_cell]][dim]);
+ const auto cell_stencil = stencil[cell_id];
+ for (size_t i_cell = 0; i_cell < cell_stencil.size(); ++i_cell) {
+ for (size_t dim = 0; dim < Dimension; ++dim) {
+ f_min[dim] = std::min(f_min[dim], f[cell_stencil[i_cell]][dim]);
+ f_max[dim] = std::max(f_max[dim], f[cell_stencil[i_cell]][dim]);
+ }
}
- }
- Rd f_bar_min = fj;
- Rd f_bar_max = fj;
+ Rd f_bar_min = fj;
+ Rd f_bar_max = fj;
- for(size_t i_cell = 0; i_cell < cell_stencil.size(); ++i_cell){
- const CellId cell_k_id = cell_stencil[i_cell];
- const Rd f_xk = DPk_fh[cell_id](xj[cell_k_id]);
+ for (size_t i_cell = 0; i_cell < cell_stencil.size(); ++i_cell) {
+ const CellId cell_k_id = cell_stencil[i_cell];
+ const Rd f_xk = DPk_fh[cell_id](xj[cell_k_id]);
- for(size_t dim = 0; dim < Dimension; ++dim){
- f_bar_min[dim] = std::min(f_bar_min[dim], f_xk[dim]);
- f_bar_max[dim] = std::max(f_bar_max[dim], f_xk[dim]);
+ for (size_t dim = 0; dim < Dimension; ++dim) {
+ f_bar_min[dim] = std::min(f_bar_min[dim], f_xk[dim]);
+ f_bar_max[dim] = std::max(f_bar_max[dim], f_xk[dim]);
+ }
}
- }
- const double eps = 1E-14;
- Rd coef1;
- for(size_t dim = 0; dim < Dimension; ++dim){
- coef1[dim] = 1;
- if (std::abs(f_bar_max[dim] - fj[dim]) > eps) {
- coef1[dim] = (f_max[dim] - fj[dim]) / ((f_bar_max[dim] - fj[dim]));
+ const double eps = 1E-14;
+ Rd coef1;
+ for (size_t dim = 0; dim < Dimension; ++dim) {
+ coef1[dim] = 1;
+ if (std::abs(f_bar_max[dim] - fj[dim]) > eps) {
+ coef1[dim] = (f_max[dim] - fj[dim]) / ((f_bar_max[dim] - fj[dim]));
+ }
}
- }
- Rd coef2;
- for(size_t dim = 0; dim < Dimension; ++dim){
- coef2[dim] = 1;
- if (std::abs(f_bar_min[dim] - fj[dim]) > eps) {
- coef2[dim] = (f_min[dim] - fj[dim]) / ((f_bar_min[dim] - fj[dim]));
+ Rd coef2;
+ for (size_t dim = 0; dim < Dimension; ++dim) {
+ coef2[dim] = 1;
+ if (std::abs(f_bar_min[dim] - fj[dim]) > eps) {
+ coef2[dim] = (f_min[dim] - fj[dim]) / ((f_bar_min[dim] - fj[dim]));
+ }
}
- }
- double min_coef1 = coef1[0];
- double min_coef2 = coef2[0];
- for(size_t dim = 1; dim < Dimension; ++dim){
- min_coef1 = std::min(min_coef1,coef1[dim]);
- min_coef2 = std::min(min_coef2,coef2[dim]);
- }
+ double min_coef1 = coef1[0];
+ double min_coef2 = coef2[0];
+ for (size_t dim = 1; dim < Dimension; ++dim) {
+ min_coef1 = std::min(min_coef1, coef1[dim]);
+ min_coef2 = std::min(min_coef2, coef2[dim]);
+ }
- const double lambda = std::max(0., std::min(1., std::min(min_coef1, min_coef2)));
+ const double lambda = std::max(0., std::min(1., std::min(min_coef1, min_coef2)));
- auto coefficients = DPk_fh.coefficients(cell_id);
+ auto coefficients = DPk_fh.coefficients(cell_id);
- coefficients[0] = (1 - lambda) * f[cell_id] + lambda * coefficients[0];
- for (size_t i = 1; i < coefficients.size(); ++i) {
- coefficients[i] *= lambda;
- }
- });
+ coefficients[0] = (1 - lambda) * f[cell_id] + lambda * coefficients[0];
+ for (size_t i = 1; i < coefficients.size(); ++i) {
+ coefficients[i] *= lambda;
+ }
+ });
}
- CellValue<double>
+ CellValue<double>
_scalar_limiter_coefficient(const MeshType& mesh,
const DiscreteFunctionP0<const double>& f,
const DiscreteFunctionDPk<Dimension, const double>& DPk_fh) const
@@ -530,80 +536,86 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
MeshData<MeshType>& mesh_data = MeshDataManager::instance().getMeshData(mesh);
StencilManager::BoundaryDescriptorList symmetry_boundary_descriptor_list;
StencilDescriptor stencil_descriptor{1, StencilDescriptor::ConnectionType::by_nodes};
- auto stencil = StencilManager::instance().getCellToCellStencilArray(mesh.connectivity(), stencil_descriptor, symmetry_boundary_descriptor_list);
+ auto stencil = StencilManager::instance().getCellToCellStencilArray(mesh.connectivity(), stencil_descriptor,
+ symmetry_boundary_descriptor_list);
const auto xj = mesh_data.xj();
CellValue<double> lambda{mesh.connectivity()};
- parallel_for(mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id){
- const double fj = f[cell_id];
+ parallel_for(
+ mesh.numberOfCells(), PUGS_LAMBDA(const CellId cell_id) {
+ const double fj = f[cell_id];
- double f_min = fj;
- double f_max = fj;
+ double f_min = fj;
+ double f_max = fj;
- const auto cell_stencil = stencil[cell_id];
- for(size_t i_cell = 0; i_cell < cell_stencil.size(); ++i_cell){
- f_min = std::min(f_min, f[cell_stencil[i_cell]]);
- f_max = std::max(f_max, f[cell_stencil[i_cell]]);
- }
+ const auto cell_stencil = stencil[cell_id];
+ for (size_t i_cell = 0; i_cell < cell_stencil.size(); ++i_cell) {
+ f_min = std::min(f_min, f[cell_stencil[i_cell]]);
+ f_max = std::max(f_max, f[cell_stencil[i_cell]]);
+ }
- double f_bar_min = fj;
- double f_bar_max = fj;
+ double f_bar_min = fj;
+ double f_bar_max = fj;
- for(size_t i_cell = 0; i_cell < cell_stencil.size(); ++i_cell){
- const CellId cell_k_id = cell_stencil[i_cell];
- const double f_xk = DPk_fh[cell_id](xj[cell_k_id]);
+ for (size_t i_cell = 0; i_cell < cell_stencil.size(); ++i_cell) {
+ const CellId cell_k_id = cell_stencil[i_cell];
+ const double f_xk = DPk_fh[cell_id](xj[cell_k_id]);
- f_bar_min = std::min(f_bar_min, f_xk);
- f_bar_max = std::max(f_bar_max, f_xk);
- }
+ f_bar_min = std::min(f_bar_min, f_xk);
+ f_bar_max = std::max(f_bar_max, f_xk);
+ }
- const double eps = 1E-14;
- double coef1 = 1;
- if (std::abs(f_bar_max - fj) > eps) {
- coef1 = (f_max - fj) / ((f_bar_max - fj));
- }
+ const double eps = 1E-14;
+ double coef1 = 1;
+ if (std::abs(f_bar_max - fj) > eps) {
+ coef1 = (f_max - fj) / ((f_bar_max - fj));
+ }
- double coef2 = 1.;
- if (std::abs(f_bar_min - fj) > eps) {
- coef2 = (f_min - fj) / ((f_bar_min - fj));
- }
+ double coef2 = 1.;
+ if (std::abs(f_bar_min - fj) > eps) {
+ coef2 = (f_min - fj) / ((f_bar_min - fj));
+ }
- lambda[cell_id] = std::max(0., std::min(1., std::min(coef1, coef2)));
- });
+ lambda[cell_id] = std::max(0., std::min(1., std::min(coef1, coef2)));
+ });
return lambda;
}
- std::tuple<const std::shared_ptr<const DiscreteFunctionVariant>, const std::shared_ptr<const DiscreteFunctionVariant>, const std::shared_ptr<const DiscreteFunctionVariant>>
- _apply_NeoHook(const DiscreteScalarFunction rho,
- const DiscreteVectorFunction u,
- const DiscreteScalarFunction E,
- const DiscreteTensorFunction CG,
- const DiscreteScalarFunction chi_solid,
- const double& lambda,
- const double& mu,
- const double& gamma,
- const double& p_inf) const
+ std::tuple<const std::shared_ptr<const DiscreteFunctionVariant>,
+ const std::shared_ptr<const DiscreteFunctionVariant>,
+ const std::shared_ptr<const DiscreteFunctionVariant>>
+ _apply_NeoHook(const DiscreteScalarFunction rho,
+ const DiscreteVectorFunction u,
+ const DiscreteScalarFunction E,
+ const DiscreteTensorFunction CG,
+ const DiscreteScalarFunction chi_solid,
+ const double& lambda,
+ const double& mu,
+ const double& gamma,
+ const double& p_inf) const
{
const TinyMatrix<Dimension> I = identity;
- const DiscreteScalarFunction& eps = E - 0.5 * dot(u,u);
- const DiscreteScalarFunction& p_d = (1-chi_solid)*(gamma-1) * rho * eps - p_inf*gamma;
- const DiscreteTensorFunction& sigma_d = chi_solid * (mu / sqrt(det(CG)) * (CG - I) + lambda / sqrt(det(CG)) * log(sqrt(det(CG))) * I) - p_d * I;
- const DiscreteScalarFunction& aL_d = sqrt(chi_solid * (lambda + 2 * mu) / rho + gamma * (p_d + p_inf) / rho);
- const DiscreteScalarFunction& aT_d = sqrt(chi_solid * mu / rho);
+ const DiscreteScalarFunction& eps = E - 0.5 * dot(u, u);
+ const DiscreteScalarFunction& p_d = (1 - chi_solid) * (gamma - 1) * rho * eps - p_inf * gamma;
+ const DiscreteTensorFunction& sigma_d =
+ chi_solid * (mu / sqrt(det(CG)) * (CG - I) + lambda / sqrt(det(CG)) * log(sqrt(det(CG))) * I) - p_d * I;
+ const DiscreteScalarFunction& aL_d = sqrt(chi_solid * (lambda + 2 * mu) / rho + gamma * (p_d + p_inf) / rho);
+ const DiscreteScalarFunction& aT_d = sqrt(chi_solid * mu / rho);
- return {std::make_shared<DiscreteFunctionVariant>(sigma_d),
- std::make_shared<DiscreteFunctionVariant>(aL_d),
+ return {std::make_shared<DiscreteFunctionVariant>(sigma_d), std::make_shared<DiscreteFunctionVariant>(aL_d),
std::make_shared<DiscreteFunctionVariant>(aT_d)};
}
- std::tuple<const std::shared_ptr<const DiscreteFunctionVariant>, const std::shared_ptr<const DiscreteFunctionVariant>, const std::shared_ptr<const DiscreteFunctionVariant>>
- _apply_NeoHook2(const DiscreteScalarFunction rho,
- const DiscreteVectorFunction u,
- const DiscreteScalarFunction E,
- const DiscreteTensorFunction CG,
+ std::tuple<const std::shared_ptr<const DiscreteFunctionVariant>,
+ const std::shared_ptr<const DiscreteFunctionVariant>,
+ const std::shared_ptr<const DiscreteFunctionVariant>>
+ _apply_NeoHook2(const DiscreteScalarFunction rho,
+ const DiscreteVectorFunction u,
+ const DiscreteScalarFunction E,
+ const DiscreteTensorFunction CG,
const DiscreteScalarFunction chi_solid,
const double& mu,
const double& gamma_f,
@@ -613,18 +625,22 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
{
const TinyMatrix<Dimension> I = identity;
- const DiscreteScalarFunction& eps = E - 0.5 * dot(u,u);
- const DiscreteScalarFunction& p_d = (1-chi_solid)*((gamma_f-1) * rho * eps)+ chi_solid * ((gamma_s - 1.) * rho * eps - gamma_s * p_inf_s);
- const DiscreteTensorFunction& sigma_d = chi_solid * 2*mu/sqrt(det(CG))*(CG-1./3. * trace(CG)*I) - p_d * I;
- const DiscreteScalarFunction& aL_d = sqrt(chi_solid * (2 * mu) / rho + ((1-chi_solid)*gamma_f * p_d + chi_solid * (gamma_s*(p_d+p_inf_s))) / rho);
- const DiscreteScalarFunction& aT_d = sqrt(chi_solid * mu / rho);
+ const DiscreteScalarFunction& eps = E - 0.5 * dot(u, u);
+ const DiscreteScalarFunction& p_d =
+ (1 - chi_solid) * ((gamma_f - 1) * rho * eps) + chi_solid * ((gamma_s - 1.) * rho * eps - gamma_s * p_inf_s);
+ const DiscreteTensorFunction& sigma_d =
+ chi_solid * 2 * mu / sqrt(det(CG)) * (CG - 1. / 3. * trace(CG) * I) - p_d * I;
+ const DiscreteScalarFunction& aL_d = sqrt(
+ chi_solid * (2 * mu) / rho + ((1 - chi_solid) * gamma_f * p_d + chi_solid * (gamma_s * (p_d + p_inf_s))) / rho);
+ const DiscreteScalarFunction& aT_d = sqrt(chi_solid * mu / rho);
- return {std::make_shared<DiscreteFunctionVariant>(sigma_d),
- std::make_shared<DiscreteFunctionVariant>(aL_d),
+ return {std::make_shared<DiscreteFunctionVariant>(sigma_d), std::make_shared<DiscreteFunctionVariant>(aL_d),
std::make_shared<DiscreteFunctionVariant>(aT_d)};
}
- std::tuple<const std::shared_ptr<const DiscreteFunctionVariant>, const std::shared_ptr<const DiscreteFunctionVariant>, const std::shared_ptr<const DiscreteFunctionVariant>>
+ std::tuple<const std::shared_ptr<const DiscreteFunctionVariant>,
+ const std::shared_ptr<const DiscreteFunctionVariant>,
+ const std::shared_ptr<const DiscreteFunctionVariant>>
_apply_state_law(const size_t law,
const std::shared_ptr<const DiscreteFunctionVariant>& rho_v,
const std::shared_ptr<const DiscreteFunctionVariant>& u_v,
@@ -636,7 +652,7 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
const double& gamma_f,
const double& p_inf_f,
const double& gamma_s,
- const double& p_inf_s) const
+ const double& p_inf_s) const
{
const DiscreteScalarFunction& chi_solid_d = chi_solid_v->get<DiscreteScalarFunction>();
const DiscreteScalarFunction& rho_d = rho_v->get<DiscreteScalarFunction>();
@@ -644,14 +660,14 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
const DiscreteScalarFunction& E_d = E_v->get<DiscreteScalarFunction>();
const DiscreteTensorFunction& CG_d = CG_v->get<DiscreteTensorFunction>();
- //const std::shared_ptr<const DiscreteFunctionVariant> sigma;
- //const std::shared_ptr<const DiscreteFunctionVariant> aL;
- //const std::shared_ptr<const DiscreteFunctionVariant> aT;
+ // const std::shared_ptr<const DiscreteFunctionVariant> sigma;
+ // const std::shared_ptr<const DiscreteFunctionVariant> aL;
+ // const std::shared_ptr<const DiscreteFunctionVariant> aT;
- if(law == 1){
- return _apply_NeoHook(rho_d,u_d,E_d,CG_d,chi_solid_d,lambda,mu,gamma_f,p_inf_f);
+ if (law == 1) {
+ return _apply_NeoHook(rho_d, u_d, E_d, CG_d, chi_solid_d, lambda, mu, gamma_f, p_inf_f);
} else {
- return _apply_NeoHook2(rho_d,u_d,E_d,CG_d,chi_solid_d,mu,gamma_f,p_inf_f,gamma_s,p_inf_s);
+ return _apply_NeoHook2(rho_d, u_d, E_d, CG_d, chi_solid_d, mu, gamma_f, p_inf_f, gamma_s, p_inf_s);
}
}
@@ -668,8 +684,7 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
std::shared_ptr mesh_v = getCommonMesh({rho_v, aL_v, aT_v, u_v, sigma_v});
if (not mesh_v) {
throw NormalError("discrete functions are not defined on the same mesh");
- } // - p* identity; // - chi_solid*P_zero_air_repo*identity;
-
+ } // - p* identity; // - chi_solid*P_zero_air_repo*identity;
if (not checkDiscretizationType({rho_v, aL_v, u_v, sigma_v}, DiscreteFunctionType::P0)) {
throw NormalError("hyperelastic solver expects P0 functions");
@@ -684,118 +699,124 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
std::vector<std::shared_ptr<const IBoundaryDescriptor>> symmetry_boundary_descriptor_list;
- for(auto&& bc_descriptor : bc_descriptor_list){
- if(bc_descriptor->type() == IBoundaryConditionDescriptor::Type::symmetry){
+ for (auto&& bc_descriptor : bc_descriptor_list) {
+ if (bc_descriptor->type() == IBoundaryConditionDescriptor::Type::symmetry) {
symmetry_boundary_descriptor_list.push_back(bc_descriptor->boundaryDescriptor_shared());
}
}
PolynomialReconstructionDescriptor reconstruction_descriptor(IntegrationMethodType::cell_center, 1,
- symmetry_boundary_descriptor_list);
+ symmetry_boundary_descriptor_list);
- //CellValue<double> limiter_j{mesh.connectivity()};
- //parallel_for(
- // mesh.numberOfCells(), PUGS_LAMBDA(CellId j){
- // limiter_j[j] = 1;
- // }
+ // CellValue<double> limiter_j{mesh.connectivity()};
+ // parallel_for(
+ // mesh.numberOfCells(), PUGS_LAMBDA(CellId j){
+ // limiter_j[j] = 1;
+ // }
//);
//
- //const CellValue<const Rdxd>& sigma_j = copy(sigma.cellValues());
- //std::vector<DiscreteFunctionDPk<Dimension, double>> sigma_coef;
- //std::vector<size_t> row;
- //std::vector<size_t> col;
- //for(size_t i = 0; i < Dimension; ++i){
- // for(size_t k = i; k < Dimension; ++k){
+ // const CellValue<const Rdxd>& sigma_j = copy(sigma.cellValues());
+ // std::vector<DiscreteFunctionDPk<Dimension, double>> sigma_coef;
+ // std::vector<size_t> row;
+ // std::vector<size_t> col;
+ // for(size_t i = 0; i < Dimension; ++i){
+ // for(size_t k = i; k < Dimension; ++k){
//
- // CellValue<double> coef{mesh.connectivity()};
- // for(CellId j = 0; j <mesh.numberOfCells(); ++j){
- // coef[j] = sigma_j[j](i,k);
- // }
+ // CellValue<double> coef{mesh.connectivity()};
+ // for(CellId j = 0; j <mesh.numberOfCells(); ++j){
+ // coef[j] = sigma_j[j](i,k);
+ // }
//
- // const auto& coef_scalar_function = DiscreteScalarFunction(mesh_v, coef);
- // auto coef_v = std::make_shared<DiscreteFunctionVariant>(coef_scalar_function);
- // auto reconstruction = PolynomialReconstruction{reconstruction_descriptor}.build(coef_v);
- // auto DPk_coef = reconstruction[0]->get<DiscreteFunctionDPk<Dimension, const double>>();
+ // const auto& coef_scalar_function = DiscreteScalarFunction(mesh_v, coef);
+ // auto coef_v = std::make_shared<DiscreteFunctionVariant>(coef_scalar_function);
+ // auto reconstruction = PolynomialReconstruction{reconstruction_descriptor}.build(coef_v);
+ // auto DPk_coef = reconstruction[0]->get<DiscreteFunctionDPk<Dimension, const double>>();
//
- // const CellValue<const double>& limiter_ik = _scalar_limiter_coefficient(mesh,coef_scalar_function,DPk_coef);
- // parallel_for(
- // mesh.numberOfCells(), PUGS_LAMBDA(CellId j) {
- // limiter_j[j] = std::min(limiter_j[j],limiter_ik[j]);
- // });
+ // const CellValue<const double>& limiter_ik =
+ // _scalar_limiter_coefficient(mesh,coef_scalar_function,DPk_coef); parallel_for(
+ // mesh.numberOfCells(), PUGS_LAMBDA(CellId j) {
+ // limiter_j[j] = std::min(limiter_j[j],limiter_ik[j]);
+ // });
//
- // sigma_coef.push_back(copy(DPk_coef));
- // row.push_back(i);
- // col.push_back(k);
- // }
- //}
+ // sigma_coef.push_back(copy(DPk_coef));
+ // row.push_back(i);
+ // col.push_back(k);
+ // }
+ // }
//
- //NodeValuePerCell<Rdxd> sigma_lim{mesh.connectivity()} ;
- //const auto& cell_to_node_matrix = mesh.connectivity().cellToNodeMatrix();
+ // NodeValuePerCell<Rdxd> sigma_lim{mesh.connectivity()} ;
+ // const auto& cell_to_node_matrix = mesh.connectivity().cellToNodeMatrix();
//
- //parallel_for(
- // mesh.numberOfCells(), PUGS_LAMBDA(CellId j){
- // const auto& cell_nodes = cell_to_node_matrix[j];
- // for(size_t i = 0; i < Dimension; ++i){
- // for(size_t k = 0; k < Dimension; ++k){
- // for(size_t R = 0; R < cell_nodes.size(); ++R){
- // sigma_lim(j,R)(i,k) = 0;
- // }
- // }
- // }
- // }
+ // parallel_for(
+ // mesh.numberOfCells(), PUGS_LAMBDA(CellId j){
+ // const auto& cell_nodes = cell_to_node_matrix[j];
+ // for(size_t i = 0; i < Dimension; ++i){
+ // for(size_t k = 0; k < Dimension; ++k){
+ // for(size_t R = 0; R < cell_nodes.size(); ++R){
+ // sigma_lim(j,R)(i,k) = 0;
+ // }
+ // }
+ // }
+ // }
//);
//
- //const NodeValue<const Rd>& xr = copy(mesh.xr());
- //parallel_for(
- // mesh.numberOfCells(), PUGS_LAMBDA(CellId j) {
- // const auto& cell_nodes = cell_to_node_matrix[j];
+ // const NodeValue<const Rd>& xr = copy(mesh.xr());
+ // parallel_for(
+ // mesh.numberOfCells(), PUGS_LAMBDA(CellId j) {
+ // const auto& cell_nodes = cell_to_node_matrix[j];
//
- // for(size_t R = 0; R < cell_nodes.size(); ++R){
- // const NodeId r = cell_nodes[R];
+ // for(size_t R = 0; R < cell_nodes.size(); ++R){
+ // const NodeId r = cell_nodes[R];
//
- // for(size_t l = 0; l < sigma_coef.size(); ++l){
- // const size_t& i = row[l];
- // const size_t& k = col[l];
+ // for(size_t l = 0; l < sigma_coef.size(); ++l){
+ // const size_t& i = row[l];
+ // const size_t& k = col[l];
//
- // auto coefficients = sigma_coef[l].coefficients(j);
+ // auto coefficients = sigma_coef[l].coefficients(j);
//
- // coefficients[0] = (1-limiter_j[j])*sigma_j[j](i,k) + limiter_j[j] * coefficients[0];
- // for (size_t indice = 1; indice < coefficients.size(); ++indice) {
- // coefficients[indice] *= limiter_j[j];
- // }
- //
- // sigma_lim(j,R)(i,k) = sigma_coef[l][j](xr[r]);
- // if(i != k){
- // sigma_lim(j,R)(i,k) *= 2;
- // }
- // }
+ // coefficients[0] = (1-limiter_j[j])*sigma_j[j](i,k) + limiter_j[j] * coefficients[0];
+ // for (size_t indice = 1; indice < coefficients.size(); ++indice) {
+ // coefficients[indice] *= limiter_j[j];
+ // }
//
- // sigma_lim(j,R) += transpose(sigma_lim(j,R));
- // sigma_lim(j,R) *= 0.5;
- // }
+ // sigma_lim(j,R)(i,k) = sigma_coef[l][j](xr[r]);
+ // if(i != k){
+ // sigma_lim(j,R)(i,k) *= 2;
+ // }
+ // }
//
- // }
+ // sigma_lim(j,R) += transpose(sigma_lim(j,R));
+ // sigma_lim(j,R) *= 0.5;
+ // }
+ //
+ // }
//);
- std::cout << "Reconstruction" << "\n";
- auto reconstruction = PolynomialReconstruction{reconstruction_descriptor}.build(u_v,sigma_v);
- auto DPk_uh = reconstruction[0]->get<DiscreteFunctionDPk<Dimension, const Rd>>();
- auto DPk_sigmah = reconstruction[1]->get<DiscreteFunctionDPk<Dimension, const Rdxd>>();
- std::cout << "Reconstruction Ok ! " << "\n";
-
- std::cout << "Limitation" << "\n";
- DiscreteFunctionDPk<Dimension, Rd> u_lim = copy(DPk_uh);
- DiscreteFunctionDPk<Dimension, Rdxd> sigma_lim = copy(DPk_sigmah);
- this->_vector_limiter(mesh,u,u_lim);
- this->_tensor_limiter(mesh,sigma,sigma_lim);
- std::cout << "Limitation OK !" << "\n";
+ std::cout << "Reconstruction"
+ << "\n";
+ auto reconstruction = PolynomialReconstruction{reconstruction_descriptor}.build(u_v, sigma_v);
+ auto DPk_uh = reconstruction[0]->get<DiscreteFunctionDPk<Dimension, const Rd>>();
+ auto DPk_sigmah = reconstruction[1]->get<DiscreteFunctionDPk<Dimension, const Rdxd>>();
+ std::cout << "Reconstruction Ok ! "
+ << "\n";
+
+ std::cout << "Limitation"
+ << "\n";
+ DiscreteFunctionDPk<Dimension, Rd> u_lim = copy(DPk_uh);
+ DiscreteFunctionDPk<Dimension, Rdxd> sigma_lim = copy(DPk_sigmah);
+ this->_vector_limiter(mesh, u, u_lim);
+ this->_tensor_limiter(mesh, sigma, sigma_lim);
+ std::cout << "Limitation OK !"
+ << "\n";
NodeValuePerCell<const Rdxd> Ajr = this->_computeAjr(solver_type, mesh, rho * aL, rho * aT);
NodeValue<Rdxd> Ar = this->_computeAr(mesh, Ajr);
- std::cout << "br" << "\n";
- NodeValue<Rd> br = this->_computeBr(mesh, Ajr, u_lim, sigma_lim);
- std::cout << "br Ok !" << "\n";
+ std::cout << "br"
+ << "\n";
+ NodeValue<Rd> br = this->_computeBr(mesh, Ajr, u_lim, sigma_lim);
+ std::cout << "br Ok !"
+ << "\n";
const BoundaryConditionList bc_list = this->_getBCList(mesh, bc_descriptor_list);
this->_applyBoundaryConditions(bc_list, mesh, Ar, br);
@@ -803,10 +824,12 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
synchronize(Ar);
synchronize(br);
- NodeValue<const Rd> ur = this->_computeUr(mesh, Ar, br);
- std::cout << "Fjr" << "\n";
+ NodeValue<const Rd> ur = this->_computeUr(mesh, Ar, br);
+ std::cout << "Fjr"
+ << "\n";
NodeValuePerCell<const Rd> Fjr = this->_computeFjr(mesh, Ajr, ur, u_lim, sigma_lim);
- std::cout << "Fjr Ok ! " << "\n";
+ std::cout << "Fjr Ok ! "
+ << "\n";
return std::make_tuple(std::make_shared<const ItemValueVariant>(ur),
std::make_shared<const SubItemValuePerItemVariant>(Fjr));
@@ -920,15 +943,15 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>>
order2_apply_fluxes(const double& dt,
- const MeshType& mesh,
- const DiscreteScalarFunction& rho,
- const DiscreteVectorFunction& u,
- const DiscreteScalarFunction& E,
- const DiscreteTensorFunction& CG,
- const DiscreteScalarFunction& rho_app,
- const DiscreteTensorFunction& CG_app,
- const NodeValue<const Rd>& ur,
- const NodeValuePerCell<const Rd>& Fjr) const
+ const MeshType& mesh,
+ const DiscreteScalarFunction& rho,
+ const DiscreteVectorFunction& u,
+ const DiscreteScalarFunction& E,
+ const DiscreteTensorFunction& CG,
+ const DiscreteScalarFunction& rho_app,
+ const DiscreteTensorFunction& CG_app,
+ const NodeValue<const Rd>& ur,
+ const NodeValuePerCell<const Rd>& Fjr) const
{
const auto& cell_to_node_matrix = mesh.connectivity().cellToNodeMatrix();
@@ -991,14 +1014,14 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>>
order2_apply_fluxes(const double& dt,
- const std::shared_ptr<const DiscreteFunctionVariant>& rho,
- const std::shared_ptr<const DiscreteFunctionVariant>& u,
- const std::shared_ptr<const DiscreteFunctionVariant>& E,
- const std::shared_ptr<const DiscreteFunctionVariant>& CG,
- const std::shared_ptr<const DiscreteFunctionVariant>& rho_app,
- const std::shared_ptr<const DiscreteFunctionVariant>& CG_app,
- const std::shared_ptr<const ItemValueVariant>& ur,
- const std::shared_ptr<const SubItemValuePerItemVariant>& Fjr) const
+ const std::shared_ptr<const DiscreteFunctionVariant>& rho,
+ const std::shared_ptr<const DiscreteFunctionVariant>& u,
+ const std::shared_ptr<const DiscreteFunctionVariant>& E,
+ const std::shared_ptr<const DiscreteFunctionVariant>& CG,
+ const std::shared_ptr<const DiscreteFunctionVariant>& rho_app,
+ const std::shared_ptr<const DiscreteFunctionVariant>& CG_app,
+ const std::shared_ptr<const ItemValueVariant>& ur,
+ const std::shared_ptr<const SubItemValuePerItemVariant>& Fjr) const
{
std::shared_ptr mesh_v = getCommonMesh({rho, u, E});
if (not mesh_v) {
@@ -1009,16 +1032,16 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
throw NormalError("hyperelastic solver expects P0 functions");
}
- return this->order2_apply_fluxes(dt, //
- *mesh_v->get<MeshType>(), //
- rho->get<DiscreteScalarFunction>(), //
- u->get<DiscreteVectorFunction>(), //
- E->get<DiscreteScalarFunction>(), //
- CG->get<DiscreteTensorFunction>(), //
- rho_app->get<DiscreteScalarFunction>(), //
- CG_app->get<DiscreteTensorFunction>(), //
- ur->get<NodeValue<const Rd>>(), //
- Fjr->get<NodeValuePerCell<const Rd>>());
+ return this->order2_apply_fluxes(dt, //
+ *mesh_v->get<MeshType>(), //
+ rho->get<DiscreteScalarFunction>(), //
+ u->get<DiscreteVectorFunction>(), //
+ E->get<DiscreteScalarFunction>(), //
+ CG->get<DiscreteTensorFunction>(), //
+ rho_app->get<DiscreteScalarFunction>(), //
+ CG_app->get<DiscreteTensorFunction>(), //
+ ur->get<NodeValue<const Rd>>(), //
+ Fjr->get<NodeValuePerCell<const Rd>>());
}
std::tuple<std::shared_ptr<const MeshVariant>,
@@ -1045,33 +1068,35 @@ class Order2HyperelasticSolverHandler::Order2HyperelasticSolver final : public O
const double& gamma_s,
const double& p_inf_s) const
{
- std::shared_ptr m_app = getCommonMesh({rho, aL, aT, u, sigma});
+ std::shared_ptr m_app = getCommonMesh({rho, aL, aT, u, sigma});
std::shared_ptr<const DiscreteFunctionVariant> rho_app = rho;
std::shared_ptr<const DiscreteFunctionVariant> u_app = u;
std::shared_ptr<const DiscreteFunctionVariant> E_app = E;
std::shared_ptr<const DiscreteFunctionVariant> CG_app = CG;
auto [ur, Fjr] = compute_fluxes(solver_type, rho, aL, aT, u, sigma, bc_descriptor_list);
- //auto [m_app, rho_app, u_app, E_app, CG_app] = apply_fluxes(dt/2., rho, u, E, CG, ur, Fjr);
+ // auto [m_app, rho_app, u_app, E_app, CG_app] = apply_fluxes(dt/2., rho, u, E, CG, ur, Fjr);
- //std::shared_ptr<const DiscreteFunctionVariant> chi_solid_app = shallowCopy(m_app, chi_solid);
+ // std::shared_ptr<const DiscreteFunctionVariant> chi_solid_app = shallowCopy(m_app, chi_solid);
- //auto [sigma_app, aL_app, aT_app] = _apply_state_law(law,rho_app,u_app,E_app,CG_app,chi_solid_app,lambda,mu,gamma_f,p_inf_f,gamma_s,p_inf_s);
+ // auto [sigma_app, aL_app, aT_app] =
+ // _apply_state_law(law,rho_app,u_app,E_app,CG_app,chi_solid_app,lambda,mu,gamma_f,p_inf_f,gamma_s,p_inf_s);
- //auto [ur_app, Fjr_app] = compute_fluxes(solver_type, rho_app, aL_app,aT_app,u_app,sigma_app,bc_descriptor_list);
+ // auto [ur_app, Fjr_app] = compute_fluxes(solver_type, rho_app, aL_app,aT_app,u_app,sigma_app,bc_descriptor_list);
return apply_fluxes(dt, rho, u, E, CG, ur, Fjr);
}
- Order2HyperelasticSolver() = default;
+ Order2HyperelasticSolver() = default;
Order2HyperelasticSolver(Order2HyperelasticSolver&&) = default;
- ~Order2HyperelasticSolver() = default;
+ ~Order2HyperelasticSolver() = default;
};
template <MeshConcept MeshType>
void
-Order2HyperelasticSolverHandler::Order2HyperelasticSolver<MeshType>::_applyPressureBC(const BoundaryConditionList& bc_list,
- const MeshType& mesh,
- NodeValue<Rd>& br) const
+Order2HyperelasticSolverHandler::Order2HyperelasticSolver<MeshType>::_applyPressureBC(
+ const BoundaryConditionList& bc_list,
+ const MeshType& mesh,
+ NodeValue<Rd>& br) const
{
for (const auto& boundary_condition : bc_list) {
std::visit(
@@ -1134,9 +1159,10 @@ Order2HyperelasticSolverHandler::Order2HyperelasticSolver<MeshType>::_applyPress
template <MeshConcept MeshType>
void
-Order2HyperelasticSolverHandler::Order2HyperelasticSolver<MeshType>::_applyNormalStressBC(const BoundaryConditionList& bc_list,
- const MeshType& mesh,
- NodeValue<Rd>& br) const
+Order2HyperelasticSolverHandler::Order2HyperelasticSolver<MeshType>::_applyNormalStressBC(
+ const BoundaryConditionList& bc_list,
+ const MeshType& mesh,
+ NodeValue<Rd>& br) const
{
for (const auto& boundary_condition : bc_list) {
std::visit(
@@ -1199,9 +1225,10 @@ Order2HyperelasticSolverHandler::Order2HyperelasticSolver<MeshType>::_applyNorma
template <MeshConcept MeshType>
void
-Order2HyperelasticSolverHandler::Order2HyperelasticSolver<MeshType>::_applySymmetryBC(const BoundaryConditionList& bc_list,
- NodeValue<Rdxd>& Ar,
- NodeValue<Rd>& br) const
+Order2HyperelasticSolverHandler::Order2HyperelasticSolver<MeshType>::_applySymmetryBC(
+ const BoundaryConditionList& bc_list,
+ NodeValue<Rdxd>& Ar,
+ NodeValue<Rd>& br) const
{
for (const auto& boundary_condition : bc_list) {
std::visit(
@@ -1230,9 +1257,10 @@ Order2HyperelasticSolverHandler::Order2HyperelasticSolver<MeshType>::_applySymme
template <MeshConcept MeshType>
void
-Order2HyperelasticSolverHandler::Order2HyperelasticSolver<MeshType>::_applyVelocityBC(const BoundaryConditionList& bc_list,
- NodeValue<Rdxd>& Ar,
- NodeValue<Rd>& br) const
+Order2HyperelasticSolverHandler::Order2HyperelasticSolver<MeshType>::_applyVelocityBC(
+ const BoundaryConditionList& bc_list,
+ NodeValue<Rdxd>& Ar,
+ NodeValue<Rd>& br) const
{
for (const auto& boundary_condition : bc_list) {
std::visit(
--
GitLab