MFEM  v4.1.0 Finite element discretization library
ex19p.cpp
Go to the documentation of this file.
1 // MFEM Example 19 - Parallel Version
2 //
3 // Compile with: make ex19p
4 //
5 // Sample runs:
6 // mpirun -np 2 ex19p -m ../data/beam-quad.mesh
7 // mpirun -np 2 ex19p -m ../data/beam-tri.mesh
8 // mpirun -np 2 ex19p -m ../data/beam-hex.mesh
9 // mpirun -np 2 ex19p -m ../data/beam-tet.mesh
10 // mpirun -np 2 ex19p -m ../data/beam-wedge.mesh
11 //
12 // Description: This examples solves a quasi-static incompressible nonlinear
13 // elasticity problem of the form 0 = H(x), where H is an
14 // incompressible hyperelastic model and x is a block state vector
15 // containing displacement and pressure variables. The geometry of
16 // the domain is assumed to be as follows:
17 //
18 // +---------------------+
19 // boundary --->| |<--- boundary
20 // attribute 1 | | attribute 2
21 // (fixed) +---------------------+ (fixed, nonzero)
22 //
23 // The example demonstrates the use of block nonlinear operators
24 // (the class RubberOperator defining H(x)) as well as a nonlinear
25 // Newton solver for the quasi-static problem. Each Newton step
26 // requires the inversion of a Jacobian matrix, which is done
27 // through a (preconditioned) inner solver. The specialized block
28 // preconditioner is implemented as a user-defined solver.
29 //
30 // We recommend viewing examples 2, 5, and 10 before viewing this
31 // example.
32
33 #include "mfem.hpp"
34 #include <memory>
35 #include <iostream>
36 #include <fstream>
37
38 using namespace std;
39 using namespace mfem;
40
41 // Custom block preconditioner for the Jacobian of the incompressible nonlinear
42 // elasticity operator. It has the form
43 //
44 // P^-1 = [ K^-1 0 ][ I -B^T ][ I 0 ]
45 // [ 0 I ][ 0 I ][ 0 -\gamma S^-1 ]
46 //
47 // where the original Jacobian has the form
48 //
49 // J = [ K B^T ]
50 // [ B 0 ]
51 //
52 // and K^-1 is an approximation of the inverse of the displacement part of the
53 // Jacobian and S^-1 is an approximation of the inverse of the Schur
54 // complement S = B K^-1 B^T. The Schur complement is approximated using
55 // a mass matrix of the pressure variables.
56 class JacobianPreconditioner : public Solver
57 {
58 protected:
59  // Finite element spaces for setting up preconditioner blocks
61
62  // Offsets for extracting block vector segments
63  Array<int> &block_trueOffsets;
64
65  // Jacobian for block access
66  BlockOperator *jacobian;
67
68  // Scaling factor for the pressure mass matrix in the block preconditioner
69  double gamma;
70
71  // Objects for the block preconditioner application
72  Operator *pressure_mass;
73  Solver *mass_pcg;
74  Solver *mass_prec;
75  Solver *stiff_pcg;
76  Solver *stiff_prec;
77
78 public:
79  JacobianPreconditioner(Array<ParFiniteElementSpace *> &fes,
80  Operator &mass, Array<int> &offsets);
81
82  virtual void Mult(const Vector &k, Vector &y) const;
83  virtual void SetOperator(const Operator &op);
84
85  virtual ~JacobianPreconditioner();
86 };
87
88 // After spatial discretization, the rubber model can be written as:
89 // 0 = H(x)
90 // where x is the block vector representing the deformation and pressure and
91 // H(x) is the nonlinear incompressible neo-Hookean operator.
92 class RubberOperator : public Operator
93 {
94 protected:
95  // Finite element spaces
97
98  // Block nonlinear form
99  ParBlockNonlinearForm *Hform;
100
101  // Pressure mass matrix for the preconditioner
102  Operator *pressure_mass;
103
104  // Newton solver for the hyperelastic operator
105  NewtonSolver newton_solver;
106
107  // Solver for the Jacobian solve in the Newton method
108  Solver *j_solver;
109
110  // Preconditioner for the Jacobian
111  Solver *j_prec;
112
113  // Shear modulus coefficient
114  Coefficient &mu;
115
116  // Block offsets for variable access
117  Array<int> &block_trueOffsets;
118
119 public:
120  RubberOperator(Array<ParFiniteElementSpace *> &fes, Array<Array<int> *>&ess_bdr,
121  Array<int> &block_trueOffsets, double rel_tol, double abs_tol,
122  int iter, Coefficient &mu);
123
124  // Required to use the native newton solver
125  virtual Operator &GetGradient(const Vector &xp) const;
126  virtual void Mult(const Vector &k, Vector &y) const;
127
128  // Driver for the newton solver
129  void Solve(Vector &xp) const;
130
131  virtual ~RubberOperator();
132 };
133
134 // Visualization driver
135 void visualize(ostream &out, ParMesh *mesh, ParGridFunction *deformed_nodes,
136  ParGridFunction *field, const char *field_name = NULL,
137  bool init_vis = false);
138
139 // Configuration definition functions
140 void ReferenceConfiguration(const Vector &x, Vector &y);
141 void InitialDeformation(const Vector &x, Vector &y);
142
143
144 int main(int argc, char *argv[])
145 {
146  // 1. Initialize MPI
147  int num_procs, myid;
148  MPI_Init(&argc, &argv);
149  MPI_Comm_size(MPI_COMM_WORLD, &num_procs);
150  MPI_Comm_rank(MPI_COMM_WORLD, &myid);
151
152  // 2. Parse command-line options
153  const char *mesh_file = "../data/beam-tet.mesh";
154  int ser_ref_levels = 0;
155  int par_ref_levels = 0;
156  int order = 2;
157  bool visualization = true;
158  double newton_rel_tol = 1e-4;
159  double newton_abs_tol = 1e-6;
160  int newton_iter = 500;
161  double mu = 1.0;
162
163  OptionsParser args(argc, argv);
164  args.AddOption(&mesh_file, "-m", "--mesh",
165  "Mesh file to use.");
166  args.AddOption(&ser_ref_levels, "-rs", "--refine-serial",
167  "Number of times to refine the mesh uniformly in serial.");
168  args.AddOption(&par_ref_levels, "-rp", "--refine-parallel",
169  "Number of times to refine the mesh uniformly in parallel.");
170  args.AddOption(&order, "-o", "--order",
171  "Order (degree) of the finite elements.");
172  args.AddOption(&visualization, "-vis", "--visualization", "-no-vis",
173  "--no-visualization",
174  "Enable or disable GLVis visualization.");
175  args.AddOption(&newton_rel_tol, "-rel", "--relative-tolerance",
176  "Relative tolerance for the Newton solve.");
177  args.AddOption(&newton_abs_tol, "-abs", "--absolute-tolerance",
178  "Absolute tolerance for the Newton solve.");
179  args.AddOption(&newton_iter, "-it", "--newton-iterations",
180  "Maximum iterations for the Newton solve.");
181  args.AddOption(&mu, "-mu", "--shear-modulus",
182  "Shear modulus for the neo-Hookean material.");
183  args.Parse();
184  if (!args.Good())
185  {
186  if (myid == 0)
187  {
188  args.PrintUsage(cout);
189  }
190  MPI_Finalize();
191  return 1;
192  }
193  if (myid == 0)
194  {
195  args.PrintOptions(cout);
196  }
197
198  // 3. Read the (serial) mesh from the given mesh file on all processors. We
199  // can handle triangular, quadrilateral, tetrahedral and hexahedral meshes
200  // with the same code.
201  Mesh *mesh = new Mesh(mesh_file, 1, 1);
202  int dim = mesh->Dimension();
203
204  // 4. Refine the mesh in serial to increase the resolution. In this example
205  // we do 'ser_ref_levels' of uniform refinement, where 'ser_ref_levels' is
206  // a command-line parameter.
207  for (int lev = 0; lev < ser_ref_levels; lev++)
208  {
209  mesh->UniformRefinement();
210  }
211
212  // 5. Define a parallel mesh by a partitioning of the serial mesh. Refine
213  // this mesh further in parallel to increase the resolution. Once the
214  // parallel mesh is defined, the serial mesh can be deleted.
215  ParMesh *pmesh = new ParMesh(MPI_COMM_WORLD, *mesh);
216  delete mesh;
217  for (int lev = 0; lev < par_ref_levels; lev++)
218  {
219  pmesh->UniformRefinement();
220  }
221
222  // 6. Define the shear modulus for the incompressible Neo-Hookean material
223  ConstantCoefficient c_mu(mu);
224
225  // 7. Define the finite element spaces for displacement and pressure
226  // (Taylor-Hood elements). By default, the displacement (u/x) is a second
227  // order vector field, while the pressure (p) is a linear scalar function.
228  H1_FECollection quad_coll(order, dim);
229  H1_FECollection lin_coll(order-1, dim);
230
231  ParFiniteElementSpace R_space(pmesh, &quad_coll, dim, Ordering::byVDIM);
232  ParFiniteElementSpace W_space(pmesh, &lin_coll);
233
235  spaces[0] = &R_space;
236  spaces[1] = &W_space;
237
238  HYPRE_Int glob_R_size = R_space.GlobalTrueVSize();
239  HYPRE_Int glob_W_size = W_space.GlobalTrueVSize();
240
241  // 8. Define the Dirichlet conditions (set to boundary attribute 1 and 2)
242  Array<Array<int> *> ess_bdr(2);
243
244  Array<int> ess_bdr_u(R_space.GetMesh()->bdr_attributes.Max());
245  Array<int> ess_bdr_p(W_space.GetMesh()->bdr_attributes.Max());
246
247  ess_bdr_p = 0;
248  ess_bdr_u = 0;
249  ess_bdr_u[0] = 1;
250  ess_bdr_u[1] = 1;
251
252  ess_bdr[0] = &ess_bdr_u;
253  ess_bdr[1] = &ess_bdr_p;
254
255  // 9. Print the mesh statistics
256  if (myid == 0)
257  {
258  std::cout << "***********************************************************\n";
259  std::cout << "dim(u) = " << glob_R_size << "\n";
260  std::cout << "dim(p) = " << glob_W_size << "\n";
261  std::cout << "dim(u+p) = " << glob_R_size + glob_W_size << "\n";
262  std::cout << "***********************************************************\n";
263  }
264
265  // 10. Define the block structure of the solution vector (u then p)
266  Array<int> block_trueOffsets(3);
267  block_trueOffsets[0] = 0;
268  block_trueOffsets[1] = R_space.TrueVSize();
269  block_trueOffsets[2] = W_space.TrueVSize();
270  block_trueOffsets.PartialSum();
271
272  BlockVector xp(block_trueOffsets);
273
274  // 11. Define grid functions for the current configuration, reference
275  // configuration, final deformation, and pressure
276  ParGridFunction x_gf(&R_space);
277  ParGridFunction x_ref(&R_space);
278  ParGridFunction x_def(&R_space);
279  ParGridFunction p_gf(&W_space);
280
283
284  x_gf.ProjectCoefficient(deform);
285  x_ref.ProjectCoefficient(refconfig);
286  p_gf = 0.0;
287
288  // 12. Set up the block solution vectors
289  x_gf.GetTrueDofs(xp.GetBlock(0));
290  p_gf.GetTrueDofs(xp.GetBlock(1));
291
292  // 13. Initialize the incompressible neo-Hookean operator
293  RubberOperator oper(spaces, ess_bdr, block_trueOffsets,
294  newton_rel_tol, newton_abs_tol, newton_iter, c_mu);
295
296  // 14. Solve the Newton system
297  oper.Solve(xp);
298
299  // 15. Distribute the shared degrees of freedom
300  x_gf.Distribute(xp.GetBlock(0));
301  p_gf.Distribute(xp.GetBlock(1));
302
303  // 16. Compute the final deformation
304  subtract(x_gf, x_ref, x_def);
305
306  // 17. Visualize the results if requested
307  socketstream vis_u, vis_p;
308  if (visualization)
309  {
310  char vishost[] = "localhost";
311  int visport = 19916;
312  vis_u.open(vishost, visport);
313  vis_u.precision(8);
314  visualize(vis_u, pmesh, &x_gf, &x_def, "Deformation", true);
315  // Make sure all ranks have sent their 'u' solution before initiating
316  // another set of GLVis connections (one from each rank):
317  MPI_Barrier(pmesh->GetComm());
318  vis_p.open(vishost, visport);
319  vis_p.precision(8);
320  visualize(vis_p, pmesh, &x_gf, &p_gf, "Pressure", true);
321  }
322
323  // 18. Save the displaced mesh, the final deformation, and the pressure
324  {
325  GridFunction *nodes = &x_gf;
326  int owns_nodes = 0;
327  pmesh->SwapNodes(nodes, owns_nodes);
328
329  ostringstream mesh_name, pressure_name, deformation_name;
330  mesh_name << "mesh." << setfill('0') << setw(6) << myid;
331  pressure_name << "pressure." << setfill('0') << setw(6) << myid;
332  deformation_name << "deformation." << setfill('0') << setw(6) << myid;
333
334  ofstream mesh_ofs(mesh_name.str().c_str());
335  mesh_ofs.precision(8);
336  pmesh->Print(mesh_ofs);
337
338  ofstream pressure_ofs(pressure_name.str().c_str());
339  pressure_ofs.precision(8);
340  p_gf.Save(pressure_ofs);
341
342  ofstream deformation_ofs(deformation_name.str().c_str());
343  deformation_ofs.precision(8);
344  x_def.Save(deformation_ofs);
345  }
346
347  // 19. Free the used memory
348  delete pmesh;
349
350  MPI_Finalize();
351
352  return 0;
353 }
354
355
356 JacobianPreconditioner::JacobianPreconditioner(Array<ParFiniteElementSpace *>
357  &fes,
358  Operator &mass,
359  Array<int> &offsets)
360  : Solver(offsets[2]), block_trueOffsets(offsets), pressure_mass(&mass)
361 {
362  fes.Copy(spaces);
363
364  gamma = 0.00001;
365
366  // The mass matrix and preconditioner do not change every Newton cycle, so
367  // we only need to define them once
368  HypreBoomerAMG *mass_prec_amg = new HypreBoomerAMG();
369  mass_prec_amg->SetPrintLevel(0);
370
371  mass_prec = mass_prec_amg;
372
373  CGSolver *mass_pcg_iter = new CGSolver(spaces[0]->GetComm());
374  mass_pcg_iter->SetRelTol(1e-12);
375  mass_pcg_iter->SetAbsTol(1e-12);
376  mass_pcg_iter->SetMaxIter(200);
377  mass_pcg_iter->SetPrintLevel(0);
378  mass_pcg_iter->SetPreconditioner(*mass_prec);
379  mass_pcg_iter->SetOperator(*pressure_mass);
380  mass_pcg_iter->iterative_mode = false;
381
382  mass_pcg = mass_pcg_iter;
383
384  // The stiffness matrix does change every Newton cycle, so we will define it
385  // during SetOperator
386  stiff_pcg = NULL;
387  stiff_prec = NULL;
388 }
389
390 void JacobianPreconditioner::Mult(const Vector &k, Vector &y) const
391 {
392  // Extract the blocks from the input and output vectors
393  Vector disp_in(k.GetData() + block_trueOffsets[0],
394  block_trueOffsets[1]-block_trueOffsets[0]);
395  Vector pres_in(k.GetData() + block_trueOffsets[1],
396  block_trueOffsets[2]-block_trueOffsets[1]);
397
398  Vector disp_out(y.GetData() + block_trueOffsets[0],
399  block_trueOffsets[1]-block_trueOffsets[0]);
400  Vector pres_out(y.GetData() + block_trueOffsets[1],
401  block_trueOffsets[2]-block_trueOffsets[1]);
402
403  Vector temp(block_trueOffsets[1]-block_trueOffsets[0]);
404  Vector temp2(block_trueOffsets[1]-block_trueOffsets[0]);
405
406  // Perform the block elimination for the preconditioner
407  mass_pcg->Mult(pres_in, pres_out);
408  pres_out *= -gamma;
409
410  jacobian->GetBlock(0,1).Mult(pres_out, temp);
411  subtract(disp_in, temp, temp2);
412
413  stiff_pcg->Mult(temp2, disp_out);
414 }
415
416 void JacobianPreconditioner::SetOperator(const Operator &op)
417 {
418  jacobian = (BlockOperator *) &op;
419
420  // Initialize the stiffness preconditioner and solver
421  if (stiff_prec == NULL)
422  {
423  HypreBoomerAMG *stiff_prec_amg = new HypreBoomerAMG();
424  stiff_prec_amg->SetPrintLevel(0);
425  stiff_prec_amg->SetElasticityOptions(spaces[0]);
426
427  stiff_prec = stiff_prec_amg;
428
429  GMRESSolver *stiff_pcg_iter = new GMRESSolver(spaces[0]->GetComm());
430  stiff_pcg_iter->SetRelTol(1e-8);
431  stiff_pcg_iter->SetAbsTol(1e-8);
432  stiff_pcg_iter->SetMaxIter(200);
433  stiff_pcg_iter->SetPrintLevel(0);
434  stiff_pcg_iter->SetPreconditioner(*stiff_prec);
435  stiff_pcg_iter->iterative_mode = false;
436
437  stiff_pcg = stiff_pcg_iter;
438  }
439
440  // At each Newton cycle, compute the new stiffness AMG preconditioner by
441  // updating the iterative solver which, in turn, updates its preconditioner
442  stiff_pcg->SetOperator(jacobian->GetBlock(0,0));
443 }
444
445 JacobianPreconditioner::~JacobianPreconditioner()
446 {
447  delete mass_pcg;
448  delete mass_prec;
449  delete stiff_prec;
450  delete stiff_pcg;
451 }
452
453
454 RubberOperator::RubberOperator(Array<ParFiniteElementSpace *> &fes,
455  Array<Array<int> *> &ess_bdr,
456  Array<int> &trueOffsets,
457  double rel_tol,
458  double abs_tol,
459  int iter,
460  Coefficient &c_mu)
461  : Operator(fes[0]->TrueVSize() + fes[1]->TrueVSize()),
462  newton_solver(fes[0]->GetComm()), mu(c_mu), block_trueOffsets(trueOffsets)
463 {
464  Array<Vector *> rhs(2);
465  rhs = NULL; // Set all entries in the array
466
467  fes.Copy(spaces);
468
469  // Define the block nonlinear form
470  Hform = new ParBlockNonlinearForm(spaces);
471
472  // Add the incompressible neo-Hookean integrator
473  Hform->AddDomainIntegrator(new IncompressibleNeoHookeanIntegrator(mu));
474
475  // Set the essential boundary conditions
476  Hform->SetEssentialBC(ess_bdr, rhs);
477
478  // Compute the pressure mass stiffness matrix
479  ParBilinearForm *a = new ParBilinearForm(spaces[1]);
480  ConstantCoefficient one(1.0);
481  OperatorHandle mass(Operator::Hypre_ParCSR);
482  a->AddDomainIntegrator(new MassIntegrator(one));
483  a->Assemble();
484  a->Finalize();
485  a->ParallelAssemble(mass);
486  delete a;
487
488  mass.SetOperatorOwner(false);
489  pressure_mass = mass.Ptr();
490
491  // Initialize the Jacobian preconditioner
492  JacobianPreconditioner *jac_prec =
493  new JacobianPreconditioner(fes, *pressure_mass, block_trueOffsets);
494  j_prec = jac_prec;
495
496  // Set up the Jacobian solver
497  GMRESSolver *j_gmres = new GMRESSolver(spaces[0]->GetComm());
498  j_gmres->iterative_mode = false;
499  j_gmres->SetRelTol(1e-12);
500  j_gmres->SetAbsTol(1e-12);
501  j_gmres->SetMaxIter(300);
502  j_gmres->SetPrintLevel(0);
503  j_gmres->SetPreconditioner(*j_prec);
504  j_solver = j_gmres;
505
506  // Set the newton solve parameters
507  newton_solver.iterative_mode = true;
508  newton_solver.SetSolver(*j_solver);
509  newton_solver.SetOperator(*this);
510  newton_solver.SetPrintLevel(1);
511  newton_solver.SetRelTol(rel_tol);
512  newton_solver.SetAbsTol(abs_tol);
513  newton_solver.SetMaxIter(iter);
514 }
515
516 // Solve the Newton system
517 void RubberOperator::Solve(Vector &xp) const
518 {
519  Vector zero;
520  newton_solver.Mult(zero, xp);
521  MFEM_VERIFY(newton_solver.GetConverged(),
522  "Newton Solver did not converge.");
523 }
524
525 // compute: y = H(x,p)
526 void RubberOperator::Mult(const Vector &k, Vector &y) const
527 {
528  Hform->Mult(k, y);
529 }
530
531 // Compute the Jacobian from the nonlinear form
532 Operator &RubberOperator::GetGradient(const Vector &xp) const
533 {
534  return Hform->GetGradient(xp);
535 }
536
537 RubberOperator::~RubberOperator()
538 {
539  delete Hform;
540  delete pressure_mass;
541  delete j_solver;
542  delete j_prec;
543 }
544
545
546 // Inline visualization
547 void visualize(ostream &out, ParMesh *mesh, ParGridFunction *deformed_nodes,
548  ParGridFunction *field, const char *field_name, bool init_vis)
549 {
550  if (!out)
551  {
552  return;
553  }
554
555  GridFunction *nodes = deformed_nodes;
556  int owns_nodes = 0;
557
558  mesh->SwapNodes(nodes, owns_nodes);
559
560  out << "parallel " << mesh->GetNRanks() << " " << mesh->GetMyRank() << "\n";
561  out << "solution\n" << *mesh << *field;
562
563  mesh->SwapNodes(nodes, owns_nodes);
564
565  if (init_vis)
566  {
567  out << "window_size 800 800\n";
568  out << "window_title '" << field_name << "'\n";
569  if (mesh->SpaceDimension() == 2)
570  {
571  out << "view 0 0\n"; // view from top
572  out << "keys jlA\n"; // turn off perspective and light, +anti-aliasing
573  }
574  out << "keys cmA\n"; // show colorbar and mesh, +anti-aliasing
575  out << "autoscale value\n"; // update value-range; keep mesh-extents fixed
576  }
577  out << flush;
578 }
579
580 void ReferenceConfiguration(const Vector &x, Vector &y)
581 {
582  // Set the reference, stress free, configuration
583  y = x;
584 }
585
586 void InitialDeformation(const Vector &x, Vector &y)
587 {
588  // Set the initial configuration. Having this different from the reference
589  // configuration can help convergence
590  y = x;
591  y[1] = x[1] + 0.25*x[0];
592 }
void visualize(ostream &out, Mesh *mesh, GridFunction *deformed_nodes, GridFunction *field, const char *field_name=NULL, bool init_vis=false)
Definition: ex10.cpp:379
void InitialDeformation(const Vector &x, Vector &y)
Definition: ex10.cpp:591
Conjugate gradient method.
Definition: solvers.hpp:152
void ReferenceConfiguration(const Vector &x, Vector &y)
Definition: ex19.cpp:526
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:27
A class to handle Vectors in a block fashion.
Definition: blockvector.hpp:30
Subclass constant coefficient.
Definition: coefficient.hpp:67
virtual Operator & GetGradient(const Vector &x) const
Evaluate the gradient operator at the point x. The default behavior in class Operator is to generate ...
Definition: operator.hpp:89
void SwapNodes(GridFunction *&nodes, int &own_nodes_)
Definition: mesh.cpp:6241
void Mult(const Table &A, const Table &B, Table &C)
C = A * B (as boolean matrices)
Definition: table.cpp:476
Pointer to an Operator of a specified type.
Definition: handle.hpp:33
void Copy(Array &copy) const
Create a copy of the current array.
Definition: array.hpp:812
virtual void Save(std::ostream &out) const
Definition: pgridfunc.cpp:485
Abstract parallel finite element space.
Definition: pfespace.hpp:28
virtual void ProjectCoefficient(Coefficient &coeff)
Definition: pgridfunc.cpp:301
bool iterative_mode
If true, use the second argument of Mult() as an initial guess.
Definition: operator.hpp:504
int main(int argc, char *argv[])
Definition: ex1.cpp:62
double * GetData() const
Return a pointer to the beginning of the Vector data.
Definition: vector.hpp:166
HypreParMatrix * ParallelAssemble()
Returns the matrix assembled on the true dofs, i.e. P^t A P.
int GetNRanks() const
Definition: pmesh.hpp:231
The BoomerAMG solver in hypre.
Definition: hypre.hpp:951
void SetPrintLevel(int print_lvl)
Definition: solvers.cpp:70
virtual void SetOperator(const Operator &op)
Also calls SetOperator for the preconditioner if there is one.
Definition: solvers.cpp:97
Mesh * GetMesh() const
Returns the mesh.
Definition: fespace.hpp:295
void UniformRefinement(int i, const DSTable &, int *, int *, int *)
Definition: mesh.cpp:7982
void SetPrintLevel(int print_level)
Definition: hypre.hpp:992
void SetMaxIter(int max_it)
Definition: solvers.hpp:65
T Max() const
Find the maximal element in the array, using the comparison operator &lt; for class T.
Definition: array.cpp:68
HYPRE_Int GlobalTrueVSize() const
Definition: pfespace.hpp:251
void Assemble(int skip_zeros=1)
Assemble the local matrix.
Newton&#39;s method for solving F(x)=b for a given operator F.
Definition: solvers.hpp:300
virtual void Print(std::ostream &out=mfem::out) const
Definition: pmesh.cpp:4102
void SetElasticityOptions(ParFiniteElementSpace *fespace)
Definition: hypre.cpp:2989
int Dimension() const
Definition: mesh.hpp:744
void PrintUsage(std::ostream &out) const
Definition: optparser.cpp:434
int SpaceDimension() const
Definition: mesh.hpp:745
void SetAbsTol(double atol)
Definition: solvers.hpp:64
Array< int > bdr_attributes
A list of all unique boundary attributes used by the Mesh.
Definition: mesh.hpp:189
int GetMyRank() const
Definition: pmesh.hpp:232
void SetRelTol(double rtol)
Definition: solvers.hpp:63
MPI_Comm GetComm() const
Definition: pmesh.hpp:230
Base class Coefficient that may optionally depend on time.
Definition: coefficient.hpp:31
void subtract(const Vector &x, const Vector &y, Vector &z)
Definition: vector.cpp:385
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)
Definition: optparser.hpp:76
void PartialSum()
Partial Sum.
Definition: array.cpp:103
GMRES method.
Definition: solvers.hpp:184
void Distribute(const Vector *tv)
Definition: pgridfunc.cpp:132
A class representing a general parallel block nonlinear operator defined on the Cartesian product of ...
double a
Definition: lissajous.cpp:41
virtual void Finalize(int skip_zeros=1)
Finalizes the matrix initialization.
int dim
Definition: ex24.cpp:43
void AddDomainIntegrator(BilinearFormIntegrator *bfi)
Adds new Domain Integrator. Assumes ownership of bfi.
void PrintOptions(std::ostream &out) const
Definition: optparser.cpp:304
HypreParVector * GetTrueDofs() const
Returns the true dofs in a new HypreParVector.
Definition: pgridfunc.cpp:143
Class for parallel bilinear form.
int open(const char hostname[], int port)
virtual void SetOperator(const Operator &op)
Also calls SetOperator for the preconditioner if there is one.
Definition: solvers.hpp:166
Vector data type.
Definition: vector.hpp:48
virtual void SetPreconditioner(Solver &pr)
This should be called before SetOperator.
Definition: solvers.cpp:91
Arbitrary order H1-conforming (continuous) finite elements.
Definition: fe_coll.hpp:83
Base class for solvers.
Definition: operator.hpp:500
Class for parallel grid function.
Definition: pgridfunc.hpp:32
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Definition: globals.hpp:66
Abstract operator.
Definition: operator.hpp:24
A class to handle Block systems in a matrix-free implementation.
Class for parallel meshes.
Definition: pmesh.hpp:32
Vector & GetBlock(int i)
Get the i-th vector in the block.
Definition: blockvector.hpp:84
bool Good() const
Definition: optparser.hpp:125