85 void SetParameters(
const Vector &
u);
87 virtual ~ConductionOperator();
92int main(
int argc,
char *argv[])
95 const char *mesh_file =
"../data/star.mesh";
98 int ode_solver_type = 3;
103 bool visualization =
true;
108 cout.precision(precision);
111 args.
AddOption(&mesh_file,
"-m",
"--mesh",
112 "Mesh file to use.");
113 args.
AddOption(&ref_levels,
"-r",
"--refine",
114 "Number of times to refine the mesh uniformly.");
116 "Order (degree) of the finite elements.");
117 args.
AddOption(&ode_solver_type,
"-s",
"--ode-solver",
118 "ODE solver: 1 - Backward Euler, 2 - SDIRK2, 3 - SDIRK3,\n\t"
119 "\t 11 - Forward Euler, 12 - RK2, 13 - RK3 SSP, 14 - RK4.");
120 args.
AddOption(&t_final,
"-tf",
"--t-final",
121 "Final time; start time is 0.");
122 args.
AddOption(&dt,
"-dt",
"--time-step",
125 "Alpha coefficient.");
127 "Kappa coefficient offset.");
128 args.
AddOption(&visualization,
"-vis",
"--visualization",
"-no-vis",
129 "--no-visualization",
130 "Enable or disable GLVis visualization.");
131 args.
AddOption(&visit,
"-visit",
"--visit-datafiles",
"-no-visit",
132 "--no-visit-datafiles",
133 "Save data files for VisIt (visit.llnl.gov) visualization.");
134 args.
AddOption(&vis_steps,
"-vs",
"--visualization-steps",
135 "Visualize every n-th timestep.");
146 Mesh *mesh =
new Mesh(mesh_file, 1, 1);
153 switch (ode_solver_type)
161 case 12: ode_solver =
new RK2Solver(0.5);
break;
163 case 14: ode_solver =
new RK4Solver;
break;
170 cout <<
"Unknown ODE solver type: " << ode_solver_type <<
'\n';
178 for (
int lev = 0; lev < ref_levels; lev++)
189 cout <<
"Number of temperature unknowns: " << fe_size << endl;
205 ofstream omesh(
"ex16.mesh");
206 omesh.precision(precision);
208 ofstream osol(
"ex16-init.gf");
209 osol.precision(precision);
230 cout <<
"Unable to connect to GLVis server at "
232 visualization =
false;
233 cout <<
"GLVis visualization disabled.\n";
237 sout.precision(precision);
238 sout <<
"solution\n" << *mesh << u_gf;
241 cout <<
"GLVis visualization paused."
242 <<
" Press space (in the GLVis window) to resume it.\n";
248 ode_solver->
Init(oper);
251 bool last_step =
false;
252 for (
int ti = 1; !last_step; ti++)
254 if (
t + dt >= t_final - dt/2)
259 ode_solver->
Step(
u,
t, dt);
261 if (last_step || (ti % vis_steps) == 0)
263 cout <<
"step " << ti <<
", t = " <<
t << endl;
268 sout <<
"solution\n" << *mesh << u_gf << flush;
278 oper.SetParameters(
u);
284 ofstream osol(
"ex16-final.gf");
285 osol.precision(precision);
299 M(NULL), K(NULL), T(NULL), current_dt(0.0), z(height)
301 const real_t rel_tol = 1e-8;
306 M->FormSystemMatrix(ess_tdof_list, Mmat);
308 M_solver.iterative_mode =
false;
309 M_solver.SetRelTol(rel_tol);
310 M_solver.SetAbsTol(0.0);
311 M_solver.SetMaxIter(30);
312 M_solver.SetPrintLevel(0);
313 M_solver.SetPreconditioner(M_prec);
314 M_solver.SetOperator(Mmat);
319 T_solver.iterative_mode =
false;
320 T_solver.SetRelTol(rel_tol);
321 T_solver.SetAbsTol(0.0);
322 T_solver.SetMaxIter(100);
323 T_solver.SetPrintLevel(0);
324 T_solver.SetPreconditioner(T_prec);
329void ConductionOperator::Mult(
const Vector &
u,
Vector &du_dt)
const
336 M_solver.
Mult(z, du_dt);
339void ConductionOperator::ImplicitSolve(
const real_t dt,
347 T =
Add(1.0, Mmat, dt, Kmat);
351 MFEM_VERIFY(dt == current_dt,
"");
354 T_solver.
Mult(z, du_dt);
357void ConductionOperator::SetParameters(
const Vector &
u)
360 u_alpha_gf.SetFromTrueDofs(
u);
361 for (
int i = 0; i < u_alpha_gf.Size(); i++)
363 u_alpha_gf(i) = kappa + alpha*u_alpha_gf(i);
378ConductionOperator::~ConductionOperator()
Backward Euler ODE solver. L-stable.
Conjugate gradient method.
virtual void SetOperator(const Operator &op)
Also calls SetOperator for the preconditioner if there is one.
virtual void Mult(const Vector &b, Vector &x) const
Iterative solution of the linear system using the Conjugate Gradient method.
Data type for scaled Jacobi-type smoother of sparse matrix.
void SetCycle(int c)
Set time cycle (for time-dependent simulations)
void SetTime(real_t t)
Set physical time (for time-dependent simulations)
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
virtual int GetTrueVSize() const
Return the number of vector true (conforming) dofs.
The classical forward Euler method.
A general function coefficient.
Coefficient defined by a GridFunction. This coefficient is mesh dependent.
Class for grid function - Vector with associated FE space.
virtual void Save(std::ostream &out) const
Save the GridFunction to an output stream.
virtual void SetFromTrueDofs(const Vector &tv)
Set the GridFunction from the given true-dof vector.
virtual void ProjectCoefficient(Coefficient &coeff)
Project coeff Coefficient to this GridFunction. The projection computation depends on the choice of t...
void GetTrueDofs(Vector &tv) const
Extract the true-dofs from the GridFunction.
Arbitrary order H1-conforming (continuous) finite elements.
Implicit midpoint method. A-stable, not L-stable.
virtual void Print(std::ostream &os=mfem::out, const std::string &comments="") const
int Dimension() const
Dimension of the reference space used within the elements.
void UniformRefinement(int i, const DSTable &, int *, int *, int *)
Abstract class for solving systems of ODEs: dx/dt = f(x,t)
virtual void Init(TimeDependentOperator &f_)
Associate a TimeDependentOperator with the ODE solver.
virtual void Step(Vector &x, real_t &t, real_t &dt)=0
Perform a time step from time t [in] to time t [out] based on the requested step size dt [in].
void Parse()
Parse the command-line options. Note that this function expects all the options provided through the ...
void PrintUsage(std::ostream &out) const
Print the usage message.
void PrintOptions(std::ostream &out) const
Print the options.
void AddOption(bool *var, const char *enable_short_name, const char *enable_long_name, const char *disable_short_name, const char *disable_long_name, const char *description, bool required=false)
Add a boolean option and set 'var' to receive the value. Enable/disable tags are used to set the bool...
bool Good() const
Return true if the command line options were parsed successfully.
Third-order, strong stability preserving (SSP) Runge-Kutta method.
The classical explicit forth-order Runge-Kutta method, RK4.
virtual void Mult(const Vector &x, Vector &y) const
Matrix vector multiplication.
Base abstract class for first order time dependent operators.
void Neg()
(*this) = -(*this)
real_t Norml2() const
Returns the l2 norm of the vector.
Data collection with VisIt I/O routines.
virtual void Save() override
Save the collection and a VisIt root file.
virtual void RegisterField(const std::string &field_name, GridFunction *gf) override
Add a grid function to the collection and update the root file.
int open(const char hostname[], int port)
Open the socket stream on 'port' at 'hostname'.
real_t InitialTemperature(const Vector &x)
real_t u(const Vector &xvec)
std::function< real_t(const Vector &)> f(real_t mass_coeff)
void Add(const DenseMatrix &A, const DenseMatrix &B, real_t alpha, DenseMatrix &C)
C = A + alpha*B.