14 #if !defined(KRATOS_TRILINOS_STRUCTURAL_MESHMOVING_STRATEGY)
15 #define KRATOS_TRILINOS_STRUCTURAL_MESHMOVING_STRATEGY
59 template <
class TSparseSpace,
86 typename TLinearSolver::Pointer pNewLinearSolver,
88 bool ReformDofSetAtEachStep =
false,
89 bool ComputeReactions =
false,
99 KRATOS_ERROR <<
"StructuralMeshMovingPart already existing when constructing TrilinosLaplacianMeshMovingStrategy";
102 m_reform_dof_set_at_each_step = ReformDofSetAtEachStep;
103 m_compute_reactions = ComputeReactions;
106 m_time_order = TimeOrder;
107 bool calculate_norm_dx_flag =
false;
116 typename SchemeType::Pointer pscheme =
typename SchemeType::Pointer(
120 BuilderSolverTypePointer builderSolver = BuilderSolverTypePointer(
124 mstrategy =
typename BaseType::Pointer(
130 m_reform_dof_set_at_each_step,
131 calculate_norm_dx_flag));
133 mstrategy->SetEchoLevel(m_echo_level);
157 i != (*mpmesh_model_part).NodesEnd();
160 (
i)->
X() = (
i)->X0();
161 (
i)->
Y() = (
i)->Y0();
162 (
i)->
Z() = (
i)->Z0();
169 if (m_reform_dof_set_at_each_step ==
true)
240 typename BaseType::Pointer mstrategy;
244 bool m_reform_dof_set_at_each_step;
245 bool m_compute_reactions;
246 bool m_calculate_mesh_velocities;
256 void GenerateMeshPart()
277 MeshElems.push_back(rElem.
Create(
278 it->Id(), it->pGetGeometry(), it->pGetProperties()));
282 CommunicatorGeneration.Execute();
Current class provides an implementation for the base builder and solving operations.
Definition: builder_and_solver.h:64
The Commmunicator class manages communication for distributed ModelPart instances.
Definition: communicator.h:67
virtual const DataCommunicator & GetDataCommunicator() const
Definition: communicator.cpp:340
Base class for all Elements.
Definition: element.h:60
virtual Pointer Create(IndexType NewId, NodesArrayType const &ThisNodes, PropertiesType::Pointer pProperties) const
It creates a new element pointer.
Definition: element.h:202
Implicit solving strategy base class This is the base class from which we will derive all the implici...
Definition: implicit_solving_strategy.h:61
static const TComponentType & Get(const std::string &rName)
Retrieves a component with the specified name.
Definition: kratos_components.h:114
ModelPart & CreateModelPart(const std::string &ModelPartName, IndexType NewBufferSize=1)
This method creates a new model part contained in the current Model with a given name and buffer size...
Definition: model.cpp:37
void DeleteModelPart(const std::string &ModelPartName)
This method deletes a modelpart with a given name.
Definition: model.cpp:64
bool HasModelPart(const std::string &rFullModelPartName) const
This method checks if a certain a model part exists given a certain name.
Definition: model.cpp:178
This class aims to manage meshes for multi-physics simulations.
Definition: model_part.h:77
MeshType::ElementsContainerType ElementsContainerType
Element container. A vector set of Elements with their Id's as key.
Definition: model_part.h:168
Communicator & GetCommunicator()
Definition: model_part.h:1821
std::string & Name()
Definition: model_part.h:1811
MeshType::NodeIterator NodeIterator
Definition: model_part.h:134
ElementsContainerType & Elements(IndexType ThisIndex=0)
Definition: model_part.h:1189
ElementIterator ElementsEnd(IndexType ThisIndex=0)
Definition: model_part.h:1179
NodesContainerType & Nodes(IndexType ThisIndex=0)
Definition: model_part.h:507
Model & GetModel()
Definition: model_part.h:323
This function recomputes the communication plan for MPI.
Definition: parallel_fill_communicator.h:56
This class provides the implementation of a static scheme.
Definition: residualbased_incrementalupdate_static_scheme.h:57
This is a very simple strategy to solve linearly the problem.
Definition: residualbased_linear_strategy.h:64
This class provides the implementation of the basic tasks that are needed by the solution strategy.
Definition: scheme.h:56
ModelPart & GetModelPart()
Operations to get the pointer to the model.
Definition: solving_strategy.h:350
Current class provides an implementation for trilinos builder and solving operations.
Definition: trilinos_block_builder_and_solver.h:85
Short class definition.
Definition: trilinos_structural_meshmoving_strategy.h:65
double Solve() override
The problem of interest is solved.
Definition: trilinos_structural_meshmoving_strategy.h:151
KRATOS_CLASS_POINTER_DEFINITION(TrilinosStructuralMeshMovingStrategy)
ImplicitSolvingStrategy< TSparseSpace, TDenseSpace, TLinearSolver > BaseType
Definition: trilinos_structural_meshmoving_strategy.h:75
void Initialize() override
Definition: trilinos_structural_meshmoving_strategy.h:148
Scheme< TSparseSpace, TDenseSpace > SchemeType
Definition: trilinos_structural_meshmoving_strategy.h:73
TrilinosStructuralMeshMovingStrategy(Epetra_MpiComm &Communicator, ModelPart &model_part, typename TLinearSolver::Pointer pNewLinearSolver, int TimeOrder=2, bool ReformDofSetAtEachStep=false, bool ComputeReactions=false, bool CalculateMeshVelocities=true, int EchoLevel=0)
Definition: trilinos_structural_meshmoving_strategy.h:84
virtual ~TrilinosStructuralMeshMovingStrategy()
Definition: trilinos_structural_meshmoving_strategy.h:140
#define KRATOS_CATCH(MoreInfo)
Definition: define.h:110
#define KRATOS_TRY
Definition: define.h:109
#define KRATOS_ERROR
Definition: exception.h:161
static int EchoLevel
Definition: co_sim_EMPIRE_API.h:42
void CalculateMeshVelocities(ModelPart &rModelPart, const TimeDiscretization::BDF1 &rBDF)
Definition: mesh_velocity_calculation.cpp:56
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
model_part
Definition: face_heat.py:14
integer i
Definition: TensorModule.f:17