Skip to content
GitLab
Commit f6a7c7aa authored by GILLES Sebastien's avatar GILLES Sebastien
Browse files

#718 FSI: clean-up (actually there were no bugs related to half sum: it was... 

#718 FSI: clean-up (actually there were no bugs related to half sum: it was Freefems cript which was at fault).
parent 882f20f7
Hide whitespace changes
Inline Side-by-side
Data/Lua/demo_input_fsi_ei.lua
View file @ f6a7c7aa
......@@ -21,7 +21,7 @@ transient = {
-- Maximum time, if set to zero run a case.
-- Expected format: VALUE
-- Constraint: v >= 0.
timeMax = 1.5e-4
timeMax = 5.5e-4
} -- transient
......
Sources/ModelInstances/FSI_EI/SolidVariationalFormulationPolicy/Hyperelasticity/VariationalFormulation.hpp
View file @ f6a7c7aa
......@@ -108,14 +108,11 @@ namespace HappyHeart
* \param[in] mpi Object that manages mpi informations.
* \param[in] time_manager Object in charge of time management.
* \param[in] god_of_dof God of dof into which the formulation works.
* \param[in] main_felt_space Finite element space upon which the variational formulation apply.
* \param[in] neumann_felt_space Finite element space upon which the Neumann condition apply.
* \param[in] felt_space Finite element space upon which the variational formulation apply.
* \param[in] solid_displacement Unknown considered; it should be a solid displacement.
* \param[in] dof_source_policy_args Supplementary arguments that might be required by the DofSourcePolicyT.
*
*/
explicit VariationalFormulationHyperElasticity(const Wrappers::Mpi& mpi,
const FEltSpace& main_felt_space,
const FEltSpace& felt_space,
const Unknown& solid_displacement,
const NumberingSubset& numbering_subset,
const time_manager_type& time_manager,
......@@ -123,7 +120,6 @@ namespace HappyHeart
DirichletBoundaryCondition::vector_shared_ptr&& boundary_condition_list,
double factor,
const GlobalVector& dof_source);
//! Destructor.
~VariationalFormulationHyperElasticity() = default;
......
Sources/ModelInstances/FSI_EI/SolidVariationalFormulationPolicy/Hyperelasticity/VariationalFormulation.hxx
View file @ f6a7c7aa
......@@ -257,15 +257,17 @@ namespace HappyHeart
inline void VariationalFormulationHyperElasticity<LawPolicyT, TimeSchemeT>
::ComputeResidual()
{
std::cout << "COMPUTE RESIDUAL"<< std::endl;
const auto& numbering_subset = GetNumberingSubset();
this->GetNonCstSystemRhs(numbering_subset).ZeroEntries(__FILE__, __LINE__);
auto& rhs = this->GetNonCstSystemRhs(numbering_subset);
auto& vectors_and_matrices = this->GetNonCstVectorsAndMatrices();
auto& displacement_previous_time_iteration = vectors_and_matrices.GetNonCstDisplacementPreviousTimeIteration();
auto& velocity_previous_time_iteration = vectors_and_matrices.GetNonCstVelocityPreviousTimeIteration();
auto& rhs = this->GetNonCstSystemRhs(numbering_subset);
rhs.ZeroEntries(__FILE__, __LINE__);
const auto& displacement_previous_time_iteration =
vectors_and_matrices.GetDisplacementPreviousTimeIteration();
const auto& velocity_previous_time_iteration =
vectors_and_matrices.GetVelocityPreviousTimeIteration();
std::lock_guard<std::mutex> lock(this->GetMutex());
......@@ -273,7 +275,8 @@ namespace HappyHeart
helper_vector.Copy(vectors_and_matrices.GetEvaluationState(), __FILE__, __LINE__);
Wrappers::Petsc::AXPY(-1., displacement_previous_time_iteration, helper_vector, __FILE__, __LINE__);
Wrappers::Petsc::AXPY(- this->GetTimeManager().GetTimeStep(), velocity_previous_time_iteration, helper_vector,
Wrappers::Petsc::AXPY(- this->GetTimeManager().GetTimeStep(), velocity_previous_time_iteration,
helper_vector,
__FILE__, __LINE__);
auto& dynamic_contribution_to_rhs = this->GetNonCstHelperGlobalVector<1>();
......@@ -287,10 +290,6 @@ namespace HappyHeart
dof_source_parent::AddToRhs(rhs);
// vectors_and_matrices.GetDisplacementPreviousTimeIteration().View(this->MpiHappyHeart(), __FILE__, __LINE__);
this->template ApplyEssentialBoundaryCondition<OperatorNS::Nature::linear>(numbering_subset,
numbering_subset);
}
......@@ -304,57 +303,22 @@ namespace HappyHeart
inline void VariationalFormulationHyperElasticity<LawPolicyT, TimeSchemeT>
::ComputeTangent()
{
std::cout << "COMPUTE TANGENT"<< std::endl;
const auto& numbering_subset = GetNumberingSubset();
this->GetNonCstSystemMatrix(numbering_subset, numbering_subset).ZeroEntries(__FILE__, __LINE__);
auto& system_matrix = this->GetNonCstSystemMatrix(numbering_subset, numbering_subset);
const auto& vectors_and_matrices = this->GetVectorsAndMatrices();
auto& system_matrix = this->GetNonCstSystemMatrix(numbering_subset, numbering_subset);
// Now all elements are settled to build the system.
{
system_matrix.Copy(vectors_and_matrices.GetMassPerSquareTime(), __FILE__, __LINE__);
static bool is_first = true;
if (is_first)
{
is_first = false;
vectors_and_matrices.GetMassPerSquareTime().View(this->MpiHappyHeart(), "/Users/sebastien/Desktop/mass.m",
__FILE__, __LINE__,
PETSC_VIEWER_ASCII_MATLAB);
vectors_and_matrices.GetMatrixNewStiffness().View(this->MpiHappyHeart(), "/Users/sebastien/Desktop/new_stiffness.m",
__FILE__, __LINE__,
PETSC_VIEWER_ASCII_MATLAB);
// vectors_and_matrices.GetMatrixStiffnessPreviousTimeIteration().View(this->MpiHappyHeart(), "/Users/sebastien/Desktop/StiffnessPreviousTimeIteration.m",
// __FILE__, __LINE__,
// PETSC_VIEWER_ASCII_MATLAB);
}
HyperelasticityNS::Private::ComputeDynamicTimeSchemeTangentContribution(vectors_and_matrices, system_matrix);
HyperelasticityNS::Private::ComputeDynamicTimeSchemeTangentContribution(vectors_and_matrices,
system_matrix);
}
this->template ApplyEssentialBoundaryCondition<OperatorNS::Nature::bilinear>(numbering_subset,
numbering_subset);
static bool first = true;
if (first)
{
std::cout << "IMPLICIT SOLID MATRIX = " << std::endl;
this->GetSystemMatrix(numbering_subset, numbering_subset).View(this->MpiHappyHeart(), "/Users/sebastien/Desktop/ImplicitSolid.m",
__FILE__, __LINE__,
PETSC_VIEWER_ASCII_MATLAB
);
first = false;
}
}
......@@ -367,7 +331,6 @@ namespace HappyHeart
void VariationalFormulationHyperElasticity<LawPolicyT, TimeSchemeT>
::UpdateVectorsAndMatrices(const Vec& a_evaluation_state)
{
std::cout << "UPDATE VM" << std::endl;
auto& vm = this->GetNonCstVectorsAndMatrices();
vm.UpdateEvaluationState(a_evaluation_state);
......
Sources/ModelInstances/Hyperelasticity/TimeSchemes/HalfSum/HalfSum.cpp
View file @ f6a7c7aa
......@@ -75,7 +75,7 @@ namespace HappyHeart
vm.GetVectorNewStiffness(),
system_rhs,
__FILE__, __LINE__);
Wrappers::Petsc::AXPY(1.,
vm.GetVectorStiffnessPreviousTimeIteration(),
system_rhs,
......@@ -83,7 +83,6 @@ namespace HappyHeart
}
void TimeSchemeDependantUpdateForNextTimeStep(VectorsAndMatrices<HyperelasticityNS::TimeScheme::half_sum>& vm)
{
// Swap current stiffness and new stiffness, because:
......
Sources/ModelInstances/Hyperelasticity/TimeSchemes/HalfSum/HalfSum.hxx
View file @ f6a7c7aa
......@@ -31,8 +31,6 @@ namespace HappyHeart
VectorsAndMatrices<HyperelasticityNS::TimeScheme::half_sum>& vm,
const Domain& domain)
{
std::cout << "ASSEMBLE STIFFNESS" << std::endl;
auto& matrix_stiffness = vm.GetNonCstMatrixNewStiffness();
auto& vector_stiffness = vm.GetNonCstVectorNewStiffness();
......@@ -62,6 +60,8 @@ namespace HappyHeart
GlobalMatrixWithCoefficient mat(matrix_stiffness_previous_time_iteration, 0.5);
GlobalVectorWithCoefficient vec(vector_stiffness_previous_time_iteration, 1.);
vm.GetNonCstDisplacementPreviousTimeIteration().ZeroEntries(__FILE__, __LINE__); // a safety: should already be zero!
the_operator.Assemble(std::make_tuple(std::ref(mat), std::ref(vec)),
vm.GetDisplacementPreviousTimeIteration(),
domain);
......
Sources/ModelInstances/Hyperelasticity/TimeSchemes/Midpoint/Midpoint.hxx
View file @ f6a7c7aa
......@@ -37,7 +37,6 @@ namespace HappyHeart
matrix_stiffness.ZeroEntries(__FILE__, __LINE__);
vector_stiffness.ZeroEntries(__FILE__, __LINE__);
GlobalMatrixWithCoefficient mat(matrix_stiffness, 1.);
GlobalVectorWithCoefficient vec(vector_stiffness, 1.);
......
Sources/ModelInstances/Hyperelasticity/VariationalFormulationHyperElasticity.hxx
View file @ f6a7c7aa
......@@ -396,7 +396,7 @@ namespace HappyHeart
inline void VariationalFormulationHyperElasticity<LawPolicyT, TimeSchemeT, VolumicIndexT, SurfacicIndexT, DofSourcePolicyT>
::ComputeResidual()
{
const auto& numbering_subset = GetNumberingSubset();
const auto& numbering_subset = GetNumberingSubset();
this->GetNonCstSystemRhs(numbering_subset).ZeroEntries(__FILE__, __LINE__);
switch(this->GetTimeManager().GetStaticOrDynamic())
......
Sources/Operators/ConformInterpolatorInstances/P1_to_P2.cpp
View file @ f6a7c7aa
......@@ -43,16 +43,15 @@ namespace HappyHeart
const auto& target_unknown_storage = interpolation_data.GetTargetData().GetExtendedUnknownList();
// \todo #714 Deactivate this to actually authorize P1 -> P1.
// assert(std::all_of(target_unknown_storage.cbegin(),
// target_unknown_storage.cend(),
// [](const auto& extended_unknown_ptr)
// {
// if (!extended_unknown_ptr)
// return true; // some unknowns may be dropped.
//
// return extended_unknown_ptr->GetShapeFunctionLabel() == "P2";
// }));
assert(std::all_of(target_unknown_storage.cbegin(),
target_unknown_storage.cend(),
[](const auto& extended_unknown_ptr)
{
if (!extended_unknown_ptr)
return true; // some unknowns may be dropped.
return extended_unknown_ptr->GetShapeFunctionLabel() == "P2";
}));
#endif // NDEBUG
}
......
Sources/Operators/GlobalVariationalOperatorInstances/NonlinearForm/GradOnGradientBasedHyperelasticityTensor.hxx
View file @ f6a7c7aa
......@@ -76,8 +76,8 @@ namespace HappyHeart
std::tuple<const GlobalVector&>&& additional_arguments) const
{
return std::make_tuple(parent::ExtractFEltValuesFrom(local_operator.GetElementaryData().GetRefFElt(this->GetNthUnknown()),
local_felt_space,
std::get<0>(additional_arguments)));
local_felt_space,
std::get<0>(additional_arguments)));
}
......
Sources/Operators/LocalVariationalOperatorInstances/NonlinearForm/GradOnGradientBasedHyperelasticityTensor.hxx
View file @ f6a7c7aa
......@@ -140,8 +140,6 @@ namespace HappyHeart
{
auto& elementary_data = GetNonCstElementaryData();
// LocalVector pressure; // Not used yet; see #7.
OperatorNS::Nature what_must_be_built = OperatorNS::Nature::nonlinear; // both matrix and vector must be built.
// see #533.
......
Sources/Operators/LocalVariationalOperatorInstances/NonlinearForm/Hyperelasticity/Private/Helper.cpp
View file @ f6a7c7aa
......@@ -92,7 +92,6 @@ namespace HappyHeart
De(2, 2) = 1. + gradient_component_dispy[1];
De(2, 3) = gradient_component_dispy[0];
return CauchyGreenTensor<2>(gradient_component_dispx, gradient_component_dispy);
}
......
Supports Markdown
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