libMesh
Loading...
Searching...
No Matches
Functions
fem_system_ex2.C File Reference

Go to the source code of this file.

Functions

void setup (EquationSystems &systems, Mesh &mesh, GetPot &args)
 
void run_timestepping (EquationSystems &systems, GetPot &args)
 
int main (int argc, char **argv)
 

Function Documentation

◆ main()

int main ( int  argc,
char **  argv 
)

Definition at line 157 of file fem_system_ex2.C.

158{
159 // Initialize libMesh and any dependent libraries
160 LibMeshInit init(argc, argv);
161
162 // This example requires a linear solver package.
163 libmesh_example_requires(libMesh::default_solver_package() != INVALID_SOLVER_PACKAGE,
164 "--enable-petsc, --enable-trilinos, or --enable-eigen");
165
166 // Trilinos gives us an inverted element on this one...
167 libmesh_example_requires(libMesh::default_solver_package() != TRILINOS_SOLVERS, "--enable-petsc");
168
169 // Threaded assembly doesn't currently work with the moving mesh
170 // code.
171 // We'll skip this example for now.
172 if (libMesh::n_threads() > 1)
173 {
174 libMesh::out << "We skip fem_system_ex2 when using threads.\n"
175 << std::endl;
176 return 0;
177 }
178
179 // read simulation parameters from file
180 GetPot args = GetPot("solid.in");
181
182 // But allow the command line to override
183 args.parse_command_line(argc, argv);
184
185 // Create System and Mesh
186 int dim = args("mesh/generation/dimension", 3);
187 libmesh_example_requires(dim <= LIBMESH_DIM, "3D support");
188
189 // Create a mesh distributed across the default MPI communicator.
190 Mesh mesh(init.comm(), cast_int<unsigned char>(dim));
191
192 EquationSystems systems(mesh);
193
194 // Create and set systems up
195 setup(systems, mesh, args);
196
197 // run the systems
198 run_timestepping(systems, args);
199
200 out << "Finished calculations" << std::endl;
201 return 0;
202}
unsigned int dim
This is the EquationSystems class.
The LibMeshInit class, when constructed, initializes the dependent libraries (e.g.
Definition libmesh.h:92
The Mesh class is a thin wrapper, around the ReplicatedMesh class by default.
Definition mesh.h:51
void run_timestepping(EquationSystems &systems, GetPot &args)
void setup(EquationSystems &systems, Mesh &mesh, GetPot &args)
MeshBase & mesh
void init(triangulateio &t)
Initializes the fields of t to nullptr/0 as necessary.
SolverPackage default_solver_package()
Definition libmesh.C:1064
OStreamProxy out
unsigned int n_threads()

References libMesh::default_solver_package(), dim, libMesh::INVALID_SOLVER_PACKAGE, main(), mesh, libMesh::n_threads(), libMesh::out, run_timestepping(), setup(), and libMesh::TRILINOS_SOLVERS.

◆ run_timestepping()

void run_timestepping ( EquationSystems systems,
GetPot &  args 
)

Definition at line 93 of file fem_system_ex2.C.

94{
95 TransientExplicitSystem & aux_system = systems.get_system<TransientExplicitSystem>("auxiliary");
96
97 SolidSystem & solid_system = systems.get_system<SolidSystem>("solid");
98
99#ifdef LIBMESH_HAVE_VTK
100 std::unique_ptr<VTKIO> io = std::make_unique<VTKIO>(systems.get_mesh());
101#endif
102
103 Real duration = args("duration", 1.0);
104
105 for (unsigned int t_step = 0; t_step < duration/solid_system.deltat; t_step++)
106 {
107 // Progress in current phase [0..1]
108 Real progress = t_step * solid_system.deltat / duration;
109 systems.parameters.set<Real>("progress") = progress;
110 systems.parameters.set<unsigned int>("step") = t_step;
111
112 // Update message
113
114 out << "===== Time Step " << std::setw(4) << t_step;
115 out << " (" << std::fixed << std::setprecision(2) << std::setw(6) << (progress*100.) << "%)";
116 out << ", time = " << std::setw(7) << solid_system.time;
117 out << " =====" << std::endl;
118
119 // Advance timestep in auxiliary system
120 aux_system.current_local_solution->close();
121 aux_system.old_local_solution->close();
122 *aux_system.older_local_solution = *aux_system.old_local_solution;
123 *aux_system.old_local_solution = *aux_system.current_local_solution;
124
125 out << "Solving Solid" << std::endl;
126 solid_system.solve();
127 aux_system.current_local_solution->close();
128 *aux_system.solution = *aux_system.current_local_solution;
129 aux_system.solution->close();
130
131 // Carry out the adaptive mesh refinement/coarsening
132 out << "Doing a reinit of the equation systems" << std::endl;
133 systems.reinit();
134
135#ifdef LIBMESH_HAVE_VTK
136 if (t_step % args("output/frequency", 1) == 0)
137 {
138 std::stringstream file_name;
139 file_name << args("results_directory", "./")
140 << "fem_"
141 << std::setw(6)
142 << std::setfill('0')
143 << t_step
144 << ".pvtu";
145
146 io->write_equation_systems(file_name.str(), systems);
147 }
148#endif
149 // Advance to the next timestep in a transient problem
150 out << "Advancing to next step" << std::endl;
151 solid_system.time_solver->advance_timestep();
152 }
153}
virtual void reinit()
Handle any mesh changes and reinitialize all the systems on the updated mesh.
const MeshBase & get_mesh() const
Parameters parameters
Data structure holding arbitrary parameters.
const T_sys & get_system(std::string_view name) const
virtual void close()=0
Calls the NumericVector's internal assembly routines, ensuring that the values are consistent across ...
T & set(const std::string &)
Definition parameters.h:494
Manages storage and variables for transient systems.
NumericVector< Number > * older_local_solution
All the values I need to compute my contribution to the simulation at hand.
NumericVector< Number > * old_local_solution
All the values I need to compute my contribution to the simulation at hand.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

References libMesh::NumericVector< T >::close(), libMesh::EquationSystems::get_mesh(), libMesh::EquationSystems::get_system(), libMesh::TransientSystem< Base >::old_local_solution, libMesh::TransientSystem< Base >::older_local_solution, libMesh::out, libMesh::EquationSystems::parameters, libMesh::Real, libMesh::EquationSystems::reinit(), and libMesh::Parameters::set().

Referenced by main().

◆ setup()

void setup ( EquationSystems systems,
Mesh mesh,
GetPot &  args 
)

Definition at line 49 of file fem_system_ex2.C.

52{
53 const unsigned int dim = mesh.mesh_dimension();
54 // We currently invert tensors with the assumption that they're 3x3
55 libmesh_assert (dim == 3);
56
57 // Generating Mesh
58 ElemType eltype = Utility::string_to_enum<ElemType>(args("mesh/generation/element_type", "hex8"));
59 int nx = args("mesh/generation/num_elem", 4, 0);
60 int ny = args("mesh/generation/num_elem", 4, 1);
61 int nz = dim > 2 ? args("mesh/generation/num_elem", 4, 2) : 0;
62 double origx = args("mesh/generation/origin", -1.0, 0);
63 double origy = args("mesh/generation/origin", -1.0, 1);
64 double origz = args("mesh/generation/origin", 0.0, 2);
65 double sizex = args("mesh/generation/size", 2.0, 0);
66 double sizey = args("mesh/generation/size", 2.0, 1);
67 double sizez = args("mesh/generation/size", 2.0, 2);
69 origx, origx+sizex, origy, origy+sizey, origz, origz+sizez, eltype);
71
72 // Creating Systems
73 SolidSystem & imms = systems.add_system<SolidSystem> ("solid");
74 imms.args = args;
75
76 // Build up auxiliary system
77 ExplicitSystem & aux_sys = systems.add_system<TransientExplicitSystem>("auxiliary");
78
79 // Initialize the system
80 systems.parameters.set<unsigned int>("phase") = 0;
81 systems.init();
82
83 imms.save_initial_mesh();
84
85 // Fill global solution vector from local ones
86 aux_sys.current_local_solution->close();
87 *aux_sys.solution = *aux_sys.current_local_solution;
88 aux_sys.solution->close();
89}
void save_initial_mesh()
virtual void init()
Initialize all the systems.
virtual System & add_system(std::string_view system_type, std::string_view name)
Add the system of type system_type named name to the systems array.
Manages consistently variables, degrees of freedom, and coefficient vectors for explicit systems.
unsigned int mesh_dimension() const
Definition mesh_base.C:430
void print_info(std::ostream &os=libMesh::out, const unsigned int verbosity=0, const bool global=true) const
Prints relevant information about the mesh.
Definition mesh_base.C:1755
std::unique_ptr< NumericVector< Number > > current_local_solution
All the values I need to compute my contribution to the simulation at hand.
Definition system.h:1667
std::unique_ptr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition system.h:1655
void build_cube(UnstructuredMesh &mesh, const unsigned int nx=0, const unsigned int ny=0, const unsigned int nz=0, const Real xmin=0., const Real xmax=1., const Real ymin=0., const Real ymax=1., const Real zmin=0., const Real zmax=1., const ElemType type=INVALID_ELEM, const bool gauss_lobatto_grid=false)
Builds a (elements) cube.
ElemType
Defines an enum for geometric element types.
libmesh_assert(ctx)

References libMesh::EquationSystems::add_system(), SolidSystem::args, libMesh::MeshTools::Generation::build_cube(), libMesh::System::current_local_solution, dim, libMesh::EquationSystems::init(), libMesh::libmesh_assert(), mesh, libMesh::MeshBase::mesh_dimension(), libMesh::EquationSystems::parameters, libMesh::MeshBase::print_info(), SolidSystem::save_initial_mesh(), libMesh::Parameters::set(), and libMesh::System::solution.

Referenced by main().