Skip to content

GitLab

Commit 2b7fc640 authored by GILLES Sebastien's avatar GILLES Sebastien
Browse files

#859 Petsc: also provide the newly introduced parameter to MatMatMult and MatMatMatMult.

parent 787f66b0
Hide whitespace changes
Inline Side-by-side
Showing with 36 additions and 65 deletions
Data/Lua/demo_input_poromechanics.lua
View file @ 2b7fc640
......@@ -21,7 +21,7 @@ transient = {
-- Maximum time, if set to zero run a static case.
-- Expected format: VALUE
-- Constraint: v >= 0.
timeMax = 1.e-3
timeMax = 2.e-3
} -- transient
......
Sources/ModelInstances/UnderDevelopment/Poromechanics/ImplicitStepFluidVariationalFormulation.cpp
View file @ 2b7fc640
......@@ -36,6 +36,12 @@ namespace HappyHeart
decltype(auto) numbering_subset =
god_of_dof.GetNumberingSubset(EnumUnderlyingType(NumberingSubsetIndex::fluid_mass_velocity_pressure));
static bool is_first_call = true;
const auto do_reuse_matrix = is_first_call
? Wrappers::Petsc::DoReuseMatrix::no
: Wrappers::Petsc::DoReuseMatrix::yes;
// \todo #820 NOT YET TESTED FOR ONE OF THE CONTRIBUTION; WAIT SECOND ITERATION TO SEE IF IT IS OK.
// auto& mass_matrix = GetNonCstMassMatrix();
......@@ -101,12 +107,14 @@ namespace HappyHeart
Wrappers::Petsc::MatMatMult(fluid_mass_matrix_on_solid_mesh,
inverse_interp_mass_from_solid_mesh_->GetInterpolationMatrix(),
intermediate,
__FILE__, __LINE__);
__FILE__, __LINE__,
do_reuse_matrix);
Wrappers::Petsc::MatTransposeMatMult(inverse_interp_mass_from_solid_mesh_->GetInterpolationMatrix(),
intermediate,
t11,
__FILE__, __LINE__);
__FILE__, __LINE__,
do_reuse_matrix);
// Wrappers::Petsc::PtAP(fluid_mass_matrix_on_solid_mesh,
// GetInterpMassFromSolidMesh().GetInterpolationMatrix(),
......@@ -382,6 +390,8 @@ namespace HappyHeart
// parent::SolveLinear<IsFactorized::no>(numbering_subset, numbering_subset);
// parent::WriteSolution(time_manager, numbering_subset);
is_first_call = false;
}
......
Sources/ModelInstances/UnderDevelopment/Poromechanics/ImplicitStepFluidVariationalFormulationInit.cpp
View file @ 2b7fc640
......@@ -263,61 +263,8 @@ namespace HappyHeart
reduce_to_mass_->Init();
}
{
// auto& interpolator = GetNonCstFluidMassInterpolationMatrix();
//
// Wrappers::Petsc::MatMatMult(GetExpandMass().GetInterpolationMatrix(),
// GetInterpMassFromSolidMesh().GetInterpolationMatrix(),
// interpolator,
// __FILE__, __LINE__);
//
//
// auto& inverse_interpolator = GetNonCstInverseFluidMassInterpolationMatrix();
//
// Wrappers::Petsc::MatMatMult(GetInterpMassFromSolidMesh().GetInterpolationMatrix(),
// GetReduceToMass().GetInterpolationMatrix(),
// inverse_interpolator,
// __FILE__, __LINE__);
//
}
{
// Allocate the pattern of a matrix from the result of the product of two other matrices.
const auto& fluid_mass_matrix_on_solid_mesh = GetFluidMassMatrixOnSolidMesh();
auto& t11 = GetNonCstFluidMassMatrixOnFluidMesh();
auto& intermediate = GetNonCstWorkMatrixFluidMassMatrixOnFluidMesh();
// decltype(auto) interp_matrix = GetInterpMassFromSolidMesh().GetInterpolationMatrix();
//
// decltype(auto) size1 = fluid_mass_matrix_on_solid_mesh.GetProcessorWiseSize(__FILE__, __LINE__);
// decltype(auto) size2 = interp_matrix.GetProcessorWiseSize(__FILE__, __LINE__);
// // decltype(auto) size3 = intermediate.GetProcessorWiseSize(__FILE__, __LINE__);
//
// Utilities::PrintTuple(size1, std::cout, ", ", MpiHappyHeart().GetRankPreffix() + " Matrix1 -> [");
// Utilities::PrintTuple(size2, std::cout, ", ", MpiHappyHeart().GetRankPreffix() + " Matrix2 -> [");
// // Utilities::PrintTuple(size3, std::cout, ", ", MpiHappyHeart().GetRankPreffix() + " Matrix3 -> [");
//
// Wrappers::Petsc::MatMatMultSymbolic(interp_matrix,
// fluid_mass_matrix_on_solid_mesh,
// intermediate,
// __FILE__, __LINE__);
//
// Wrappers::Petsc::MatMatMultNumeric(interp_matrix,
// fluid_mass_matrix_on_solid_mesh,
// intermediate,
// __FILE__, __LINE__);
//
// Wrappers::Petsc::MatTransposeMatMultSymbolic(interp_matrix,
// intermediate,
// t11,
// __FILE__, __LINE__);
}
{
// \todo #820 Only for dev purposes (comparison with freefem output).
......
Sources/ThirdParty/Wrappers/Petsc/Matrix/MatrixOperations.hpp
View file @ 2b7fc640
......@@ -162,6 +162,8 @@ namespace HappyHeart
* One can also reuse a pattern for several matrices, but so far I do not need that in HappyHeart so
* I equally use the default setting.
*
* \copydetails doxygen_hide_petsc_do_reuse_matrix_arg
*
*/
template
<
......@@ -177,7 +179,8 @@ namespace HappyHeart
MatMatMult(const MatrixT& matrix1,
const MatrixT& matrix2,
MatrixU& matrix3,
const char* invoking_file, int invoking_line);
const char* invoking_file, int invoking_line,
DoReuseMatrix do_reuse_matrix = DoReuseMatrix::no);
/*!
......@@ -243,7 +246,7 @@ namespace HappyHeart
* If the operation is performed many times with each time the same non zero pattern for the matrices,
* rather use MatMatMatMultSymbolic/MatMatMatMultNumeric to improve efficiency.
*
* \copydetails
* \copydetails doxygen_hide_petsc_do_reuse_matrix_arg
*/
template
<
......@@ -260,7 +263,8 @@ namespace HappyHeart
const MatrixT& matrix2,
const MatrixT& matrix3,
MatrixU& matrix4,
const char* invoking_file, int invoking_line);
const char* invoking_file, int invoking_line,
DoReuseMatrix do_reuse_matrix = DoReuseMatrix::no);
/*!
......
Sources/ThirdParty/Wrappers/Petsc/Matrix/MatrixOperations.hxx
View file @ 2b7fc640
......@@ -65,18 +65,23 @@ namespace HappyHeart
>
MatMatMult(const MatrixT& matrix1,
const MatrixT& matrix2,
MatrixU& matrix3,
const char* invoking_file, int invoking_line)
MatrixU& out,
const char* invoking_file, int invoking_line,
DoReuseMatrix do_reuse_matrix)
{
Mat result;
if (do_reuse_matrix == DoReuseMatrix::yes)
result = out.Internal();
int error_code = ::MatMatMult(matrix1.Internal(),
matrix2.Internal(),
MAT_INITIAL_MATRIX,
PETSC_DEFAULT,
&result);
matrix3.SetFromPetscMat(result);
if (do_reuse_matrix == DoReuseMatrix::no)
out.SetFromPetscMat(result);
if (error_code)
throw ExceptionNS::Exception(error_code, "MatMatMult", invoking_file, invoking_line);
......@@ -155,11 +160,15 @@ namespace HappyHeart
MatMatMatMult(const MatrixT& matrix1,
const MatrixT& matrix2,
const MatrixT& matrix3,
MatrixU& matrix4,
const char* invoking_file, int invoking_line)
MatrixU& out,
const char* invoking_file, int invoking_line,
DoReuseMatrix do_reuse_matrix)
{
Mat result;
if (do_reuse_matrix == DoReuseMatrix::yes)
result = out.Internal();
int error_code = ::MatMatMatMult(matrix1.Internal(),
matrix2.Internal(),
matrix3.Internal(),
......@@ -167,7 +176,8 @@ namespace HappyHeart
PETSC_DEFAULT,
&result);
matrix4.SetFromPetscMat(result);
if (do_reuse_matrix == DoReuseMatrix::no)
out.SetFromPetscMat(result);
if (error_code)
throw ExceptionNS::Exception(error_code, "MatMatMatMult", invoking_file, invoking_line);
......
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