Commit 82a4e5ff authored by GILLES Sebastien's avatar GILLES Sebastien
Browse files

#531 Hyperelastic operator now also uses up the new assembling procedure....

#531 Hyperelastic operator now also uses up the new assembling procedure. Previous commit has been checked successfully: the computation is right whatever the ordering inside the tuple.
parent f989345d
......@@ -135,20 +135,29 @@ namespace HappyHeart
* for Assemble() method rather than figuring out what is a variadic method and which additional
* arguments are required (none for this specific operator).
*/
void Assemble(double matrix_coefficient, GlobalMatrix& global_matrix,
double vector_coefficient, GlobalVector& global_vector,
// void Assemble(double matrix_coefficient, GlobalMatrix& global_matrix,
// double vector_coefficient, GlobalVector& global_vector,
// const GlobalVector& previous_iteration_data,
// const Domain& domain = Domain()) const;
//
// //! Same as overload above except only a matrix is assembled.
// void Assemble(double matrix_coefficient, GlobalMatrix& global_matrix,
// const GlobalVector& previous_iteration_data,
// const Domain& domain = Domain()) const;
//
// //! Same as overload above except only a vector is assembled.
// void Assemble(double vector_coefficient, GlobalVector& global_vector,
// const GlobalVector& previous_iteration_data,
// const Domain& domain = Domain()) const;
template<class LinearAlgebraTupleT>
void Assemble(LinearAlgebraTupleT&& linear_algebra_tuple,
const GlobalVector& previous_iteration_data,
const Domain& domain = Domain()) const;
//! Same as overload above except only a matrix is assembled.
void Assemble(double matrix_coefficient, GlobalMatrix& global_matrix,
const GlobalVector& previous_iteration_data,
const Domain& domain = Domain()) const;
//! Same as overload above except only a vector is assembled.
void Assemble(double vector_coefficient, GlobalVector& global_vector,
const GlobalVector& previous_iteration_data,
const Domain& domain = Domain()) const;
private:
......
......@@ -52,41 +52,53 @@ namespace HappyHeart
}
template<class LawPolicyT>
inline void
GradOnGradientBasedHyperelasticityTensor<LawPolicyT>
::Assemble(double matrix_coefficient, GlobalMatrix& global_matrix,
double vector_coefficient, GlobalVector& global_vector,
const GlobalVector& input_vector,
const Domain& domain) const
{
return Parent::AssembleWithVariadicArguments(matrix_coefficient, global_matrix,
vector_coefficient, global_vector,
domain, input_vector);
}
template<class LawPolicyT>
inline void
GradOnGradientBasedHyperelasticityTensor<LawPolicyT>
::Assemble(double matrix_coefficient, GlobalMatrix& global_matrix,
const GlobalVector& input_vector,
const Domain& domain) const
{
return Parent::AssembleWithVariadicArguments(matrix_coefficient, global_matrix,
domain, input_vector);
}
// template<class LawPolicyT>
// inline void
// GradOnGradientBasedHyperelasticityTensor<LawPolicyT>
// ::Assemble(double matrix_coefficient, GlobalMatrix& global_matrix,
// double vector_coefficient, GlobalVector& global_vector,
// const GlobalVector& input_vector,
// const Domain& domain) const
// {
// return Parent::AssembleWithVariadicArguments(matrix_coefficient, global_matrix,
// vector_coefficient, global_vector,
// domain, input_vector);
// }
//
//
// template<class LawPolicyT>
// inline void
// GradOnGradientBasedHyperelasticityTensor<LawPolicyT>
// ::Assemble(double matrix_coefficient, GlobalMatrix& global_matrix,
// const GlobalVector& input_vector,
// const Domain& domain) const
// {
// return Parent::AssembleWithVariadicArguments(matrix_coefficient, global_matrix,
// domain, input_vector);
// }
//
//
// template<class LawPolicyT>
// inline void
// GradOnGradientBasedHyperelasticityTensor<LawPolicyT>
// ::Assemble(double vector_coefficient, GlobalVector& global_vector,
// const GlobalVector& input_vector,
// const Domain& domain) const
// {
// return Parent::AssembleWithVariadicArguments(vector_coefficient, global_vector,
// domain, input_vector);
// }
//
template<class LawPolicyT>
template<class LinearAlgebraTupleT>
inline void
GradOnGradientBasedHyperelasticityTensor<LawPolicyT>
::Assemble(double vector_coefficient, GlobalVector& global_vector,
::Assemble(LinearAlgebraTupleT&& linear_algebra_tuple,
const GlobalVector& input_vector,
const Domain& domain) const
{
return Parent::AssembleWithVariadicArguments(vector_coefficient, global_vector,
domain, input_vector);
return Parent::Assemble531(std::move(linear_algebra_tuple), domain, input_vector);
}
......
......@@ -122,17 +122,13 @@ namespace HappyHeart
* \tparam MatrixVectorNatureT Depending on this parameter, the local matrix and/or the local vector is
* computed by the method.
*
* \param[in] matrix_coefficient If the local matrix is computed, this coefficient is applied.
* \param[in] vector_coefficient If the local vector is computed, this coefficient is applied.
* \param[in] local_displacement Last known displacement for the current finite element.
* \internal This parameter is computed by
* GlobalVariationalOperatorNS::GradOnGradientBasedHyperelasticityTensor<LawPolicyT>::SupplArgumentsComputeEltArray() \
* method.
*/
template<MatrixVectorNature MatrixVectorNatureT>
void ComputeEltArray(double matrix_coefficient,
double vector_coefficient,
const std::vector<double>& local_displacement);
void ComputeEltArray(const std::vector<double>& local_displacement);
private:
......@@ -145,9 +141,7 @@ namespace HappyHeart
* \brief Part of ComputeEltArray once all the dimension depending operations have been performed.
*/
template<MatrixVectorNature MatrixVectorNatureT>
void ComputeAtQuadraturePoint(double matrix_coefficient,
double vector_coefficient,
const Private::ShapeFunctionAtQuadraturePoint& infos_at_quad_pt,
void ComputeAtQuadraturePoint(const Private::ShapeFunctionAtQuadraturePoint& infos_at_quad_pt,
const LocalMatrix& tangent_matrix,
const LocalVector& rhs_part,
ElementaryData<MatrixVectorNature::vector_and_matrix>& elementary_data) const;
......
......@@ -142,9 +142,7 @@ namespace HappyHeart
template <class LawPolicyT>
template<MatrixVectorNature MatrixVectorNatureT>
void GradOnGradientBasedHyperelasticityTensor<LawPolicyT>
::ComputeEltArray(double matrix_coefficient,
double vector_coefficient,
const std::vector<double>& local_displacement)
::ComputeEltArray(const std::vector<double>& local_displacement)
{
auto& elementary_data = GetNonCstElementaryData();
......@@ -189,9 +187,7 @@ namespace HappyHeart
// Complete the tangent matrix with the linear part!
Seldon::Add(1., linear_part, tangent_matrix);
ComputeAtQuadraturePoint<MatrixVectorNatureT>(matrix_coefficient,
vector_coefficient,
infos_at_quad_pt,
ComputeAtQuadraturePoint<MatrixVectorNatureT>(infos_at_quad_pt,
tangent_matrix,
rhs_part,
elementary_data);
......@@ -204,9 +200,7 @@ namespace HappyHeart
template <class LawPolicyT>
template<MatrixVectorNature MatrixVectorNatureT>
void GradOnGradientBasedHyperelasticityTensor<LawPolicyT>
::ComputeAtQuadraturePoint(const double matrix_coefficient,
const double vector_coefficient,
const Private::ShapeFunctionAtQuadraturePoint& infos_at_quad_pt,
::ComputeAtQuadraturePoint(const Private::ShapeFunctionAtQuadraturePoint& infos_at_quad_pt,
const LocalMatrix& tangent_matrix,
const LocalVector& rhs_part,
ElementaryData<MatrixVectorNature::vector_and_matrix>& elementary_data) const
......@@ -223,9 +217,6 @@ namespace HappyHeart
const auto weight_meas = infos_at_quad_pt.GetQuadraturePoint().GetWeight()
* infos_at_quad_pt.GetJacobianDeterminant();
const auto weight_meas_with_matrix_coefficient = weight_meas * matrix_coefficient;
const auto weight_meas_with_vector_coefficient = weight_meas * vector_coefficient;
const auto& ref_felt = elementary_data.GetRefFElt(this->GetUnknown());
......@@ -261,7 +252,7 @@ namespace HappyHeart
gradient_based_block,
dPhi_mult_gradient_based_block);
Seldon::Mlt(weight_meas_with_matrix_coefficient,
Seldon::Mlt(weight_meas,
dPhi_mult_gradient_based_block,
transposed_dphi,
block_contribution);
......@@ -296,7 +287,7 @@ namespace HappyHeart
}
// Multiply trans_dPhi by this local rhs. This will yield the new contribution to vector_.
Seldon::Mlt(weight_meas_with_vector_coefficient, dPhi, local_rhs, new_contribution);
Seldon::Mlt(weight_meas, dPhi, local_rhs, new_contribution);
// Finally add this contribution to vector_.
Seldon::Add(1., new_contribution, block_vector);
......
......@@ -129,6 +129,9 @@ namespace HappyHeart
bool is_first_dynamic_call_done,
GlobalVector& irrelevant_argument,
VectorsAndMatrices<HyperelasticityNS::TimeScheme::half_sum>& vm,
GlobalMatrix& foo,
GlobalMatrix& bar,
GlobalVector& baz,
const Domain& domain = Domain());
......
......@@ -30,19 +30,29 @@ namespace HappyHeart
bool is_first_dynamic_call,
GlobalVector& into_vector,
VectorsAndMatrices<HyperelasticityNS::TimeScheme::half_sum>& vm,
GlobalMatrix& foo,
GlobalMatrix& bar,
GlobalVector& baz,
const Domain& domain)
{
static_cast<void>(into_vector);
the_operator.Assemble(0.5, vm.GetNonCstMatrixNewStiffness(),
1., vm.GetNonCstVectorNewStiffness(),
the_operator.Assemble(std::make_tuple(GlobalVectorWithCoefficient(baz, -2111.),
GlobalMatrixWithCoefficient(vm.GetNonCstMatrixNewStiffness(), 0.5),
GlobalMatrixWithCoefficient(bar, 10000.),
GlobalVectorWithCoefficient(vm.GetNonCstVectorNewStiffness(), 1.),
GlobalMatrixWithCoefficient(foo, 10000.)
),
vm.GetEvaluationState(),
domain);
if (is_first_dynamic_call)
the_operator.Assemble(0.5, vm.GetNonCstMatrixCurrentStiffness(),
1., vm.GetNonCstVectorCurrentStiffness(),
the_operator.Assemble(std::make_tuple(GlobalMatrixWithCoefficient(vm.GetNonCstMatrixCurrentStiffness(), 0.5),
GlobalVectorWithCoefficient(vm.GetNonCstVectorCurrentStiffness(), 1.)),
vm.GetCurrentDisplacement(),
domain);
}
......
......@@ -219,7 +219,7 @@ namespace HappyHeart
GlobalMatrix::unique_ptr foo = nullptr; // \todo #531 Remove this once ticket is done!
GlobalMatrix::unique_ptr bar = nullptr;
GlobalVector::unique_ptr baz = nullptr;
///@}
//! Set the volumic mass.
......
......@@ -144,6 +144,7 @@ namespace HappyHeart
foo = std::make_unique<GlobalMatrix>(system_matrix);
bar = std::make_unique<GlobalMatrix>(system_matrix);
baz = std::make_unique<GlobalVector>(system_solution);
{
std::lock_guard<std::mutex> lock(this->GetMutex());
......@@ -275,7 +276,7 @@ namespace HappyHeart
Private::AssembleStiffness(this->GetNonCstStiffnessOperator(),
is_first_dynamic_call,
this->GetNonCstSystemRhs(this->GetNumberingSubset()),
vm);
vm, *foo, *bar, *baz);
}
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment