From 52c2655bcb769bebe1212cf7d135560bc57e21de Mon Sep 17 00:00:00 2001
From: Ludovic Aubry <ludovic.aubry@cea.fr>
Date: Mon, 9 Sep 2024 13:20:10 +0200
Subject: [PATCH] Fix a case when mfront would raise an exception and leave
PetSc vectors in an unwanted state
---
MathLib/LinAlg/LinAlg.h | 16 +++++++++++++++
NumLib/ODESolver/NonlinearSolver.cpp | 6 +++++-
NumLib/ODESolver/TimeDiscretizedODESystem.cpp | 17 +++++++++++-----
.../SmallDeformationProcess.cpp | 20 ++++++++++++++-----
4 files changed, 48 insertions(+), 11 deletions(-)
diff --git a/MathLib/LinAlg/LinAlg.h b/MathLib/LinAlg/LinAlg.h
index 206395e8f90..038dd762743 100644
--- a/MathLib/LinAlg/LinAlg.h
+++ b/MathLib/LinAlg/LinAlg.h
@@ -17,6 +17,8 @@
#ifdef USE_PETSC
#include <petscsystypes.h>
+
+#include "PETSc/PETScVector.h"
#endif
namespace MathLib
@@ -204,6 +206,14 @@ void linearSysNormalize(PETScMatrix const& A, PETScMatrix& new_A,
void finalizeAssembly(PETScMatrix& A);
void finalizeAssembly(PETScVector& x);
+// Reduce operations for interprocess communications while using Petsc
+static inline int reduceMin(int val)
+{
+ int result;
+ MPI_Allreduce(&val, &result, 1, MPI_INTEGER, MPI_MIN, PETSC_COMM_WORLD);
+ return result;
+}
+
} // namespace LinAlg
} // namespace MathLib
@@ -276,6 +286,12 @@ void linearSysNormalize(EigenMatrix const& A, EigenMatrix& new_A,
void finalizeAssembly(EigenMatrix& x);
void finalizeAssembly(EigenVector& A);
+// Reduce operations for interprocess communications while using Petsc
+static inline int reduceMin(int val)
+{
+ return val;
+}
+
} // namespace LinAlg
} // namespace MathLib
diff --git a/NumLib/ODESolver/NonlinearSolver.cpp b/NumLib/ODESolver/NonlinearSolver.cpp
index 59c0886b55b..520e18442d2 100644
--- a/NumLib/ODESolver/NonlinearSolver.cpp
+++ b/NumLib/ODESolver/NonlinearSolver.cpp
@@ -356,6 +356,7 @@ NonlinearSolverStatus NonlinearSolver<NonlinearSolverTag::Newton>::solve(
BaseLib::RunTime time_assembly;
time_assembly.start();
+ int ok = 1, g_ok = 0;
try
{
sys.assemble(x, x_prev, process_id);
@@ -366,8 +367,11 @@ NonlinearSolverStatus NonlinearSolver<NonlinearSolverTag::Newton>::solve(
e.what());
error_norms_met = false;
iteration = _maxiter;
- break;
+ ok = 0;
}
+ g_ok = MathLib::LinAlg::reduceMin(ok);
+ if (g_ok == 0)
+ break;
sys.getResidual(*x[process_id], *x_prev[process_id], res);
sys.getJacobian(J);
INFO("[time] Assembly took {:g} s.", time_assembly.elapsed());
diff --git a/NumLib/ODESolver/TimeDiscretizedODESystem.cpp b/NumLib/ODESolver/TimeDiscretizedODESystem.cpp
index a536e2ff034..cfc584031a5 100644
--- a/NumLib/ODESolver/TimeDiscretizedODESystem.cpp
+++ b/NumLib/ODESolver/TimeDiscretizedODESystem.cpp
@@ -79,11 +79,18 @@ void TimeDiscretizedODESystem<ODESystemTag::FirstOrderImplicitQuasilinear,
_b->setZero();
_Jac->setZero();
-
- _ode.preAssemble(t, dt, x_curr);
- _ode.assembleWithJacobian(t, dt, x_new_timestep, x_prev, process_id, *_b,
- *_Jac);
-
+ try
+ {
+ _ode.preAssemble(t, dt, x_curr);
+ _ode.assembleWithJacobian(t, dt, x_new_timestep, x_prev, process_id,
+ *_b, *_Jac);
+ }
+ catch (AssemblyException const&)
+ {
+ MathLib::LinAlg::finalizeAssembly(*_b);
+ MathLib::LinAlg::finalizeAssembly(*_Jac);
+ throw;
+ }
MathLib::LinAlg::finalizeAssembly(*_b);
MathLib::LinAlg::finalizeAssembly(*_Jac);
}
diff --git a/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp b/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp
index 8e9ea109e51..fdb73c82536 100644
--- a/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp
+++ b/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp
@@ -14,6 +14,7 @@
#include "MeshLib/Utils/IntegrationPointWriter.h"
#include "MeshLib/Utils/getOrCreateMeshProperty.h"
+#include "NumLib/Exceptions.h"
#include "ProcessLib/Deformation/SolidMaterialInternalToSecondaryVariables.h"
#include "ProcessLib/Output/CellAverageAlgorithm.h"
#include "ProcessLib/Process.h"
@@ -162,11 +163,20 @@ void SmallDeformationProcess<DisplacementDim>::
std::vector<NumLib::LocalToGlobalIndexMap const*> dof_table = {
_local_to_global_index_map.get()};
// Call global assembler for each local assembly item.
- GlobalExecutor::executeSelectedMemberDereferenced(
- _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
- _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev,
- process_id, &b, &Jac);
-
+ try
+ {
+ GlobalExecutor::executeSelectedMemberDereferenced(
+ _global_assembler, &VectorMatrixAssembler::assembleWithJacobian,
+ _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x,
+ x_prev, process_id, &b, &Jac);
+ }
+ catch (NumLib::AssemblyException const&)
+ {
+ transformVariableFromGlobalVector(b, 0, *_local_to_global_index_map,
+ *_nodal_forces,
+ std::negate<double>());
+ throw;
+ }
transformVariableFromGlobalVector(b, 0, *_local_to_global_index_map,
*_nodal_forces, std::negate<double>());
--
GitLab