MFEM v4.7.0
Finite element discretization library
Loading...
Searching...
No Matches
ex31p.cpp
Go to the documentation of this file.
1// MFEM Example 31 - Parallel Version
2//
3// Compile with: make ex31p
4//
5// Sample runs: mpirun -np 4 ex31p -m ../data/hexagon.mesh -o 2
6// mpirun -np 4 ex31p -m ../data/star.mesh
7// mpirun -np 4 ex31p -m ../data/square-disc.mesh -o 2
8// mpirun -np 4 ex31p -m ../data/fichera.mesh -o 3 -rs 1 -rp 0
9// mpirun -np 4 ex31p -m ../data/square-disc-nurbs.mesh -o 3
10// mpirun -np 4 ex31p -m ../data/amr-quad.mesh -o 2 -rs 1
11// mpirun -np 4 ex31p -m ../data/amr-hex.mesh -rs 1
12//
13// Description: This example code solves a simple electromagnetic diffusion
14// problem corresponding to the second order definite Maxwell
15// equation curl curl E + sigma E = f with boundary condition
16// E x n = <given tangential field>. In this example sigma is an
17// anisotropic 3x3 tensor. Here, we use a given exact solution E
18// and compute the corresponding r.h.s. f. We discretize with
19// Nedelec finite elements in 1D, 2D, or 3D.
20//
21// The example demonstrates the use of restricted H(curl) finite
22// element spaces with the curl-curl and the (vector finite
23// element) mass bilinear form, as well as the computation of
24// discretization error when the exact solution is known. These
25// restricted spaces allow the solution of 1D or 2D
26// electromagnetic problems which involve 3D field vectors. Such
27// problems arise in plasma physics and crystallography.
28//
29// We recommend viewing example 3 before viewing this example.
30
31#include "mfem.hpp"
32#include <fstream>
33#include <iostream>
34
35using namespace std;
36using namespace mfem;
37
38// Exact solution, E, and r.h.s., f. See below for implementation.
39void E_exact(const Vector &, Vector &);
40void CurlE_exact(const Vector &, Vector &);
41void f_exact(const Vector &, Vector &);
43int dim;
44
45int main(int argc, char *argv[])
46{
47 // 1. Initialize MPI.
48 Mpi::Init(argc, argv);
49 int num_procs = Mpi::WorldSize();
50 int myid = Mpi::WorldRank();
52
53 // 2. Parse command-line options.
54 const char *mesh_file = "../data/inline-quad.mesh";
55 int ser_ref_levels = 2;
56 int par_ref_levels = 1;
57 int order = 1;
58 bool use_ams = true;
59 bool visualization = 1;
60
61 OptionsParser args(argc, argv);
62 args.AddOption(&mesh_file, "-m", "--mesh",
63 "Mesh file to use.");
64 args.AddOption(&ser_ref_levels, "-rs", "--refine-serial",
65 "Number of times to refine the mesh uniformly in serial.");
66 args.AddOption(&par_ref_levels, "-rp", "--refine-parallel",
67 "Number of times to refine the mesh uniformly in parallel.");
68 args.AddOption(&order, "-o", "--order",
69 "Finite element order (polynomial degree).");
70 args.AddOption(&freq, "-f", "--frequency", "Set the frequency for the exact"
71 " solution.");
72 args.AddOption(&use_ams, "-ams", "--hypre-ams", "-slu",
73 "--superlu", "Use AMS or SuperLU solver.");
74 args.AddOption(&visualization, "-vis", "--visualization", "-no-vis",
75 "--no-visualization",
76 "Enable or disable GLVis visualization.");
77 args.ParseCheck();
78
79 kappa = freq * M_PI;
80
81 // 3. Read the (serial) mesh from the given mesh file on all processors. We
82 // can handle triangular, quadrilateral, tetrahedral, hexahedral, surface
83 // and volume meshes with the same code.
84 Mesh *mesh = new Mesh(mesh_file, 1, 1);
85 dim = mesh->Dimension();
86
87 // 4. Refine the serial mesh on all processors to increase the resolution. In
88 // this example we do 'ref_levels' of uniform refinement (2 by default, or
89 // specified on the command line with -rs).
90 for (int lev = 0; lev < ser_ref_levels; lev++)
91 {
92 mesh->UniformRefinement();
93 }
94
95 // 5. Define a parallel mesh by a partitioning of the serial mesh. Refine
96 // this mesh further in parallel to increase the resolution (1 time by
97 // default, or specified on the command line with -rp). Once the parallel
98 // mesh is defined, the serial mesh can be deleted.
99 ParMesh pmesh(MPI_COMM_WORLD, *mesh);
100 delete mesh;
101 for (int lev = 0; lev < par_ref_levels; lev++)
102 {
103 pmesh.UniformRefinement();
104 }
105
106 // 6. Define a parallel finite element space on the parallel mesh. Here we
107 // use the Nedelec finite elements of the specified order.
108 FiniteElementCollection *fec = NULL;
109 if (dim == 1)
110 {
111 fec = new ND_R1D_FECollection(order, dim);
112 }
113 else if (dim == 2)
114 {
115 fec = new ND_R2D_FECollection(order, dim);
116 }
117 else
118 {
119 fec = new ND_FECollection(order, dim);
120 }
121 ParFiniteElementSpace fespace(&pmesh, fec);
122 HYPRE_Int size = fespace.GlobalTrueVSize();
123 if (Mpi::Root()) { cout << "Number of H(Curl) unknowns: " << size << endl; }
124
125 // 7. Determine the list of true (i.e. parallel conforming) essential
126 // boundary dofs. In this example, the boundary conditions are defined
127 // by marking all the boundary attributes from the mesh as essential
128 // (Dirichlet) and converting them to a list of true dofs.
129 Array<int> ess_tdof_list;
130 if (pmesh.bdr_attributes.Size())
131 {
132 Array<int> ess_bdr(pmesh.bdr_attributes.Max());
133 ess_bdr = 1;
134 fespace.GetEssentialTrueDofs(ess_bdr, ess_tdof_list);
135 }
136
137 // 8. Set up the parallel linear form b(.) which corresponds to the
138 // right-hand side of the FEM linear system, which in this case is
139 // (f,phi_i) where f is given by the function f_exact and phi_i are the
140 // basis functions in the finite element fespace.
142 ParLinearForm b(&fespace);
143 b.AddDomainIntegrator(new VectorFEDomainLFIntegrator(f));
144 b.Assemble();
145
146 // 9. Define the solution vector x as a parallel finite element grid function
147 // corresponding to fespace. Initialize x by projecting the exact
148 // solution. Note that only values from the boundary edges will be used
149 // when eliminating the non-homogeneous boundary condition to modify the
150 // r.h.s. vector b.
151 ParGridFunction sol(&fespace);
154 sol.ProjectCoefficient(E);
155
156 // 10. Set up the parallel bilinear form corresponding to the EM diffusion
157 // operator curl muinv curl + sigma I, by adding the curl-curl and the
158 // mass domain integrators.
159 DenseMatrix sigmaMat(3);
160 sigmaMat(0,0) = 2.0; sigmaMat(1,1) = 2.0; sigmaMat(2,2) = 2.0;
161 sigmaMat(0,2) = 0.0; sigmaMat(2,0) = 0.0;
162 sigmaMat(0,1) = M_SQRT1_2; sigmaMat(1,0) = M_SQRT1_2; // 1/sqrt(2) in cmath
163 sigmaMat(1,2) = M_SQRT1_2; sigmaMat(2,1) = M_SQRT1_2;
164
165 ConstantCoefficient muinv(1.0);
167 ParBilinearForm a(&fespace);
168 a.AddDomainIntegrator(new CurlCurlIntegrator(muinv));
169 a.AddDomainIntegrator(new VectorFEMassIntegrator(sigma));
170
171 // 11. Assemble the parallel bilinear form and the corresponding linear
172 // system, applying any necessary transformations such as: parallel
173 // assembly, eliminating boundary conditions, applying conforming
174 // constraints for non-conforming AMR, etc.
175 a.Assemble();
176
177 OperatorPtr A;
178 Vector B, X;
179
180 a.FormLinearSystem(ess_tdof_list, sol, b, A, X, B);
181
182 // 12. Solve the system AX=B using PCG with the AMS preconditioner from hypre
183 if (use_ams)
184 {
185 if (Mpi::Root())
186 {
187 cout << "Size of linear system: "
188 << A.As<HypreParMatrix>()->GetGlobalNumRows() << endl;
189 }
190
191 HypreAMS ams(*A.As<HypreParMatrix>(), &fespace);
192
193 HyprePCG pcg(*A.As<HypreParMatrix>());
194 pcg.SetTol(1e-12);
195 pcg.SetMaxIter(1000);
196 pcg.SetPrintLevel(2);
197 pcg.SetPreconditioner(ams);
198 pcg.Mult(B, X);
199 }
200 else
201#ifdef MFEM_USE_SUPERLU
202 {
203 if (Mpi::Root())
204 {
205 cout << "Size of linear system: "
206 << A.As<HypreParMatrix>()->GetGlobalNumRows() << endl;
207 }
208
209 SuperLURowLocMatrix A_SuperLU(*A.As<HypreParMatrix>());
210 SuperLUSolver AInv(MPI_COMM_WORLD);
211 AInv.SetOperator(A_SuperLU);
212 AInv.Mult(B,X);
213 }
214#else
215 {
216 if (Mpi::Root()) { cout << "No solvers available." << endl; }
217 return 1;
218 }
219#endif
220
221 // 13. Recover the parallel grid function corresponding to X. This is the
222 // local finite element solution on each processor.
223 a.RecoverFEMSolution(X, b, sol);
224
225 // 14. Compute and print the H(Curl) norm of the error.
226 {
227 real_t error = sol.ComputeHCurlError(&E, &CurlE);
228 if (Mpi::Root())
229 {
230 cout << "\n|| E_h - E ||_{H(Curl)} = " << error << '\n' << endl;
231 }
232 }
233
234
235 // 15. Save the refined mesh and the solution in parallel. This output can
236 // be viewed later using GLVis: "glvis -np <np> -m mesh -g sol".
237 {
238 ostringstream mesh_name, sol_name;
239 mesh_name << "mesh." << setfill('0') << setw(6) << myid;
240 sol_name << "sol." << setfill('0') << setw(6) << myid;
241
242 ofstream mesh_ofs(mesh_name.str().c_str());
243 mesh_ofs.precision(8);
244 pmesh.Print(mesh_ofs);
245
246 ofstream sol_ofs(sol_name.str().c_str());
247 sol_ofs.precision(8);
248 sol.Save(sol_ofs);
249 }
250
251 // 16. Send the solution by socket to a GLVis server.
252 if (visualization)
253 {
254 char vishost[] = "localhost";
255 int visport = 19916;
256
257 VectorGridFunctionCoefficient solCoef(&sol);
258 CurlGridFunctionCoefficient dsolCoef(&sol);
259
260 if (dim ==1)
261 {
265 socketstream dy_sock(vishost, visport);
266 socketstream dz_sock(vishost, visport);
267 x_sock.precision(8);
268 y_sock.precision(8);
269 z_sock.precision(8);
270 dy_sock.precision(8);
271 dz_sock.precision(8);
272
273 Vector xVec(3); xVec = 0.0; xVec(0) = 1;
274 Vector yVec(3); yVec = 0.0; yVec(1) = 1;
275 Vector zVec(3); zVec = 0.0; zVec(2) = 1;
276 VectorConstantCoefficient xVecCoef(xVec);
277 VectorConstantCoefficient yVecCoef(yVec);
278 VectorConstantCoefficient zVecCoef(zVec);
279
280 H1_FECollection fec_h1(order, dim);
281 L2_FECollection fec_l2(order-1, dim);
282
283 ParFiniteElementSpace fes_h1(&pmesh, &fec_h1);
284 ParFiniteElementSpace fes_l2(&pmesh, &fec_l2);
285
286 ParGridFunction xComp(&fes_l2);
287 ParGridFunction yComp(&fes_h1);
288 ParGridFunction zComp(&fes_h1);
289
290 ParGridFunction dyComp(&fes_l2);
291 ParGridFunction dzComp(&fes_l2);
292
293 InnerProductCoefficient xCoef(xVecCoef, solCoef);
294 InnerProductCoefficient yCoef(yVecCoef, solCoef);
295 InnerProductCoefficient zCoef(zVecCoef, solCoef);
296
297 xComp.ProjectCoefficient(xCoef);
298 yComp.ProjectCoefficient(yCoef);
299 zComp.ProjectCoefficient(zCoef);
300
301 x_sock << "parallel " << num_procs << " " << myid << "\n"
302 << "solution\n" << pmesh << xComp << flush
303 << "window_title 'X component'" << endl;
304 y_sock << "parallel " << num_procs << " " << myid << "\n"
305 << "solution\n" << pmesh << yComp << flush
306 << "window_geometry 403 0 400 350 "
307 << "window_title 'Y component'" << endl;
308 z_sock << "parallel " << num_procs << " " << myid << "\n"
309 << "solution\n" << pmesh << zComp << flush
310 << "window_geometry 806 0 400 350 "
311 << "window_title 'Z component'" << endl;
312
313 InnerProductCoefficient dyCoef(yVecCoef, dsolCoef);
314 InnerProductCoefficient dzCoef(zVecCoef, dsolCoef);
315
316 dyComp.ProjectCoefficient(dyCoef);
317 dzComp.ProjectCoefficient(dzCoef);
318
319 dy_sock << "parallel " << num_procs << " " << myid << "\n"
320 << "solution\n" << pmesh << dyComp << flush
321 << "window_geometry 403 375 400 350 "
322 << "window_title 'Y component of Curl'" << endl;
323 dz_sock << "parallel " << num_procs << " " << myid << "\n"
324 << "solution\n" << pmesh << dzComp << flush
325 << "window_geometry 806 375 400 350 "
326 << "window_title 'Z component of Curl'" << endl;
327 }
328 else if (dim == 2)
329 {
330 socketstream xy_sock(vishost, visport);
332 socketstream dxy_sock(vishost, visport);
333 socketstream dz_sock(vishost, visport);
334
335 DenseMatrix xyMat(2,3); xyMat = 0.0;
336 xyMat(0,0) = 1.0; xyMat(1,1) = 1.0;
337 MatrixConstantCoefficient xyMatCoef(xyMat);
338 Vector zVec(3); zVec = 0.0; zVec(2) = 1;
339 VectorConstantCoefficient zVecCoef(zVec);
340
341 MatrixVectorProductCoefficient xyCoef(xyMatCoef, solCoef);
342 InnerProductCoefficient zCoef(zVecCoef, solCoef);
343
344 H1_FECollection fec_h1(order, dim);
345 ND_FECollection fec_nd(order, dim);
346 RT_FECollection fec_rt(order-1, dim);
347 L2_FECollection fec_l2(order-1, dim);
348
349 ParFiniteElementSpace fes_h1(&pmesh, &fec_h1);
350 ParFiniteElementSpace fes_nd(&pmesh, &fec_nd);
351 ParFiniteElementSpace fes_rt(&pmesh, &fec_rt);
352 ParFiniteElementSpace fes_l2(&pmesh, &fec_l2);
353
354 ParGridFunction xyComp(&fes_nd);
355 ParGridFunction zComp(&fes_h1);
356
357 ParGridFunction dxyComp(&fes_rt);
358 ParGridFunction dzComp(&fes_l2);
359
360 xyComp.ProjectCoefficient(xyCoef);
361 zComp.ProjectCoefficient(zCoef);
362
363 xy_sock << "parallel " << num_procs << " " << myid << "\n";
364 xy_sock.precision(8);
365 xy_sock << "solution\n" << pmesh << xyComp
366 << "window_title 'XY components'\n" << flush;
367 z_sock << "parallel " << num_procs << " " << myid << "\n"
368 << "solution\n" << pmesh << zComp << flush
369 << "window_geometry 403 0 400 350 "
370 << "window_title 'Z component'" << endl;
371
372 MatrixVectorProductCoefficient dxyCoef(xyMatCoef, dsolCoef);
373 InnerProductCoefficient dzCoef(zVecCoef, dsolCoef);
374
375 dxyComp.ProjectCoefficient(dxyCoef);
376 dzComp.ProjectCoefficient(dzCoef);
377
378 dxy_sock << "parallel " << num_procs << " " << myid << "\n"
379 << "solution\n" << pmesh << dxyComp << flush
380 << "window_geometry 0 375 400 350 "
381 << "window_title 'XY components of Curl'" << endl;
382 dz_sock << "parallel " << num_procs << " " << myid << "\n"
383 << "solution\n" << pmesh << dzComp << flush
384 << "window_geometry 403 375 400 350 "
385 << "window_title 'Z component of Curl'" << endl;
386 }
387 else
388 {
389 socketstream sol_sock(vishost, visport);
390 socketstream dsol_sock(vishost, visport);
391
392 RT_FECollection fec_rt(order-1, dim);
393
394 ParFiniteElementSpace fes_rt(&pmesh, &fec_rt);
395
396 ParGridFunction dsol(&fes_rt);
397
398 dsol.ProjectCoefficient(dsolCoef);
399
400 sol_sock << "parallel " << num_procs << " " << myid << "\n";
401 sol_sock.precision(8);
402 sol_sock << "solution\n" << pmesh << sol
403 << "window_title 'Solution'" << flush << endl;
404 dsol_sock << "parallel " << num_procs << " " << myid << "\n"
405 << "solution\n" << pmesh << dsol << flush
406 << "window_geometry 0 375 400 350 "
407 << "window_title 'Curl of solution'" << endl;
408 }
409 }
410
411 // 17. Free the used memory.
412 delete fec;
413
414 return 0;
415}
416
417
418void E_exact(const Vector &x, Vector &E)
419{
420 if (dim == 1)
421 {
422 E(0) = 1.1 * sin(kappa * x(0) + 0.0 * M_PI);
423 E(1) = 1.2 * sin(kappa * x(0) + 0.4 * M_PI);
424 E(2) = 1.3 * sin(kappa * x(0) + 0.9 * M_PI);
425 }
426 else if (dim == 2)
427 {
428 E(0) = 1.1 * sin(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.0 * M_PI);
429 E(1) = 1.2 * sin(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.4 * M_PI);
430 E(2) = 1.3 * sin(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.9 * M_PI);
431 }
432 else
433 {
434 E(0) = 1.1 * sin(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.0 * M_PI);
435 E(1) = 1.2 * sin(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.4 * M_PI);
436 E(2) = 1.3 * sin(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.9 * M_PI);
437 E *= cos(kappa * x(2));
438 }
439}
440
441void CurlE_exact(const Vector &x, Vector &dE)
442{
443 if (dim == 1)
444 {
445 real_t c4 = cos(kappa * x(0) + 0.4 * M_PI);
446 real_t c9 = cos(kappa * x(0) + 0.9 * M_PI);
447
448 dE(0) = 0.0;
449 dE(1) = -1.3 * c9;
450 dE(2) = 1.2 * c4;
451 dE *= kappa;
452 }
453 else if (dim == 2)
454 {
455 real_t c0 = cos(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.0 * M_PI);
456 real_t c4 = cos(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.4 * M_PI);
457 real_t c9 = cos(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.9 * M_PI);
458
459 dE(0) = 1.3 * c9;
460 dE(1) = -1.3 * c9;
461 dE(2) = 1.2 * c4 - 1.1 * c0;
462 dE *= kappa * M_SQRT1_2;
463 }
464 else
465 {
466 real_t s0 = sin(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.0 * M_PI);
467 real_t c0 = cos(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.0 * M_PI);
468 real_t s4 = sin(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.4 * M_PI);
469 real_t c4 = cos(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.4 * M_PI);
470 real_t c9 = cos(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.9 * M_PI);
471 real_t sk = sin(kappa * x(2));
472 real_t ck = cos(kappa * x(2));
473
474 dE(0) = 1.2 * s4 * sk + 1.3 * M_SQRT1_2 * c9 * ck;
475 dE(1) = -1.1 * s0 * sk - 1.3 * M_SQRT1_2 * c9 * ck;
476 dE(2) = -M_SQRT1_2 * (1.1 * c0 - 1.2 * c4) * ck;
477 dE *= kappa;
478 }
479}
480
481void f_exact(const Vector &x, Vector &f)
482{
483 if (dim == 1)
484 {
485 real_t s0 = sin(kappa * x(0) + 0.0 * M_PI);
486 real_t s4 = sin(kappa * x(0) + 0.4 * M_PI);
487 real_t s9 = sin(kappa * x(0) + 0.9 * M_PI);
488
489 f(0) = 2.2 * s0 + 1.2 * M_SQRT1_2 * s4;
490 f(1) = 1.2 * (2.0 + kappa * kappa) * s4 +
491 M_SQRT1_2 * (1.1 * s0 + 1.3 * s9);
492 f(2) = 1.3 * (2.0 + kappa * kappa) * s9 + 1.2 * M_SQRT1_2 * s4;
493 }
494 else if (dim == 2)
495 {
496 real_t s0 = sin(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.0 * M_PI);
497 real_t s4 = sin(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.4 * M_PI);
498 real_t s9 = sin(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.9 * M_PI);
499
500 f(0) = 0.55 * (4.0 + kappa * kappa) * s0 +
501 0.6 * (M_SQRT2 - kappa * kappa) * s4;
502 f(1) = 0.55 * (M_SQRT2 - kappa * kappa) * s0 +
503 0.6 * (4.0 + kappa * kappa) * s4 +
504 0.65 * M_SQRT2 * s9;
505 f(2) = 0.6 * M_SQRT2 * s4 + 1.3 * (2.0 + kappa * kappa) * s9;
506 }
507 else
508 {
509 real_t s0 = sin(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.0 * M_PI);
510 real_t c0 = cos(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.0 * M_PI);
511 real_t s4 = sin(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.4 * M_PI);
512 real_t c4 = cos(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.4 * M_PI);
513 real_t s9 = sin(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.9 * M_PI);
514 real_t c9 = cos(kappa * M_SQRT1_2 * (x(0) + x(1)) + 0.9 * M_PI);
515 real_t sk = sin(kappa * x(2));
516 real_t ck = cos(kappa * x(2));
517
518 f(0) = 0.55 * (4.0 + 3.0 * kappa * kappa) * s0 * ck +
519 0.6 * (M_SQRT2 - kappa * kappa) * s4 * ck -
520 0.65 * M_SQRT2 * kappa * kappa * c9 * sk;
521
522 f(1) = 0.55 * (M_SQRT2 - kappa * kappa) * s0 * ck +
523 0.6 * (4.0 + 3.0 * kappa * kappa) * s4 * ck +
524 0.65 * M_SQRT2 * s9 * ck -
525 0.65 * M_SQRT2 * kappa * kappa * c9 * sk;
526
527 f(2) = 0.6 * M_SQRT2 * s4 * ck -
528 M_SQRT2 * kappa * kappa * (0.55 * c0 + 0.6 * c4) * sk
529 + 1.3 * (2.0 + kappa * kappa) * s9 * ck;
530 }
531}
T Max() const
Find the maximal element in the array, using the comparison operator < for class T.
Definition array.cpp:68
int Size() const
Return the logical size of the array.
Definition array.hpp:144
A coefficient that is constant across space and time.
Integrator for for Nedelec elements.
Vector coefficient defined as the Curl of a vector GridFunction.
Data type dense matrix using column-major storage.
Definition densemat.hpp:24
Collection of finite elements from the same family in multiple dimensions. This class is used to matc...
Definition fe_coll.hpp:27
Arbitrary order H1-conforming (continuous) finite elements.
Definition fe_coll.hpp:260
The Auxiliary-space Maxwell Solver in hypre.
Definition hypre.hpp:1845
PCG solver in hypre.
Definition hypre.hpp:1275
void SetPrintLevel(int print_lvl)
Definition hypre.cpp:4156
void SetPreconditioner(HypreSolver &precond)
Set the hypre solver to be used as a preconditioner.
Definition hypre.cpp:4161
virtual void Mult(const HypreParVector &b, HypreParVector &x) const
Solve Ax=b with hypre's PCG.
Definition hypre.cpp:4184
void SetMaxIter(int max_iter)
Definition hypre.cpp:4146
void SetTol(real_t tol)
Definition hypre.cpp:4136
Wrapper for hypre's ParCSR matrix class.
Definition hypre.hpp:388
static void Init()
Initialize hypre by calling HYPRE_Init() and set default options. After calling Hypre::Init(),...
Definition hypre.hpp:74
Scalar coefficient defined as the inner product of two vector coefficients.
Arbitrary order "L2-conforming" discontinuous finite elements.
Definition fe_coll.hpp:330
A matrix coefficient that is constant in space and time.
Vector coefficient defined as a product of a matrix coefficient and a vector coefficient.
Mesh data type.
Definition mesh.hpp:56
Array< int > bdr_attributes
A list of all unique boundary attributes used by the Mesh.
Definition mesh.hpp:282
int Dimension() const
Dimension of the reference space used within the elements.
Definition mesh.hpp:1160
void UniformRefinement(int i, const DSTable &, int *, int *, int *)
Definition mesh.cpp:10970
static bool Root()
Return true if the rank in MPI_COMM_WORLD is zero.
static int WorldRank()
Return the MPI rank in MPI_COMM_WORLD.
static int WorldSize()
Return the size of MPI_COMM_WORLD.
static void Init(int &argc, char **&argv, int required=default_thread_required, int *provided=nullptr)
Singleton creation with Mpi::Init(argc, argv).
Arbitrary order H(curl)-conforming Nedelec finite elements.
Definition fe_coll.hpp:465
Arbitrary order 3D H(curl)-conforming Nedelec finite elements in 1D.
Definition fe_coll.hpp:520
Arbitrary order 3D H(curl)-conforming Nedelec finite elements in 2D.
Definition fe_coll.hpp:584
Pointer to an Operator of a specified type.
Definition handle.hpp:34
OpType * As() const
Return the Operator pointer statically cast to a specified OpType. Similar to the method Get().
Definition handle.hpp:104
void ParseCheck(std::ostream &out=mfem::out)
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...
Definition optparser.hpp:82
Class for parallel bilinear form.
Abstract parallel finite element space.
Definition pfespace.hpp:29
void GetEssentialTrueDofs(const Array< int > &bdr_attr_is_ess, Array< int > &ess_tdof_list, int component=-1) const override
HYPRE_BigInt GlobalTrueVSize() const
Definition pfespace.hpp:285
Class for parallel grid function.
Definition pgridfunc.hpp:33
void Save(std::ostream &out) const override
void ProjectCoefficient(Coefficient &coeff) override
Project coeff Coefficient to this GridFunction. The projection computation depends on the choice of t...
real_t ComputeHCurlError(VectorCoefficient *exsol, VectorCoefficient *excurl, const IntegrationRule *irs[]=NULL) const override
Returns the error measured H(curl)-norm for ND elements.
Class for parallel linear form.
Class for parallel meshes.
Definition pmesh.hpp:34
void Print(std::ostream &out=mfem::out, const std::string &comments="") const override
Definition pmesh.cpp:4801
Arbitrary order H(div)-conforming Raviart-Thomas finite elements.
Definition fe_coll.hpp:386
void Mult(const Vector &x, Vector &y) const
Factor and solve the linear system .
Definition superlu.cpp:587
void SetOperator(const Operator &op)
Set the operator/matrix.
Definition superlu.cpp:475
Vector coefficient that is constant in space and time.
for VectorFiniteElements (Nedelec, Raviart-Thomas)
Definition lininteg.hpp:347
A general vector function coefficient.
Vector coefficient defined by a vector GridFunction.
Vector data type.
Definition vector.hpp:80
real_t sigma(const Vector &x)
Definition maxwell.cpp:102
real_t kappa
Definition ex31p.cpp:42
void CurlE_exact(const Vector &, Vector &)
Definition ex31p.cpp:441
int dim
Definition ex31p.cpp:43
void E_exact(const Vector &, Vector &)
Definition ex31p.cpp:418
void f_exact(const Vector &, Vector &)
Definition ex31p.cpp:481
real_t freq
Definition ex31p.cpp:42
int main()
real_t b
Definition lissajous.cpp:42
real_t a
Definition lissajous.cpp:41
const int visport
float real_t
Definition config.hpp:43
std::function< real_t(const Vector &)> f(real_t mass_coeff)
Definition lor_mms.hpp:30
const char vishost[]