Model v3 - third iteration of the model redesign

applications/DEM_application/custom_strategies/strategies/explicit_solver_strategy.h

@@ -70,11 +70,11 @@ namespace Kratos {
 
         ~ExplicitSolverSettings() {
         }
-        ModelPart::Pointer r_model_part;
-        ModelPart::Pointer contact_model_part;
-        ModelPart::Pointer fem_model_part;
-        ModelPart::Pointer cluster_model_part;
-        ModelPart::Pointer inlet_model_part;
+        ModelPart* r_model_part;
+        ModelPart* contact_model_part;
+        ModelPart* fem_model_part;
+        ModelPart* cluster_model_part;
+        ModelPart* inlet_model_part;
     };
 
     class KRATOS_API(DEM_APPLICATION) ExplicitSolverStrategy {

applications/FluidDynamicsApplication/custom_elements/vms.h

@@ -301,7 +301,8 @@ class VMS : public Element
         this->AddMomentumRHS(rRightHandSideVector, Density, N, Area);
 
         // For OSS: Add projection of residuals to RHS
-        if (rCurrentProcessInfo[OSS_SWITCH] == 1)
+        const ProcessInfo& r_const_process_info = rCurrentProcessInfo;
+        if (r_const_process_info[OSS_SWITCH] == 1)
         {
             array_1d<double, 3 > AdvVel;
             this->GetAdvectiveVel(AdvVel, N);
@@ -353,7 +354,8 @@ class VMS : public Element
          These terms are not used in OSS, as they belong to the finite element
          space and cancel out with their projections.
          */
-        if (rCurrentProcessInfo[OSS_SWITCH] != 1)
+        const ProcessInfo& r_const_process_info = rCurrentProcessInfo;
+        if (r_const_process_info[OSS_SWITCH] != 1)
         {
             double ElemSize = this->ElementSize(Area);
             double Viscosity = this->EffectiveViscosity(Density,N,DN_DX,ElemSize,rCurrentProcessInfo);

applications/FluidDynamicsApplication/custom_processes/boussinesq_force_process.cpp

@@ -4,7 +4,7 @@ namespace Kratos
 {
     /* Public functions *******************************************************/
     BoussinesqForceProcess::BoussinesqForceProcess(
-        ModelPart::Pointer pModelPart,
+        ModelPart* pModelPart,
         Parameters& rParameters):
     mpModelPart(pModelPart),
     mrGravity(array_1d<double,3>(3,0.0))

applications/FluidDynamicsApplication/custom_processes/boussinesq_force_process.h

@@ -83,7 +83,7 @@ namespace Kratos
       ///@{
 
       /// Constructor
-      BoussinesqForceProcess(ModelPart::Pointer pModelPart, Parameters& rParameters);
+      BoussinesqForceProcess(ModelPart* pModelPart, Parameters& rParameters);
 
       /// Destructor.
       ~BoussinesqForceProcess() override;

applications/FluidDynamicsApplication/custom_processes/spalart_allmaras_turbulence_model.h

@@ -27,6 +27,8 @@
 
 // Project includes
 #include "includes/define.h"
+#include "includes/kernel.h"
+#include "containers/model.h"
 #include "processes/process.h"
 #include "includes/cfd_variables.h"
 #include "solving_strategies/strategies/solving_strategy.h"
@@ -103,7 +105,13 @@ class SpalartAllmarasTurbulenceModel : public Process
         unsigned int MaxIter,
         bool ReformDofSet,
         unsigned int TimeOrder)
-        : mr_model_part(rModelPart), mdomain_size(DomainSize), mtol(NonLinearTol), mmax_it(MaxIter), mtime_order(TimeOrder),madapt_for_fractional_step(false)
+        : mr_model_part(rModelPart),
+        mrspalart_model_part(Kernel::GetModel().CreateModelPart("SpalartModelPart")),
+        mdomain_size(DomainSize),
+        mtol(NonLinearTol),
+        mmax_it(MaxIter),
+        mtime_order(TimeOrder),
+        madapt_for_fractional_step(false)
     {
         //************************************************************************************************
         //check that the variables needed are in the model part
@@ -129,11 +137,11 @@ class SpalartAllmarasTurbulenceModel : public Process
 
         //************************************************************************************************
         //construct a new auxiliary model part
-        mspalart_model_part.GetNodalSolutionStepVariablesList() = mr_model_part.GetNodalSolutionStepVariablesList();
-        mspalart_model_part.SetBufferSize(3);
-        mspalart_model_part.Nodes() = mr_model_part.Nodes();
-        mspalart_model_part.SetProcessInfo(mr_model_part.pGetProcessInfo());
-        mspalart_model_part.SetProperties(mr_model_part.pProperties());
+        mrspalart_model_part.GetNodalSolutionStepVariablesList() = mr_model_part.GetNodalSolutionStepVariablesList();
+        mrspalart_model_part.SetBufferSize(3);
+        mrspalart_model_part.Nodes() = mr_model_part.Nodes();
+        mrspalart_model_part.SetProcessInfo(mr_model_part.pGetProcessInfo());
+        mrspalart_model_part.SetProperties(mr_model_part.pProperties());
 
         std::string ElementName;
         if (DomainSize == 2)
@@ -148,7 +156,7 @@ class SpalartAllmarasTurbulenceModel : public Process
         {
             Properties::Pointer properties = iii->pGetProperties();
             Element::Pointer p_element = rReferenceElement.Create(iii->Id(), iii->GetGeometry(), properties);
-            mspalart_model_part.Elements().push_back(p_element);
+            mrspalart_model_part.Elements().push_back(p_element);
         }
 
         // pointer types for the solution strategy construcion
@@ -173,7 +181,7 @@ class SpalartAllmarasTurbulenceModel : public Process
         bool CalculateReactions = false;
         bool MoveMesh = false;
 
-        mpSolutionStrategy = StrategyPointerType( new ResidualBasedNewtonRaphsonStrategy<TSparseSpace, TDenseSpace, TLinearSolver>(mspalart_model_part,pScheme,pLinearSolver,pConvCriteria,pBuildAndSolver,MaxIter,CalculateReactions,ReformDofSet,MoveMesh));
+        mpSolutionStrategy = StrategyPointerType( new ResidualBasedNewtonRaphsonStrategy<TSparseSpace, TDenseSpace, TLinearSolver>(mrspalart_model_part,pScheme,pLinearSolver,pConvCriteria,pBuildAndSolver,MaxIter,CalculateReactions,ReformDofSet,MoveMesh));
         mpSolutionStrategy->SetEchoLevel(0);
         mpSolutionStrategy->Check();
     }
@@ -182,6 +190,7 @@ class SpalartAllmarasTurbulenceModel : public Process
 
     ~SpalartAllmarasTurbulenceModel() override
     {
+        Kernel::GetModel().DeleteModelPart("SpalartModelPart");
     }
 
 
@@ -201,22 +210,22 @@ class SpalartAllmarasTurbulenceModel : public Process
 
         if(madapt_for_fractional_step == true)
         {
-            if (!(mspalart_model_part.NodesBegin()->SolutionStepsDataHas(FRACT_VEL)))
+            if (!(mrspalart_model_part.NodesBegin()->SolutionStepsDataHas(FRACT_VEL)))
                 KRATOS_THROW_ERROR(std::logic_error, "Variable is not in the model part:", FRACT_VEL);
 
             #pragma omp parallel for
-            for (int i = 0; i < static_cast<int>(mspalart_model_part.Nodes().size()); i++)
+            for (int i = 0; i < static_cast<int>(mrspalart_model_part.Nodes().size()); i++)
             {
-                ModelPart::NodesContainerType::iterator it = mspalart_model_part.NodesBegin() + i;
+                ModelPart::NodesContainerType::iterator it = mrspalart_model_part.NodesBegin() + i;
                 it->FastGetSolutionStepValue(VELOCITY) = it->FastGetSolutionStepValue(FRACT_VEL);
             }
         }
 
         AuxSolve();
 
         //update viscosity on the nodes
-        for (ModelPart::NodeIterator i = mspalart_model_part.NodesBegin();
-                i != mspalart_model_part.NodesEnd(); ++i)
+        for (ModelPart::NodeIterator i = mrspalart_model_part.NodesBegin();
+                i != mrspalart_model_part.NodesEnd(); ++i)
         {
             double molecular_viscosity = i->FastGetSolutionStepValue(MOLECULAR_VISCOSITY);
             double turbulent_viscosity = i->FastGetSolutionStepValue(TURBULENT_VISCOSITY);
@@ -264,11 +273,11 @@ class SpalartAllmarasTurbulenceModel : public Process
     {
         KRATOS_TRY;
 
-        mspalart_model_part.GetProcessInfo()[C_DES] = CDES;
+        mrspalart_model_part.GetProcessInfo()[C_DES] = CDES;
 /*
         //update viscosity on the nodes
-        for (ModelPart::NodeIterator i = mspalart_model_part.NodesBegin();
-                i != mspalart_model_part.NodesEnd(); ++i)
+        for (ModelPart::NodeIterator i = mrspalart_model_part.NodesBegin();
+                i != mrspalart_model_part.NodesEnd(); ++i)
         {
             double distance = i->FastGetSolutionStepValue(DISTANCE);
             const array_1d<double,3>& xc = i->Coordinates();
@@ -353,7 +362,7 @@ class SpalartAllmarasTurbulenceModel : public Process
     ///@{
 
     ModelPart& mr_model_part;
-    ModelPart mspalart_model_part;
+    ModelPart& mrspalart_model_part;
     unsigned int mdomain_size;
     double mtol;
     unsigned int mmax_it;
@@ -390,7 +399,7 @@ class SpalartAllmarasTurbulenceModel : public Process
         :
         Process(),
         mr_model_part(rModelPart),
-        mspalart_model_part()
+        mrspalart_model_part(Kernel::GetModel().CreateModelPart("SpalartModelPart"))
     {}
 
     ///@}
@@ -422,7 +431,7 @@ class SpalartAllmarasTurbulenceModel : public Process
         KRATOS_TRY
 
         //calculate the BDF coefficients
-        ProcessInfo& rCurrentProcessInfo = mspalart_model_part.GetProcessInfo();
+        ProcessInfo& rCurrentProcessInfo = mrspalart_model_part.GetProcessInfo();
         double Dt = rCurrentProcessInfo[DELTA_TIME];
 
         if (mtime_order == 2)
@@ -509,8 +518,8 @@ class SpalartAllmarasTurbulenceModel : public Process
 
 
 
-        for (ModelPart::NodeIterator i = mspalart_model_part.NodesBegin();
-                i != mspalart_model_part.NodesEnd(); ++i)
+        for (ModelPart::NodeIterator i = mrspalart_model_part.NodesBegin();
+                i != mrspalart_model_part.NodesEnd(); ++i)
         {
             norm += pow(i->FastGetSolutionStepValue(TURBULENT_VISCOSITY), 2);
         }
@@ -526,11 +535,11 @@ class SpalartAllmarasTurbulenceModel : public Process
     {
         KRATOS_TRY;
 
-        ProcessInfo& rCurrentProcessInfo = mspalart_model_part.GetProcessInfo();
+        ProcessInfo& rCurrentProcessInfo = mrspalart_model_part.GetProcessInfo();
 
         //first of all set to zero the nodal variables to be updated nodally
-        for (ModelPart::NodeIterator i = mspalart_model_part.NodesBegin();
-                i != mspalart_model_part.NodesEnd(); ++i)
+        for (ModelPart::NodeIterator i = mrspalart_model_part.NodesBegin();
+                i != mrspalart_model_part.NodesEnd(); ++i)
         {
             (i)->FastGetSolutionStepValue(TEMP_CONV_PROJ) = 0.00;
             (i)->FastGetSolutionStepValue(NODAL_AREA) = 0.00;
@@ -540,20 +549,20 @@ class SpalartAllmarasTurbulenceModel : public Process
         //and the determination of the nodal area
 
 
-        for (ModelPart::ElementIterator i = mspalart_model_part.ElementsBegin();
-                i != mspalart_model_part.ElementsEnd(); ++i)
+        for (ModelPart::ElementIterator i = mrspalart_model_part.ElementsBegin();
+                i != mrspalart_model_part.ElementsEnd(); ++i)
         {
             (i)->InitializeSolutionStep(rCurrentProcessInfo);
         }
 
-        Communicator& rComm = mspalart_model_part.GetCommunicator();
+        Communicator& rComm = mrspalart_model_part.GetCommunicator();
 
         rComm.AssembleCurrentData(NODAL_AREA);
         rComm.AssembleCurrentData(TEMP_CONV_PROJ);
 
         // Obtain nodal projection of the residual
-        for (ModelPart::NodeIterator i = mspalart_model_part.NodesBegin();
-                i != mspalart_model_part.NodesEnd(); ++i)
+        for (ModelPart::NodeIterator i = mrspalart_model_part.NodesBegin();
+                i != mrspalart_model_part.NodesEnd(); ++i)
         {
             const double NodalArea = i->FastGetSolutionStepValue(NODAL_AREA);
             if(NodalArea > 0.0)

applications/FluidDynamicsApplication/custom_processes/stokes_initialization_process.h

@@ -11,6 +11,8 @@
 
 // Project includes
 #include "includes/define.h"
+#include "includes/kernel.h"
+#include "containers/model.h"
 #include "includes/model_part.h"
 #include "processes/process.h"
 #include "solving_strategies/strategies/solving_strategy.h"
@@ -66,24 +68,24 @@ class  StokesInitializationProcess : public Process
     ///@name Life Cycle
     ///@{
 
-    StokesInitializationProcess(const ModelPart::Pointer pModelPart,
+    StokesInitializationProcess(ModelPart& rModelPart,
                                 typename TLinearSolver::Pointer pLinearSolver,
                                 unsigned int DomainSize,
                                 const Variable<int>& PeriodicPairIndicesVar):
         Process(),
-        mpReferenceModelPart(pModelPart),
+        mrReferenceModelPart(rModelPart),
         mpLinearSolver(pLinearSolver),
         mDomainSize(DomainSize)
     {
         KRATOS_TRY;
 
         // Initialize new model part (same nodes, new elements, no conditions)
-        mpStokesModelPart = ModelPart::Pointer(new ModelPart("StokesModelPart"));
-        mpStokesModelPart->GetNodalSolutionStepVariablesList() = mpReferenceModelPart->GetNodalSolutionStepVariablesList();
+        mpStokesModelPart = &(Kernel::GetModel().CreateModelPart("StokesModelPart"));
+        mpStokesModelPart->GetNodalSolutionStepVariablesList() = mrReferenceModelPart.GetNodalSolutionStepVariablesList();
         mpStokesModelPart->SetBufferSize(1);
-        mpStokesModelPart->SetNodes( mpReferenceModelPart->pNodes() );
-        mpStokesModelPart->SetProcessInfo(mpReferenceModelPart->pGetProcessInfo());
-        mpStokesModelPart->SetProperties(mpReferenceModelPart->pProperties());
+        mpStokesModelPart->SetNodes( mrReferenceModelPart.pNodes() );
+        mpStokesModelPart->SetProcessInfo(mrReferenceModelPart.pGetProcessInfo());
+        mpStokesModelPart->SetProperties(mrReferenceModelPart.pProperties());
 
         // Retrieve Stokes element model
         std::string ElementName;
@@ -95,7 +97,7 @@ class  StokesInitializationProcess : public Process
         const Element& rReferenceElement = KratosComponents<Element>::Get(ElementName);
 
         // Generate Stokes elements
-        for (ModelPart::ElementsContainerType::iterator itElem = mpReferenceModelPart->ElementsBegin(); itElem != mpReferenceModelPart->ElementsEnd(); itElem++)
+        for (ModelPart::ElementsContainerType::const_iterator itElem = mrReferenceModelPart.ElementsBegin(); itElem != mrReferenceModelPart.ElementsEnd(); itElem++)
         {
             Element::Pointer pElem = rReferenceElement.Create(itElem->Id(), itElem->GetGeometry(), itElem->pGetProperties() );
             mpStokesModelPart->Elements().push_back(pElem);
@@ -232,13 +234,13 @@ class  StokesInitializationProcess : public Process
     ///@name Protected member Variables
     ///@{
 
-    const ModelPart::Pointer mpReferenceModelPart;
+    ModelPart& mrReferenceModelPart;
 
     typename TLinearSolver::Pointer mpLinearSolver;
 
     unsigned int mDomainSize;
 
-    ModelPart::Pointer mpStokesModelPart;
+    ModelPart* mpStokesModelPart;
 
     typename SolvingStrategy<TSparseSpace, TDenseSpace, TLinearSolver>::Pointer mpSolutionStrategy;
 
@@ -249,12 +251,12 @@ class  StokesInitializationProcess : public Process
     ///@{
 
     /// Protected constructor to be used by derived classes
-    StokesInitializationProcess(const ModelPart::Pointer pModelPart,
+    StokesInitializationProcess(ModelPart& rModelPart,
                                 typename TLinearSolver::Pointer pLinearSolver,
                                 unsigned int DomainSize,
                                 const StokesInitializationProcess* pThis):
         Process(),
-        mpReferenceModelPart(pModelPart),
+        mrReferenceModelPart(rModelPart),
         mpLinearSolver(pLinearSolver),
         mDomainSize(DomainSize)
     {}

applications/FluidDynamicsApplication/custom_python/add_custom_processes_to_python.cpp

@@ -61,13 +61,13 @@ void AddCustomProcessesToPython(pybind11::module& m)
 
     class_<StokesInitializationProcess< SparseSpaceType, LocalSpaceType, LinearSolverType >, StokesInitializationProcess< SparseSpaceType, LocalSpaceType, LinearSolverType >::Pointer, Process>
     (m,"StokesInitializationProcess")
-    .def(init<ModelPart::Pointer, LinearSolverType::Pointer, unsigned int, const Kratos::Variable<int>& >())
+    .def(init<ModelPart&, LinearSolverType::Pointer, unsigned int, const Kratos::Variable<int>& >())
     .def("SetConditions",&StokesInitializationProcess<SparseSpaceType, LocalSpaceType, LinearSolverType>::SetConditions)
     ;
 
     class_<BoussinesqForceProcess, BoussinesqForceProcess::Pointer, Process>
     (m,"BoussinesqForceProcess")
-    .def(init<ModelPart::Pointer, Parameters& >())
+    .def(init<ModelPart*, Parameters& >())
     ;
 
     class_<WindkesselModel, WindkesselModel::Pointer, Process>

applications/FluidDynamicsApplication/python_scripts/adjoint_fluid_analysis.py

@@ -56,8 +56,6 @@ def Initialize(self):
             self.solver.ImportModelPart()
             self.solver.AddDofs()
 
-        self.model.AddModelPart(self.main_model_part)
-
         # this should let eventual derived stages modify the model after reading.
         self.ModifyInitialProperties()
         self.ModifyInitialGeometry()
@@ -226,6 +224,5 @@ def _SetUpRestart(self):
         primal_simulation = FluidDynamicsAnalysis(model,parameters["primal_settings"])
         primal_simulation.Run()
 
-    adjoint_model = Kratos.Model()
-    adjoint_simulation = AdjointFluidAnalysis(adjoint_model,parameters["adjoint_settings"])
+    adjoint_simulation = AdjointFluidAnalysis(model,parameters["adjoint_settings"])
     adjoint_simulation.Run()

applications/FluidDynamicsApplication/python_scripts/fluid_dynamics_analysis.py

@@ -50,8 +50,6 @@ def Initialize(self):
             self.solver.ImportModelPart()
             self.solver.AddDofs()
 
-        self.model.AddModelPart(self.main_model_part)
-
         # this should let eventual derived stages modify the model after reading.
         self.ModifyInitialProperties()
         self.ModifyInitialGeometry()

applications/FluidDynamicsApplication/python_scripts/python_solvers_wrapper_fluid.py

@@ -4,7 +4,7 @@
 
 def CreateSolver(main_model_part, custom_settings):
 
-    if (type(main_model_part) != KratosMultiphysics.ModelPart):
+    if (type(main_model_part) != KratosMultiphysics.ModelPartInterface):
         raise Exception("input is expected to be provided as a Kratos ModelPart object")
 
     if (type(custom_settings) != KratosMultiphysics.Parameters):