MFEM  v4.6.0 Finite element discretization library
ex34p.cpp
Go to the documentation of this file.
1 // MFEM Example 34 - Parallel Version
2 //
3 // Compile with: make ex34p
4 //
5 // Sample runs: mpirun -np 4 ex34p -o 2
6 // mpirun -np 4 ex34p -o 2 -hex -pa
7 //
8 // Device sample runs:
9 // mpirun -np 4 ex34p -o 2 -hex -pa -d cuda
10 // mpirun -np 4 ex34p -o 2 -no-pa -d cuda
11 //
12 // Description: This example code solves a simple magnetostatic problem
13 // curl curl A = J where the current density J is computed on a
14 // subset of the domain as J = -sigma grad phi. We discretize the
15 // vector potential with Nedelec finite elements, the scalar
16 // potential with Lagrange finite elements, and the current
17 // density with Raviart-Thomas finite elements.
18 //
19 // The example demonstrates the use of a SubMesh to compute the
20 // scalar potential and its associated current density which is
21 // then transferred to the original mesh and used as a source
22 // function.
23 //
24 // Note that this example takes certain liberties with the
25 // current density which is not necessarily divergence free
26 // as it should be. This was done to focus on the use of the
27 // SubMesh to transfer information between a full mesh and a
28 // sub-domain. A more rigorous implementation might employ an
29 // H(div) saddle point solver to obtain a divergence free J on
30 // the SubMesh. It would then also need to ensure that the r.h.s.
31 // of curl curl A = J does in fact lie in the range of the weak
32 // curl operator by performing a divergence cleaning procedure
33 // before the solve. After divergence cleaning the delta
34 // parameter would probably not be needed.
35 //
36 // This example is designed to make use of a specific mesh which
37 // has a known configuration of elements and boundary attributes.
38 // Other meshes could be used but extra care would be required to
39 // properly define the SubMesh and the necessary boundaries.
40 //
41 // We recommend viewing examples 1 and 3 before viewing this
42 // example.
43
44 #include "mfem.hpp"
45 #include <fstream>
46 #include <iostream>
47
48 using namespace std;
49 using namespace mfem;
50
51 void ComputeCurrentDensityOnSubMesh(int order,
52  const Array<int> &phi0_attr,
53  const Array<int> &phi1_attr,
54  const Array<int> &jn_zero_attr,
55  ParGridFunction &j_cond);
56
57 int main(int argc, char *argv[])
58 {
59  // 1. Initialize MPI and HYPRE.
60  Mpi::Init(argc, argv);
61  int num_procs = Mpi::WorldSize();
62  int myid = Mpi::WorldRank();
63  Hypre::Init();
64
65  // 2. Parse command-line options.
66  const char *mesh_file = "../data/fichera-mixed.mesh";
67  Array<int> cond_attr;
68  Array<int> submesh_elems;
69  Array<int> sym_plane_attr;
70  Array<int> phi0_attr;
71  Array<int> phi1_attr;
72  Array<int> jn_zero_attr;
73  int ser_ref_levels = 1;
74  int par_ref_levels = 1;
75  int order = 1;
76  double delta_const = 1e-6;
77  bool mixed = true;
78  bool static_cond = false;
79  bool pa = false;
80  const char *device_config = "cpu";
81  bool visualization = true;
82 #ifdef MFEM_USE_AMGX
83  bool useAmgX = false;
84 #endif
85
86  OptionsParser args(argc, argv);
88  "Number of times to refine the mesh uniformly in serial.");
90  "Number of times to refine the mesh uniformly in parallel.");
92  "Finite element order (polynomial degree).");
94  "Magnetic Conductivity");
96  "--no-static-condensation", "Enable static condensation.");
98  "--hex-mesh", "Mixed mesh of hexahedral mesh.");
100  "--no-partial-assembly", "Enable Partial Assembly.");
102  "Device configuration string, see Device::Configure().");
104  "--no-visualization",
105  "Enable or disable GLVis visualization.");
106 #ifdef MFEM_USE_AMGX
108  "--no-useAmgX",
109  "Enable or disable AmgX in MatrixFreeAMS.");
110 #endif
111
112  args.Parse();
113  if (!args.Good())
114  {
115  if (myid == 0)
116  {
117  args.PrintUsage(cout);
118  }
119  return 1;
120  }
121  if (myid == 0)
122  {
123  args.PrintOptions(cout);
124  }
125
126  if (!mixed || pa)
127  {
128  mesh_file = "../data/fichera.mesh";
129  }
130
131  if (submesh_elems.Size() == 0)
132  {
133  if (strcmp(mesh_file, "../data/fichera-mixed.mesh") == 0)
134  {
135  submesh_elems.SetSize(5);
136  submesh_elems[0] = 0;
137  submesh_elems[1] = 2;
138  submesh_elems[2] = 3;
139  submesh_elems[3] = 4;
140  submesh_elems[4] = 9;
141  }
142  else if (strcmp(mesh_file, "../data/fichera.mesh") == 0)
143  {
144  submesh_elems.SetSize(7);
145  submesh_elems[0] = 10;
146  submesh_elems[1] = 14;
147  submesh_elems[2] = 34;
148  submesh_elems[3] = 36;
149  submesh_elems[4] = 37;
150  submesh_elems[5] = 38;
151  submesh_elems[6] = 39;
152  }
153  }
154  if (sym_plane_attr.Size() == 0)
155  {
156  if (strcmp(mesh_file, "../data/fichera-mixed.mesh") == 0 ||
157  strcmp(mesh_file, "../data/fichera.mesh") == 0)
158  {
159  sym_plane_attr.SetSize(8);
160  sym_plane_attr[0] = 9;
161  sym_plane_attr[1] = 10;
162  sym_plane_attr[2] = 11;
163  sym_plane_attr[3] = 12;
164  sym_plane_attr[4] = 13;
165  sym_plane_attr[5] = 14;
166  sym_plane_attr[6] = 15;
167  sym_plane_attr[7] = 16;
168  }
169  }
170  if (phi0_attr.Size() == 0)
171  {
172  if (strcmp(mesh_file, "../data/fichera-mixed.mesh") == 0 ||
173  strcmp(mesh_file, "../data/fichera.mesh") == 0)
174  {
175  phi0_attr.Append(2);
176  }
177  }
178  if (phi1_attr.Size() == 0)
179  {
180  if (strcmp(mesh_file, "../data/fichera-mixed.mesh") == 0 ||
181  strcmp(mesh_file, "../data/fichera.mesh") == 0)
182  {
183  phi1_attr.Append(23);
184  }
185  }
186  if (jn_zero_attr.Size() == 0)
187  {
188  if (strcmp(mesh_file, "../data/fichera-mixed.mesh") == 0 ||
189  strcmp(mesh_file, "../data/fichera.mesh") == 0)
190  {
191  jn_zero_attr.Append(25);
192  }
193  for (int i=0; i<sym_plane_attr.Size(); i++)
194  {
195  jn_zero_attr.Append(sym_plane_attr[i]);
196  }
197  }
198
199  // 3. Enable hardware devices such as GPUs, and programming models such as
200  // CUDA, OCCA, RAJA and OpenMP based on command line options.
201  Device device(device_config);
202  if (myid == 0) { device.Print(); }
203
204  // 4. Read the (serial) mesh from the given mesh file on all processors. We
205  // can handle triangular, quadrilateral, tetrahedral, hexahedral, surface
206  // and volume meshes with the same code.
207  Mesh *mesh = new Mesh(mesh_file, 1, 1);
208  int dim = mesh->Dimension();
209
210  if (!mixed || pa)
211  {
212  mesh->UniformRefinement();
213
214  if (ser_ref_levels > 0)
215  {
216  ser_ref_levels--;
217  }
218  else
219  {
220  par_ref_levels--;
221  }
222  }
223
224  int submesh_attr = -1;
225  if (cond_attr.Size() == 0 && submesh_elems.Size() > 0)
226  {
227  int max_attr = mesh->attributes.Max();
228  submesh_attr = max_attr + 1;
229
230  for (int i=0; i<submesh_elems.Size(); i++)
231  {
232  mesh->SetAttribute(submesh_elems[i], submesh_attr);
233  }
234  mesh->SetAttributes();
235
236  if (cond_attr.Size() == 0)
237  {
238  cond_attr.Append(submesh_attr);
239  }
240  }
241
242  // 5. Refine the serial mesh on all processors to increase the resolution. In
243  // this example we do 'ref_levels' of uniform refinement.
244  {
245  int ref_levels = ser_ref_levels;
246  for (int l = 0; l < ref_levels; l++)
247  {
248  mesh->UniformRefinement();
249  }
250  }
251
252  // 6. Define a parallel mesh by a partitioning of the serial mesh. Refine
253  // this mesh further in parallel to increase the resolution. Once the
254  // parallel mesh is defined, the serial mesh can be deleted.
255  ParMesh pmesh(MPI_COMM_WORLD, *mesh);
256  delete mesh;
257  {
258  for (int l = 0; l < par_ref_levels; l++)
259  {
260  pmesh.UniformRefinement();
261  }
262  }
263
264  // 6b. Extract a submesh covering a portion of the domain
265  ParSubMesh pmesh_cond(ParSubMesh::CreateFromDomain(pmesh, cond_attr));
266
267  // 7. Define a suitable finite element space on the SubMesh and compute
268  // the current density as an H(div) field.
269  RT_FECollection fec_cond_rt(order - 1, dim);
270  ParFiniteElementSpace fes_cond_rt(&pmesh_cond, &fec_cond_rt);
271  ParGridFunction j_cond(&fes_cond_rt);
272
273  ComputeCurrentDensityOnSubMesh(order, phi0_attr, phi1_attr, jn_zero_attr,
274  j_cond);
275
276  // 7a. Save the SubMesh and associated current density in parallel. This
277  // output can be viewed later using GLVis:
278  // "glvis -np <np> -m cond_mesh -g cond_j"
279  {
280  ostringstream mesh_name, cond_name;
281  mesh_name << "cond_mesh." << setfill('0') << setw(6) << myid;
282  cond_name << "cond_j." << setfill('0') << setw(6) << myid;
283
284  ofstream mesh_ofs(mesh_name.str().c_str());
285  mesh_ofs.precision(8);
286  pmesh_cond.Print(mesh_ofs);
287
288  ofstream cond_ofs(cond_name.str().c_str());
289  cond_ofs.precision(8);
290  j_cond.Save(cond_ofs);
291  }
292  // 7b. Send the current density, computed on the SubMesh, to a GLVis server.
293  if (visualization)
294  {
295  char vishost[] = "localhost";
296  int visport = 19916;
297  socketstream port_sock(vishost, visport);
298  port_sock << "parallel " << num_procs << " " << myid << "\n";
299  port_sock.precision(8);
300  port_sock << "solution\n" << pmesh_cond << j_cond
301  << "window_title 'Conductor J'"
302  << "window_geometry 400 0 400 350" << flush;
303  }
304
305  // 8. Define a parallel finite element space on the full mesh. Here we use
306  // the H(curl) finite elements for the vector potential and H(div) for the
307  // current density.
308  ND_FECollection fec_nd(order, dim);
309  RT_FECollection fec_rt(order - 1, dim);
310  ParFiniteElementSpace fespace_nd(&pmesh, &fec_nd);
311  ParFiniteElementSpace fespace_rt(&pmesh, &fec_rt);
312
313  ParGridFunction j_full(&fespace_rt);
314  j_full = 0.0;
315  pmesh_cond.Transfer(j_cond, j_full);
316
317  // 8a. Send the transferred current density to a GLVis server.
318  if (visualization)
319  {
320  char vishost[] = "localhost";
321  int visport = 19916;
322  socketstream sol_sock(vishost, visport);
323  sol_sock << "parallel " << num_procs << " " << myid << "\n";
324  sol_sock.precision(8);
325  sol_sock << "solution\n" << pmesh << j_full
326  << "window_title 'J Full'"
327  << "window_geometry 400 430 400 350" << flush;
328  }
329
330  // 9. Determine the list of true (i.e. parallel conforming) essential
331  // boundary dofs. In this example, the boundary conditions are defined
332  // by marking all the boundary attributes except for those on a symmetry
333  // plane as essential (Dirichlet) and converting them to a list of
334  // true dofs.
335  Array<int> ess_tdof_list;
336  Array<int> ess_bdr;
337  if (pmesh.bdr_attributes.Size())
338  {
339  ess_bdr.SetSize(pmesh.bdr_attributes.Max());
340  ess_bdr = 1;
341  for (int i=0; i<sym_plane_attr.Size(); i++)
342  {
343  ess_bdr[sym_plane_attr[i]-1] = 0;
344  }
345  fespace_nd.GetEssentialTrueDofs(ess_bdr, ess_tdof_list);
346  }
347
348  // 10. Set up the parallel linear form b(.) which corresponds to the
349  // right-hand side of the FEM linear system, which in this case is
350  // (J,W_i) where J is given by the function H(div) field transferred
351  // from the SubMesh and W_i are the basis functions in the finite
352  // element fespace.
353  VectorGridFunctionCoefficient jCoef(&j_full);
354  ParLinearForm b(&fespace_nd);
356  b.Assemble();
357
358  // 11. Define the solution vector x as a parallel finite element grid
359  // function corresponding to fespace. Initialize x to zero.
360  ParGridFunction x(&fespace_nd);
361  x = 0.0;
362
363  // 12. Set up the parallel bilinear form corresponding to the EM diffusion
364  // operator curl muinv curl + delta I, by adding the curl-curl and the
365  // mass domain integrators. For standard magnetostatics equations choose
366  // delta << 1. Larger values of delta should make the linear system
367  // easier to solve at the expense of resembling a diffusive quasistatic
368  // magnetic field. A reasonable balance must be found whenever the mesh
369  // or problem setup is altered.
370  ConstantCoefficient muinv(1.0);
371  ConstantCoefficient delta(delta_const);
372  ParBilinearForm a(&fespace_nd);
373  if (pa) { a.SetAssemblyLevel(AssemblyLevel::PARTIAL); }
376
377  // 13. Assemble the parallel bilinear form and the corresponding linear
378  // system, applying any necessary transformations such as: parallel
379  // assembly, eliminating boundary conditions, applying conforming
380  // constraints for non-conforming AMR, static condensation, etc.
381  if (static_cond) { a.EnableStaticCondensation(); }
382  a.Assemble();
383
384  OperatorPtr A;
385  Vector B, X;
386  a.FormLinearSystem(ess_tdof_list, x, b, A, X, B);
387
388  if (myid == 0)
389  {
390  cout << "\nSolving for magnetic vector potential "
391  << "using CG with AMS" << endl;
392  }
393
394  // 14. Solve the system AX=B using PCG with an AMS preconditioner.
395  if (pa)
396  {
397 #ifdef MFEM_USE_AMGX
398  MatrixFreeAMS ams(a, *A, fespace_nd, &muinv, &delta, NULL, ess_bdr,
399  useAmgX);
400 #else
401  MatrixFreeAMS ams(a, *A, fespace_nd, &muinv, &delta, NULL, ess_bdr);
402 #endif
403  CGSolver cg(MPI_COMM_WORLD);
404  cg.SetRelTol(1e-12);
405  cg.SetMaxIter(1000);
406  cg.SetPrintLevel(1);
407  cg.SetOperator(*A);
408  cg.SetPreconditioner(ams);
409  cg.Mult(B, X);
410  }
411  else
412  {
413  if (myid == 0)
414  {
415  cout << "Size of linear system: "
416  << A.As<HypreParMatrix>()->GetGlobalNumRows() << endl;
417  }
418
419  ParFiniteElementSpace *prec_fespace =
420  (a.StaticCondensationIsEnabled() ? a.SCParFESpace() : &fespace_nd);
421  HypreAMS ams(*A.As<HypreParMatrix>(), prec_fespace);
422  HyprePCG pcg(*A.As<HypreParMatrix>());
423  pcg.SetTol(1e-12);
424  pcg.SetMaxIter(500);
425  pcg.SetPrintLevel(2);
426  pcg.SetPreconditioner(ams);
427  pcg.Mult(B, X);
428  }
429
430  // 15. Recover the parallel grid function corresponding to X. This is the
431  // local finite element solution on each processor.
432  a.RecoverFEMSolution(X, b, x);
433
434  // 16. Save the refined mesh and the solution in parallel. This output can
435  // be viewed later using GLVis: "glvis -np <np> -m mesh -g sol".
436  {
437  ostringstream mesh_name, sol_name;
438  mesh_name << "mesh." << setfill('0') << setw(6) << myid;
439  sol_name << "sol." << setfill('0') << setw(6) << myid;
440
441  ofstream mesh_ofs(mesh_name.str().c_str());
442  mesh_ofs.precision(8);
443  pmesh.Print(mesh_ofs);
444
445  ofstream sol_ofs(sol_name.str().c_str());
446  sol_ofs.precision(8);
447  x.Save(sol_ofs);
448  }
449
450  // 17. Send the solution by socket to a GLVis server.
451  if (visualization)
452  {
453  char vishost[] = "localhost";
454  int visport = 19916;
455  socketstream sol_sock(vishost, visport);
456  sol_sock << "parallel " << num_procs << " " << myid << "\n";
457  sol_sock.precision(8);
458  sol_sock << "solution\n" << pmesh << x
459  << "window_title 'Vector Potential'"
460  << "window_geometry 800 0 400 350" << flush;
461  }
462
463  // 18. Compute the magnetic flux as the curl of the solution
464  ParDiscreteLinearOperator curl(&fespace_nd, &fespace_rt);
466  curl.Assemble();
467  curl.Finalize();
468
469  ParGridFunction dx(&fespace_rt);
470  curl.Mult(x, dx);
471
472  // 19. Save the curl of the solution in parallel. This output can
473  // be viewed later using GLVis: "glvis -np <np> -m mesh -g dsol".
474  {
475  ostringstream dsol_name;
476  dsol_name << "dsol." << setfill('0') << setw(6) << myid;
477
478  ofstream dsol_ofs(dsol_name.str().c_str());
479  dsol_ofs.precision(8);
480  dx.Save(dsol_ofs);
481  }
482
483  // 20. Send the curl of the solution by socket to a GLVis server.
484  if (visualization)
485  {
486  char vishost[] = "localhost";
487  int visport = 19916;
488  socketstream sol_sock(vishost, visport);
489  sol_sock << "parallel " << num_procs << " " << myid << "\n";
490  sol_sock.precision(8);
491  sol_sock << "solution\n" << pmesh << dx
492  << "window_title 'Magnetic Flux'"
493  << "window_geometry 1200 0 400 350" << flush;
494  }
495
496  // 21. Clean exit
497  return 0;
498 }
499
501  const Array<int> &phi0_attr,
502  const Array<int> &phi1_attr,
503  const Array<int> &jn_zero_attr,
504  ParGridFunction &j_cond)
505 {
506  // Extract the finite element space and mesh on which j_cond is defined
507  ParFiniteElementSpace &fes_cond_rt = *j_cond.ParFESpace();
508  ParMesh &pmesh_cond = *fes_cond_rt.GetParMesh();
509  int myid = fes_cond_rt.GetMyRank();
510  int dim = pmesh_cond.Dimension();
511
512  // Define a parallel finite element space on the SubMesh. Here we use the
513  // H1 finite elements for the electrostatic potential.
514  H1_FECollection fec_h1(order, dim);
515  ParFiniteElementSpace fes_cond_h1(&pmesh_cond, &fec_h1);
516
517  // Define the conductivity coefficient and the boundaries associated with the
518  // fixed potentials phi0 and phi1 which will drive the current.
519  ConstantCoefficient sigmaCoef(1.0);
520  Array<int> ess_bdr_phi(pmesh_cond.bdr_attributes.Max());
521  Array<int> ess_bdr_j(pmesh_cond.bdr_attributes.Max());
522  Array<int> ess_bdr_tdof_phi;
523  ess_bdr_phi = 0;
524  ess_bdr_j = 0;
525  for (int i=0; i<phi0_attr.Size(); i++)
526  {
527  ess_bdr_phi[phi0_attr[i]-1] = 1;
528  }
529  for (int i=0; i<phi1_attr.Size(); i++)
530  {
531  ess_bdr_phi[phi1_attr[i]-1] = 1;
532  }
533  for (int i=0; i<jn_zero_attr.Size(); i++)
534  {
535  ess_bdr_j[jn_zero_attr[i]-1] = 1;
536  }
537  fes_cond_h1.GetEssentialTrueDofs(ess_bdr_phi, ess_bdr_tdof_phi);
538
539  // Setup the bilinear form corresponding to -Div(sigma Grad phi)
540  ParBilinearForm a_h1(&fes_cond_h1);
542  a_h1.Assemble();
543
544  // Set the r.h.s. to zero
545  ParLinearForm b_h1(&fes_cond_h1);
546  b_h1 = 0.0;
547
548  // Setup the boundary conditions on phi
549  ConstantCoefficient one(1.0);
550  ConstantCoefficient zero(0.0);
551  ParGridFunction phi_h1(&fes_cond_h1);
552  phi_h1 = 0.0;
553
554  Array<int> bdr0(pmesh_cond.bdr_attributes.Max()); bdr0 = 0;
555  for (int i=0; i<phi0_attr.Size(); i++)
556  {
557  bdr0[phi0_attr[i]-1] = 1;
558  }
559  phi_h1.ProjectBdrCoefficient(zero, bdr0);
560
561  Array<int> bdr1(pmesh_cond.bdr_attributes.Max()); bdr1 = 0;
562  for (int i=0; i<phi1_attr.Size(); i++)
563  {
564  bdr1[phi1_attr[i]-1] = 1;
565  }
566  phi_h1.ProjectBdrCoefficient(one, bdr1);
567
568  // Solve the linear system using algebraic multigrid
569  {
570  if (myid == 0)
571  {
572  cout << "\nSolving for electric potential "
573  << "using CG with AMG" << endl;
574  }
575  OperatorPtr A;
576  Vector B, X;
577  a_h1.FormLinearSystem(ess_bdr_tdof_phi, phi_h1, b_h1, A, X, B);
578
579  HypreBoomerAMG prec;
580  CGSolver cg(MPI_COMM_WORLD);
581  cg.SetRelTol(1e-12);
582  cg.SetMaxIter(2000);
583  cg.SetPrintLevel(1);
584  cg.SetPreconditioner(prec);
585  cg.SetOperator(*A);
586  cg.Mult(B, X);
587  a_h1.RecoverFEMSolution(X, b_h1, phi_h1);
588  }
589  {
590  int num_procs = fes_cond_h1.GetNRanks();
591  char vishost[] = "localhost";
592  int visport = 19916;
593  socketstream port_sock(vishost, visport);
594  port_sock << "parallel " << num_procs << " " << myid << "\n";
595  port_sock.precision(8);
596  port_sock << "solution\n" << pmesh_cond << phi_h1
597  << "window_title 'Conductor Potential'"
598  << "window_geometry 0 0 400 350" << flush;
599  }
600
601  // Solve for the current density J = -sigma Grad phi with boundary conditions
602  // J.n = 0 on the walls of the conductor but not on the ports where phi=0 and
603  // phi=1.
604
605  // J will be computed in H(div) so we need an RT mass matrix
606  ParBilinearForm m_rt(&fes_cond_rt);
608  m_rt.Assemble();
609
610  // Assemble the (sigma Grad phi) operator
611  ParMixedBilinearForm d_h1(&fes_cond_h1, &fes_cond_rt);
613  d_h1.Assemble();
614
615  // Compute the r.h.s, b_rt = sigma E = -sigma Grad phi
616  ParLinearForm b_rt(&fes_cond_rt);
617  d_h1.Mult(phi_h1, b_rt);
618  b_rt *= -1.0;
619
620  // Apply the necessary boundary conditions and solve for J in H(div)
621  HYPRE_BigInt glb_size_rt = fes_cond_rt.GlobalTrueVSize();
622  if (myid == 0)
623  {
624  cout << "\nSolving for current density in H(Div) "
625  << "using diagonally scaled CG" << endl;
626  cout << "Size of linear system: "
627  << glb_size_rt << endl;
628  }
629  Array<int> ess_bdr_tdof_rt;
630  OperatorPtr M;
631  Vector B, X;
632
633  fes_cond_rt.GetEssentialTrueDofs(ess_bdr_j, ess_bdr_tdof_rt);
634
635  j_cond = 0.0;
636  m_rt.FormLinearSystem(ess_bdr_tdof_rt, j_cond, b_rt, M, X, B);
637
638  HypreDiagScale prec;
639
640  CGSolver cg(MPI_COMM_WORLD);
641  cg.SetRelTol(1e-12);
642  cg.SetMaxIter(2000);
643  cg.SetPrintLevel(1);
644  cg.SetPreconditioner(prec);
645  cg.SetOperator(*M);
646  cg.Mult(B, X);
647  m_rt.RecoverFEMSolution(X, b_rt, j_cond);
648 }
void SetTol(double tol)
Definition: hypre.cpp:3990
virtual void GetEssentialTrueDofs(const Array< int > &bdr_attr_is_ess, Array< int > &ess_tdof_list, int component=-1)
Definition: pfespace.cpp:1031
int visport
Definition: solvers.hpp:493
The Auxiliary-space Maxwell Solver in hypre.
Definition: hypre.hpp:1743
A coefficient that is constant across space and time.
Definition: coefficient.hpp:84
virtual void FormLinearSystem(const Array< int > &ess_tdof_list, Vector &x, Vector &b, OperatorHandle &A, Vector &X, Vector &B, int copy_interior=0)
Form the linear system A X = B, corresponding to this bilinear form and the linear form b(...
void PrintOptions(std::ostream &out) const
Print the options.
Definition: optparser.cpp:331
int Dimension() const
Dimension of the reference space used within the elements.
Definition: mesh.hpp:1020
void PrintUsage(std::ostream &out) const
Print the usage message.
Definition: optparser.cpp:462
virtual void Mult(const HypreParVector &b, HypreParVector &x) const
Solve Ax=b with hypre&#39;s PCG.
Definition: hypre.cpp:4038
Integrator for (curl u, curl v) for Nedelec elements.
Pointer to an Operator of a specified type.
Definition: handle.hpp:33
virtual void Mult(const Vector &b, Vector &x) const
Operator application: y=A(x).
Definition: solvers.cpp:718
T Max() const
Find the maximal element in the array, using the comparison operator < for class T.
Definition: array.cpp:68
int main(int argc, char *argv[])
Definition: ex34p.cpp:57
void Print(std::ostream &out=mfem::out)
Print the configuration of the MFEM virtual device object.
Definition: device.cpp:279
bool Good() const
Return true if the command line options were parsed successfully.
Definition: optparser.hpp:159
Abstract parallel finite element space.
Definition: pfespace.hpp:28
STL namespace.
virtual void Finalize(int skip_zeros=1)
Finalizes the matrix initialization.
virtual void Mult(const Vector &x, Vector &y) const
Matrix multiplication: .
void ComputeCurrentDensityOnSubMesh(int order, const Array< int > &phi0_attr, const Array< int > &phi1_attr, const Array< int > &jn_zero_attr, ParGridFunction &j_cond)
Definition: ex34p.cpp:500
void SetPrintLevel(int print_lvl)
Definition: hypre.cpp:4010
Subdomain representation of a topological parent in another ParMesh.
Definition: psubmesh.hpp:51
Adds a domain interpolator. Assumes ownership of di.
Adds a domain integrator. Assumes ownership of bfi.
The BoomerAMG solver in hypre.
Definition: hypre.hpp:1590
Class for parallel linear form.
Definition: plinearform.hpp:26
An auxiliary Maxwell solver for a high-order curl-curl system without high-order assembly.
Definition: auxiliary.hpp:167
virtual void SetPrintLevel(int print_lvl)
Legacy method to set the level of verbosity of the solver output.
Definition: solvers.cpp:71
void Parse()
Parse the command-line options. Note that this function expects all the options provided through the ...
Definition: optparser.cpp:151
char vishost[]
virtual void SetAttributes()
Determine the sets of unique attribute values in domain and boundary elements.
Definition: mesh.cpp:1572
int Append(const T &el)
Append element &#39;el&#39; to array, resize if necessary.
Definition: array.hpp:759
Jacobi preconditioner in hypre.
Definition: hypre.hpp:1421
double b
Definition: lissajous.cpp:42
void UniformRefinement(int i, const DSTable &, int *, int *, int *)
Definition: mesh.cpp:10232
ParMesh * GetParMesh() const
Definition: pfespace.hpp:273
void SetMaxIter(int max_it)
Definition: solvers.hpp:201
double delta
Definition: lissajous.cpp:43
void Assemble(int skip_zeros=1)
void Assemble(int skip_zeros=1)
Assemble the local matrix.
HYPRE_BigInt GlobalTrueVSize() const
Definition: pfespace.hpp:281
void SetMaxIter(int max_iter)
Definition: hypre.cpp:4000
Arbitrary order H(div)-conforming Raviart-Thomas finite elements.
Definition: fe_coll.hpp:380
Array< int > bdr_attributes
A list of all unique boundary attributes used by the Mesh.
Definition: mesh.hpp:275
void SetRelTol(double rtol)
Definition: solvers.hpp:199
PCG solver in hypre.
Definition: hypre.hpp:1215
void ProjectBdrCoefficient(Coefficient *coeff[], VectorCoefficient *vcoeff, Array< int > &attr)
Definition: pgridfunc.cpp:669
ParFiniteElementSpace * ParFESpace() const
Definition: pgridfunc.hpp:113
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 &#39;var&#39; to receive the value. Enable/disable tags are used to set the bool...
Definition: optparser.hpp:82
void SetSize(int nsize)
Change the logical size of the array, keep existing entries.
Definition: array.hpp:687
HYPRE_Int HYPRE_BigInt
virtual void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x)
virtual void Save(std::ostream &out) const
Definition: pgridfunc.cpp:909
double a
Definition: lissajous.cpp:41
Class for parallel bilinear form using different test and trial FE spaces.
OpType * As() const
Return the Operator pointer statically cast to a specified OpType. Similar to the method Get()...
Definition: handle.hpp:104
int dim
Definition: ex24.cpp:53
void SetPreconditioner(HypreSolver &precond)
Set the hypre solver to be used as a preconditioner.
Definition: hypre.cpp:4015
Adds new Domain Integrator. Assumes ownership of bfi.
Class for parallel bilinear form.
int Size() const
Return the logical size of the array.
Definition: array.hpp:141
for VectorFiniteElements (Nedelec, Raviart-Thomas)
Definition: lininteg.hpp:346
virtual void SetOperator(const Operator &op)
Also calls SetOperator for the preconditioner if there is one.
Definition: solvers.hpp:507
void SetAttribute(int i, int attr)
Set the attribute of element i.
Definition: mesh.hpp:1193
Arbitrary order H(curl)-conforming Nedelec finite elements.
Definition: fe_coll.hpp:454
Vector data type.
Definition: vector.hpp:58
Vector coefficient defined by a vector GridFunction.
virtual void SetPreconditioner(Solver &pr)
This should be called before SetOperator.
Definition: solvers.cpp:173
Arbitrary order H1-conforming (continuous) finite elements.
Definition: fe_coll.hpp:259
void Print(std::ostream &out=mfem::out) const override
Definition: pmesh.cpp:4825
Class for parallel grid function.
Definition: pgridfunc.hpp:32
static void Transfer(const ParGridFunction &src, ParGridFunction &dst)
Transfer the dofs of a ParGridFunction.
Definition: psubmesh.cpp:932
The MFEM Device class abstracts hardware devices such as GPUs, as well as programming models such as ...
Definition: device.hpp:121
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:343
virtual void Assemble(int skip_zeros=1)
Construct the internal matrix representation of the discrete linear operator.
Class for parallel meshes.
Definition: pmesh.hpp:32
Array< int > attributes
A list of all unique element attributes used by the Mesh.
Definition: mesh.hpp:273