Skip to content

Commit

Permalink
[GeoMechanicsApplication] Clean up strategies (#12792)
Browse files Browse the repository at this point in the history 
Changes include:
- Removed a redundant overridden member.
- Removed unused function parameters (all related to linear solvers) from the constructors of the strategies. The base class of all our strategies (`GeoMechanicsNewtonRaphsonStrategy`) used to receive a pointer to a linear solver, but that parameter was never used. Consequently, there is no point in supplying one. In fact, requesting one but not using it is very confusing. The constructors of the derived classes simply forwarded the linear solver pointers to the base class constructor. Also there, the corresponding parameters have been removed from the constructors.
- Moved a few member functions and data members that were defined in class `GeoMechanicsNewtonRaphsonStrategy` to one of the derived classes (`GeoMechanicsRammArcLengthStrategy`), since only there they were actually being used. As a result, the moved members could be made `private` rather than `protected`.
- Replaced a local implementation of the L2 norm calculation by calling Boost's `norm_2` function. There is one thing we need to keep in mind here: the local implementation used to be done in parallel (although it was using deprecated functionality to achieve that). In other words, there might be a performance impact here. The changes made also fix three code smells found by SonarQube.
  • Loading branch information
@avdg81
avdg81 authored Oct 28, 2024
1 parent 98d2dbf commit 3ecd8ee
Show file tree
Hide file tree
Showing 7 changed files with 50 additions and 105 deletions.
6 changes: 3 additions & 3 deletions applications/GeoMechanicsApplication/custom_python/add_custom_strategies_to_python.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -110,18 +110,18 @@ void AddCustomStrategiesToPython(pybind11::module& m)

py::class_<GeoMechanicsNewtonRaphsonStrategyType, typename GeoMechanicsNewtonRaphsonStrategyType::Pointer, BaseSolvingStrategyType>(
m, "GeoMechanicsNewtonRaphsonStrategy")
.def(py::init<ModelPart&, BaseSchemeType::Pointer, LinearSolverType::Pointer, ConvergenceCriteriaType::Pointer,
.def(py::init<ModelPart&, BaseSchemeType::Pointer, ConvergenceCriteriaType::Pointer,
BuilderAndSolverType::Pointer, Parameters&, int, bool, bool, bool>());

py::class_<GeoMechanicsNewtonRaphsonErosionProcessStrategyType,
typename GeoMechanicsNewtonRaphsonErosionProcessStrategyType::Pointer, BaseSolvingStrategyType>(
m, "GeoMechanicsNewtonRaphsonErosionProcessStrategy")
.def(py::init<ModelPart&, BaseSchemeType::Pointer, LinearSolverType::Pointer, ConvergenceCriteriaType::Pointer,
.def(py::init<ModelPart&, BaseSchemeType::Pointer, ConvergenceCriteriaType::Pointer,
BuilderAndSolverType::Pointer, Parameters&, int, bool, bool, bool>());

py::class_<GeoMechanicsRammArcLengthStrategyType, typename GeoMechanicsRammArcLengthStrategyType::Pointer, BaseSolvingStrategyType>(
m, "GeoMechanicsRammArcLengthStrategy")
.def(py::init<ModelPart&, BaseSchemeType::Pointer, LinearSolverType::Pointer, ConvergenceCriteriaType::Pointer,
.def(py::init<ModelPart&, BaseSchemeType::Pointer, ConvergenceCriteriaType::Pointer,
BuilderAndSolverType::Pointer, Parameters&, int, bool, bool, bool>())
.def("UpdateLoads", &GeoMechanicsRammArcLengthStrategyType::UpdateLoads);

Expand Down
12 changes: 1 addition & 11 deletions ...on/custom_strategies/strategies/geo_mechanics_newton_raphson_erosion_process_strategy.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -54,7 +54,6 @@ class GeoMechanicsNewtonRaphsonErosionProcessStrategy

GeoMechanicsNewtonRaphsonErosionProcessStrategy(ModelPart& model_part,
typename TSchemeType::Pointer pScheme,
typename TLinearSolver::Pointer pNewLinearSolver,
typename TConvergenceCriteriaType::Pointer pNewConvergenceCriteria,
typename TBuilderAndSolverType::Pointer pNewBuilderAndSolver,
Parameters& rParameters,
Expand All @@ -63,16 +62,7 @@ class GeoMechanicsNewtonRaphsonErosionProcessStrategy
bool ReformDofSetAtEachStep = false,
bool MoveMeshFlag = false)
: GeoMechanicsNewtonRaphsonStrategy<TSparseSpace, TDenseSpace, TLinearSolver>(
model_part,
pScheme,
pNewLinearSolver,
pNewConvergenceCriteria,
pNewBuilderAndSolver,
rParameters,
MaxIterations,
CalculateReactions,
ReformDofSetAtEachStep,
MoveMeshFlag)
model_part, pScheme, pNewConvergenceCriteria, pNewBuilderAndSolver, rParameters, MaxIterations, CalculateReactions, ReformDofSetAtEachStep, MoveMeshFlag)
{
rank = model_part.GetCommunicator().MyPID();
mPipingIterations = rParameters["max_piping_iterations"].GetInt();
Expand Down
72 changes: 3 additions & 69 deletions ...chanicsApplication/custom_strategies/strategies/geo_mechanics_newton_raphson_strategy.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -52,16 +52,14 @@ class GeoMechanicsNewtonRaphsonStrategy
* @brief Default constructor
* @param rModelPart The model part of the problem
* @param pScheme The integration scheme
* @param pNewLinearSolver The linear solver employed
* @param pNewConvergenceCriteria The convergence criteria employed
* @param MaxIterations The maximum number of iterations
* @param CalculateReactions The flag for the reaction calculation
* @param ReformDofSetAtEachStep The flag that allows to compute the modification of the DOF
* @param MoveMeshFlag The flag that allows to move the mesh
*/
GeoMechanicsNewtonRaphsonStrategy(ModelPart& model_part,
typename TSchemeType::Pointer pScheme,
typename TLinearSolver::Pointer pNewLinearSolver,
GeoMechanicsNewtonRaphsonStrategy(ModelPart& model_part,
typename TSchemeType::Pointer pScheme,
typename TConvergenceCriteriaType::Pointer pNewConvergenceCriteria,
typename TBuilderAndSolverType::Pointer pNewBuilderAndSolver,
Parameters& rParameters,
Expand All @@ -70,15 +68,7 @@ class GeoMechanicsNewtonRaphsonStrategy
bool ReformDofSetAtEachStep = false,
bool MoveMeshFlag = false)
: ResidualBasedNewtonRaphsonStrategy<TSparseSpace, TDenseSpace, TLinearSolver>(
model_part,
pScheme,
/*pNewLinearSolver,*/
pNewConvergenceCriteria,
pNewBuilderAndSolver,
MaxIterations,
CalculateReactions,
ReformDofSetAtEachStep,
MoveMeshFlag)
model_part, pScheme, pNewConvergenceCriteria, pNewBuilderAndSolver, MaxIterations, CalculateReactions, ReformDofSetAtEachStep, MoveMeshFlag)
{
// only include validation with c++11 since raw_literals do not exist in c++03
Parameters default_parameters(R"(
Expand All @@ -100,62 +90,6 @@ class GeoMechanicsNewtonRaphsonStrategy

// Validate against defaults -- this also ensures no type mismatch
rParameters.ValidateAndAssignDefaults(default_parameters);

// Set Load SubModelParts and Variable names
if (rParameters["loads_sub_model_part_list"].size() > 0) {
mSubModelPartList.resize(rParameters["loads_sub_model_part_list"].size());
mVariableNames.resize(rParameters["loads_variable_list"].size());

if (mSubModelPartList.size() != mVariableNames.size())
KRATOS_ERROR << "For each SubModelPart there must be a corresponding nodal Variable"
<< std::endl;

for (unsigned int i = 0; i < mVariableNames.size(); ++i) {
mSubModelPartList[i] =
&(model_part.GetSubModelPart(rParameters["loads_sub_model_part_list"][i].GetString()));
mVariableNames[i] = rParameters["loads_variable_list"][i].GetString();
}
}
}

protected:
/// Member Variables
std::vector<ModelPart*> mSubModelPartList; /// List of every SubModelPart associated to an external load
std::vector<std::string> mVariableNames; /// Name of the nodal variable associated to every SubModelPart

int Check() override
{
KRATOS_TRY

return MotherType::Check();

KRATOS_CATCH("")
}

double CalculateReferenceDofsNorm(DofsArrayType& rDofSet)
{
double ReferenceDofsNorm = 0.0;

int NumThreads = ParallelUtilities::GetNumThreads();
OpenMPUtils::PartitionVector DofSetPartition;
OpenMPUtils::DivideInPartitions(rDofSet.size(), NumThreads, DofSetPartition);

#pragma omp parallel reduction(+ : ReferenceDofsNorm)
{
int k = OpenMPUtils::ThisThread();

typename DofsArrayType::iterator DofsBegin = rDofSet.begin() + DofSetPartition[k];
typename DofsArrayType::iterator DofsEnd = rDofSet.begin() + DofSetPartition[k + 1];

for (typename DofsArrayType::iterator itDof = DofsBegin; itDof != DofsEnd; ++itDof) {
if (itDof->IsFree()) {
const double& temp = itDof->GetSolutionStepValue();
ReferenceDofsNorm += temp * temp;
}
}
}

return sqrt(ReferenceDofsNorm);
}
}; // Class GeoMechanicsNewtonRaphsonStrategy

Expand Down
56 changes: 40 additions & 16 deletions ...hanicsApplication/custom_strategies/strategies/geo_mechanics_ramm_arc_length_strategy.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,9 @@
// Application includes
#include "geo_mechanics_application_variables.h"

#include <boost/numeric/ublas/vector_expression.hpp>
#include <boost/range/adaptor/filtered.hpp>

namespace Kratos
{

Expand Down Expand Up @@ -50,12 +53,9 @@ class GeoMechanicsRammArcLengthStrategy
using GrandMotherType::mpDx; // Delta x of iteration i
using GrandMotherType::mpScheme;
using GrandMotherType::mReformDofSetAtEachStep;
using MotherType::mSubModelPartList;
using MotherType::mVariableNames;

GeoMechanicsRammArcLengthStrategy(ModelPart& model_part,
typename TSchemeType::Pointer pScheme,
typename TLinearSolver::Pointer pNewLinearSolver,
GeoMechanicsRammArcLengthStrategy(ModelPart& model_part,
typename TSchemeType::Pointer pScheme,
typename TConvergenceCriteriaType::Pointer pNewConvergenceCriteria,
typename TBuilderAndSolverType::Pointer pNewBuilderAndSolver,
Parameters& rParameters,
Expand All @@ -64,17 +64,24 @@ class GeoMechanicsRammArcLengthStrategy
bool ReformDofSetAtEachStep = false,
bool MoveMeshFlag = false)
: GeoMechanicsNewtonRaphsonStrategy<TSparseSpace, TDenseSpace, TLinearSolver>(
model_part,
pScheme,
pNewLinearSolver,
pNewConvergenceCriteria,
pNewBuilderAndSolver,
rParameters,
MaxIterations,
CalculateReactions,
ReformDofSetAtEachStep,
MoveMeshFlag)
model_part, pScheme, pNewConvergenceCriteria, pNewBuilderAndSolver, rParameters, MaxIterations, CalculateReactions, ReformDofSetAtEachStep, MoveMeshFlag)
{
// Set Load SubModelParts and Variable names
if (rParameters["loads_sub_model_part_list"].size() > 0) {
mSubModelPartList.resize(rParameters["loads_sub_model_part_list"].size());
mVariableNames.resize(rParameters["loads_variable_list"].size());

if (mSubModelPartList.size() != mVariableNames.size())
KRATOS_ERROR << "For each SubModelPart there must be a corresponding nodal Variable"
<< std::endl;

for (unsigned int i = 0; i < mVariableNames.size(); ++i) {
mSubModelPartList[i] =
&(model_part.GetSubModelPart(rParameters["loads_sub_model_part_list"][i].GetString()));
mVariableNames[i] = rParameters["loads_variable_list"][i].GetString();
}
}

mDesiredIterations = rParameters["desired_iterations"].GetInt();
mMaxRadiusFactor = rParameters["max_radius_factor"].GetDouble();
mMinRadiusFactor = rParameters["min_radius_factor"].GetDouble();
Expand Down Expand Up @@ -306,7 +313,7 @@ class GeoMechanicsRammArcLengthStrategy
else if (mRadius < mMinRadiusFactor * mRadius_0) mRadius = mMinRadiusFactor * mRadius_0;

// Update Norm of x
mNormxEquilibrium = this->CalculateReferenceDofsNorm(rDofSet);
mNormxEquilibrium = CalculateReferenceDofsNorm(rDofSet);
} else {
std::cout << "************ NO CONVERGENCE: restoring equilibrium path ************" << std::endl;

Expand Down Expand Up @@ -534,6 +541,23 @@ class GeoMechanicsRammArcLengthStrategy
// Save the applied Lambda factor
mLambda_old = mLambda;
}

private:
static double CalculateReferenceDofsNorm(DofsArrayType& rDofSet)
{
auto is_free_dof = [](const auto& rDof) { return rDof.IsFree(); };
auto free_dofs = rDofSet | boost::adaptors::filtered(is_free_dof);

auto free_dof_values = Vector{boost::size(free_dofs)};
auto get_dof_value = [](const auto& rDof) { return rDof.GetSolutionStepValue(); };
std::transform(std::begin(free_dofs), std::end(free_dofs), free_dof_values.begin(), get_dof_value);

return norm_2(free_dof_values);
}

std::vector<ModelPart*> mSubModelPartList; // List of every SubModelPart associated to an external load
std::vector<std::string> mVariableNames; // Name of the nodal variable associated to every SubModelPart

}; // Class GeoMechanicsRammArcLengthStrategy

} // namespace Kratos
2 changes: 1 addition & 1 deletion applications/GeoMechanicsApplication/custom_utilities/solving_strategy_factory.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -69,7 +69,7 @@ class SolvingStrategyFactory
ParametersUtilities::CopyOptionalParameters(rSolverSettings, strategy_entries);
auto result =
std::make_unique<GeoMechanicsNewtonRaphsonStrategy<TSparseSpace, TDenseSpace, TLinearSolver>>(
rModelPart, scheme, solver, criteria, builder_and_solver, strategy_parameters,
rModelPart, scheme, criteria, builder_and_solver, strategy_parameters,
max_iterations, calculate_reactions, reform_dof_set_at_each_step, move_mesh_flag);
result->SetEchoLevel(echo_level);
return result;
Expand Down
4 changes: 2 additions & 2 deletions applications/GeoMechanicsApplication/custom_workflows/dgeoflow.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -138,8 +138,8 @@ KratosExecute::GeoMechanicsNewtonRaphsonErosionProcessStrategyType::Pointer Krat

auto pSolvingStrategy =
Kratos::make_unique<GeoMechanicsNewtonRaphsonErosionProcessStrategy<SparseSpaceType, LocalSpaceType, KratosExecute::LinearSolverType>>(
rModelPart, p_scheme, p_solver, p_criteria, p_builder_and_solver, p_parameters,
MaxIterations, CalculateReactions, ReformDofSetAtEachStep, MoveMeshFlag);
rModelPart, p_scheme, p_criteria, p_builder_and_solver, p_parameters, MaxIterations,
CalculateReactions, ReformDofSetAtEachStep, MoveMeshFlag);

pSolvingStrategy->Check();
return pSolvingStrategy;
Expand Down
3 changes: 0 additions & 3 deletions applications/GeoMechanicsApplication/python_scripts/geomechanics_solver.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -489,7 +489,6 @@ def _ConstructSolver(self, builder_and_solver, strategy_type):
self.strategy_params.AddValue("loads_variable_list",self.settings["loads_variable_list"])
solving_strategy = GeoMechanicsApplication.GeoMechanicsNewtonRaphsonStrategy(self.computing_model_part,
self.scheme,
self.linear_solver,
self.convergence_criterion,
builder_and_solver,
self.strategy_params,
Expand All @@ -504,7 +503,6 @@ def _ConstructSolver(self, builder_and_solver, strategy_type):
self.strategy_params.AddValue("max_piping_iterations", self.settings["max_piping_iterations"])
solving_strategy = GeoMechanicsApplication.GeoMechanicsNewtonRaphsonErosionProcessStrategy(self.computing_model_part,
self.scheme,
self.linear_solver,
self.convergence_criterion,
builder_and_solver,
self.strategy_params,
Expand Down Expand Up @@ -545,7 +543,6 @@ def _ConstructSolver(self, builder_and_solver, strategy_type):
self.strategy_params.AddValue("loads_variable_list",self.settings["loads_variable_list"])
solving_strategy = GeoMechanicsApplication.GeoMechanicsRammArcLengthStrategy(self.computing_model_part,
self.scheme,
self.linear_solver,
self.convergence_criterion,
builder_and_solver,
self.strategy_params,
Expand Down

0 comments on commit 3ecd8ee

Please sign in to comment.