MFEM v4.8.0
Finite element discretization library
Loading...
Searching...
No Matches
tesla_solver.cpp
Go to the documentation of this file.
1// Copyright (c) 2010-2025, Lawrence Livermore National Security, LLC. Produced
2// at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3// LICENSE and NOTICE for details. LLNL-CODE-806117.
4//
5// This file is part of the MFEM library. For more information and source code
6// availability visit https://mfem.org.
7//
8// MFEM is free software; you can redistribute it and/or modify it under the
9// terms of the BSD-3 license. We welcome feedback and contributions, see file
10// CONTRIBUTING.md for details.
11
12#include "tesla_solver.hpp"
13
14#ifdef MFEM_USE_MPI
15
16using namespace std;
17
18namespace mfem
19{
20
21using namespace common;
22
23namespace electromagnetics
24{
25
27 Array<int> & kbcs,
28 Array<int> & vbcs, Vector & vbcv,
29 Coefficient & muInvCoef,
30 void (*a_bc )(const Vector&, Vector&),
31 void (*j_src)(const Vector&, Vector&),
32 void (*m_src)(const Vector&, Vector&))
33 : myid_(0),
34 num_procs_(1),
35 order_(order),
36 pmesh_(&pmesh),
37 visit_dc_(NULL),
38 H1FESpace_(NULL),
39 HCurlFESpace_(NULL),
40 HDivFESpace_(NULL),
41 curlMuInvCurl_(NULL),
42 hCurlMass_(NULL),
43 hDivHCurlMuInv_(NULL),
44 weakCurlMuInv_(NULL),
45 grad_(NULL),
46 curl_(NULL),
47 a_(NULL),
48 b_(NULL),
49 h_(NULL),
50 jr_(NULL),
51 j_(NULL),
52 k_(NULL),
53 m_(NULL),
54 bd_(NULL),
55 jd_(NULL),
56 DivFreeProj_(NULL),
57 SurfCur_(NULL),
58 muInvCoef_(&muInvCoef),
59 aBCCoef_(NULL),
60 jCoef_(NULL),
61 mCoef_(NULL),
62 a_bc_(a_bc),
63 j_src_(j_src),
64 m_src_(m_src)
65{
66 // Initialize MPI variables
67 MPI_Comm_size(pmesh_->GetComm(), &num_procs_);
68 MPI_Comm_rank(pmesh_->GetComm(), &myid_);
69
70 // Define compatible parallel finite element spaces on the parallel
71 // mesh. Here we use arbitrary order H1, Nedelec, and Raviart-Thomas finite
72 // elements.
73 H1FESpace_ = new H1_ParFESpace(pmesh_,order,pmesh_->Dimension());
74 HCurlFESpace_ = new ND_ParFESpace(pmesh_,order,pmesh_->Dimension());
75 HDivFESpace_ = new RT_ParFESpace(pmesh_,order,pmesh_->Dimension());
76
77 int irOrder = pmesh.GetTypicalElementTransformation()->OrderW()
78 + 2 * order;
80 const IntegrationRule * ir = &IntRules.Get(geom, irOrder);
81
82 // Select surface attributes for Dirichlet BCs
83 ess_bdr_.SetSize(pmesh.bdr_attributes.Max());
84 non_k_bdr_.SetSize(pmesh.bdr_attributes.Max());
85 ess_bdr_ = 1; // All outer surfaces
86 non_k_bdr_ = 1; // Surfaces without applied surface currents
87
88 for (int i=0; i<kbcs.Size(); i++)
89 {
90 non_k_bdr_[kbcs[i]-1] = 0;
91 }
92
93 // Setup various coefficients
94
95 // Vector Potential on the outer surface
96 if ( a_bc_ == NULL )
97 {
98 Vector Zero(3);
99 Zero = 0.0;
100 aBCCoef_ = new VectorConstantCoefficient(Zero);
101 }
102 else
103 {
104 aBCCoef_ = new VectorFunctionCoefficient(pmesh_->SpaceDimension(),
105 *a_bc_);
106 }
107
108 // Volume Current Density
109 if ( j_src_ != NULL )
110 {
111 jCoef_ = new VectorFunctionCoefficient(pmesh_->SpaceDimension(),
112 j_src_);
113 }
114
115 // Magnetization
116 if ( m_src_ != NULL )
117 {
118 mCoef_ = new VectorFunctionCoefficient(pmesh_->SpaceDimension(),
119 m_src_);
120 }
121
122 // Bilinear Forms
123 curlMuInvCurl_ = new ParBilinearForm(HCurlFESpace_);
124 curlMuInvCurl_->AddDomainIntegrator(new CurlCurlIntegrator(*muInvCoef_));
125
126 BilinearFormIntegrator * hCurlMassInteg = new VectorFEMassIntegrator;
127 hCurlMassInteg->SetIntRule(ir);
128 hCurlMass_ = new ParBilinearForm(HCurlFESpace_);
129 hCurlMass_->AddDomainIntegrator(hCurlMassInteg);
130
131 BilinearFormIntegrator * hDivHCurlInteg =
132 new VectorFEMassIntegrator(*muInvCoef_);
133 hDivHCurlInteg->SetIntRule(ir);
134 hDivHCurlMuInv_ = new ParMixedBilinearForm(HDivFESpace_, HCurlFESpace_);
135 hDivHCurlMuInv_->AddDomainIntegrator(hDivHCurlInteg);
136
137 // Discrete Curl operator
138 curl_ = new ParDiscreteCurlOperator(HCurlFESpace_, HDivFESpace_);
139
140 // Build grid functions
141 a_ = new ParGridFunction(HCurlFESpace_);
142 b_ = new ParGridFunction(HDivFESpace_);
143 h_ = new ParGridFunction(HCurlFESpace_);
144 bd_ = new ParGridFunction(HCurlFESpace_);
145 jd_ = new ParGridFunction(HCurlFESpace_);
146
147 if ( jCoef_ || kbcs.Size() > 0 )
148 {
149 grad_ = new ParDiscreteGradOperator(H1FESpace_, HCurlFESpace_);
150 }
151 if ( jCoef_ )
152 {
153 jr_ = new ParGridFunction(HCurlFESpace_);
154 j_ = new ParGridFunction(HCurlFESpace_);
155 DivFreeProj_ = new DivergenceFreeProjector(*H1FESpace_, *HCurlFESpace_,
156 irOrder, NULL, NULL, grad_);
157 }
158
159 if ( kbcs.Size() > 0 )
160 {
161 k_ = new ParGridFunction(HCurlFESpace_);
162
163 // Object to solve the subproblem of computing surface currents
164 SurfCur_ = new SurfaceCurrent(*H1FESpace_, *grad_,
165 kbcs, vbcs, vbcv);
166 }
167
168 if ( mCoef_ )
169 {
170 m_ = new ParGridFunction(HDivFESpace_);
171
172 weakCurlMuInv_ = new ParMixedBilinearForm(HDivFESpace_, HCurlFESpace_);
173 weakCurlMuInv_->AddDomainIntegrator(
174 new VectorFECurlIntegrator(*muInvCoef_));
175 }
176}
177
179{
180 delete jCoef_;
181 delete mCoef_;
182 delete aBCCoef_;
183
184 delete DivFreeProj_;
185 delete SurfCur_;
186
187 delete a_;
188 delete b_;
189 delete h_;
190 delete jr_;
191 delete j_;
192 delete k_;
193 delete m_;
194 delete bd_;
195 delete jd_;
196
197 delete grad_;
198 delete curl_;
199
200 delete curlMuInvCurl_;
201 delete hCurlMass_;
202 delete hDivHCurlMuInv_;
203 delete weakCurlMuInv_;
204
205 delete H1FESpace_;
206 delete HCurlFESpace_;
207 delete HDivFESpace_;
208
209 map<string,socketstream*>::iterator mit;
210 for (mit=socks_.begin(); mit!=socks_.end(); mit++)
211 {
212 delete mit->second;
213 }
214}
215
218{
219 return HCurlFESpace_->GlobalTrueVSize();
220}
221
222void
224{
225 HYPRE_BigInt size_h1 = H1FESpace_->GlobalTrueVSize();
226 HYPRE_BigInt size_nd = HCurlFESpace_->GlobalTrueVSize();
227 HYPRE_BigInt size_rt = HDivFESpace_->GlobalTrueVSize();
228 if (myid_ == 0)
229 {
230 cout << "Number of H1 unknowns: " << size_h1 << endl;
231 cout << "Number of H(Curl) unknowns: " << size_nd << endl;
232 cout << "Number of H(Div) unknowns: " << size_rt << endl;
233 }
234}
235
236void
238{
239 if (myid_ == 0) { cout << "Assembling ..." << flush; }
240
241 curlMuInvCurl_->Assemble();
242 curlMuInvCurl_->Finalize();
243
244 hDivHCurlMuInv_->Assemble();
245 hDivHCurlMuInv_->Finalize();
246
247 hCurlMass_->Assemble();
248 hCurlMass_->Finalize();
249
250 curl_->Assemble();
251 curl_->Finalize();
252
253 if ( grad_ )
254 {
255 grad_->Assemble();
256 grad_->Finalize();
257 }
258 if ( weakCurlMuInv_ )
259 {
260 weakCurlMuInv_->Assemble();
261 weakCurlMuInv_->Finalize();
262 }
263
264 if (myid_ == 0) { cout << " done." << endl; }
265}
266
267void
269{
270 if (myid_ == 0) { cout << "Updating ..." << endl; }
271
272 // Inform the spaces that the mesh has changed
273 // Note: we don't need to interpolate any GridFunctions on the new mesh
274 // so we pass 'false' to skip creation of any transformation matrices.
275 H1FESpace_->Update(false);
276 HCurlFESpace_->Update(false);
277 HDivFESpace_->Update(false);
278
279 HCurlFESpace_->GetEssentialTrueDofs(ess_bdr_, ess_bdr_tdofs_);
280
281 // Inform the grid functions that the space has changed.
282 a_->Update();
283 h_->Update();
284 b_->Update();
285 bd_->Update();
286 jd_->Update();
287 if ( jr_ ) { jr_->Update(); }
288 if ( j_ ) { j_->Update(); }
289 if ( k_ ) { k_->Update(); }
290 if ( m_ ) { m_->Update(); }
291
292 // Inform the bilinear forms that the space has changed.
293 curlMuInvCurl_->Update();
294 hCurlMass_->Update();
295 hDivHCurlMuInv_->Update();
296 if ( weakCurlMuInv_ ) { weakCurlMuInv_->Update(); }
297
298 // Inform the other objects that the space has changed.
299 curl_->Update();
300 if ( grad_ ) { grad_->Update(); }
301 if ( DivFreeProj_ ) { DivFreeProj_->Update(); }
302 if ( SurfCur_ ) { SurfCur_->Update(); }
303}
304
305void
307{
308 if (myid_ == 0) { cout << "Running solver ... " << endl; }
309
310 // Initialize the magnetic vector potential with its boundary conditions
311 *a_ = 0.0;
312
313 // Apply surface currents if available
314 if ( k_ )
315 {
316 SurfCur_->ComputeSurfaceCurrent(*k_);
317 *a_ = *k_;
318 }
319
320 // Apply uniform B boundary condition on remaining surfaces
321 a_->ProjectBdrCoefficientTangent(*aBCCoef_, non_k_bdr_);
322
323 // Initialize the RHS vector to zero
324 *jd_ = 0.0;
325
326 // Initialize the volumetric current density
327 if ( jr_ )
328 {
329 jr_->ProjectCoefficient(*jCoef_);
330
331 // Compute the discretely divergence-free portion of jr_
332 DivFreeProj_->Mult(*jr_, *j_);
333
334 // Compute the dual of j_
335 hCurlMass_->AddMult(*j_, *jd_);
336 }
337
338 // Initialize the Magnetization
339 if ( m_ )
340 {
341 m_->ProjectCoefficient(*mCoef_);
342 weakCurlMuInv_->AddMult(*m_, *jd_, mu0_);
343 }
344
345 // Apply Dirichlet BCs to matrix and right hand side and otherwise
346 // prepare the linear system
347 HypreParMatrix CurlMuInvCurl;
348 HypreParVector A(HCurlFESpace_);
349 HypreParVector RHS(HCurlFESpace_);
350
351 curlMuInvCurl_->FormLinearSystem(ess_bdr_tdofs_, *a_, *jd_, CurlMuInvCurl,
352 A, RHS);
353
354 // Define and apply a parallel PCG solver for AX=B with the AMS
355 // preconditioner from hypre.
356 HypreAMS ams(CurlMuInvCurl, HCurlFESpace_);
357 ams.SetSingularProblem();
358
359 HyprePCG pcg (CurlMuInvCurl);
360 pcg.SetTol(1e-12);
361 pcg.SetMaxIter(50);
362 pcg.SetPrintLevel(2);
363 pcg.SetPreconditioner(ams);
364 pcg.Mult(RHS, A);
365
366 // Extract the parallel grid function corresponding to the finite
367 // element approximation A. This is the local solution on each
368 // processor.
369 curlMuInvCurl_->RecoverFEMSolution(A, *jd_, *a_);
370
371 // Compute the negative Gradient of the solution vector. This is
372 // the magnetic field corresponding to the scalar potential
373 // represented by phi.
374 curl_->Mult(*a_, *b_);
375
376 // Compute magnetic field (H) from B and M
377 if (myid_ == 0) { cout << "Computing H ... " << flush; }
378
379 hDivHCurlMuInv_->Mult(*b_, *bd_);
380 if ( m_ )
381 {
382 hDivHCurlMuInv_->AddMult(*m_, *bd_, -1.0 * mu0_);
383 }
384
385 HypreParMatrix MassHCurl;
386 Vector BD, H;
387
388 Array<int> dbc_dofs_h;
389 hCurlMass_->FormLinearSystem(dbc_dofs_h, *h_, *bd_, MassHCurl, H, BD);
390
391 HyprePCG pcgM(MassHCurl);
392 pcgM.SetTol(1e-12);
393 pcgM.SetMaxIter(500);
394 pcgM.SetPrintLevel(0);
395 HypreDiagScale diagM;
396 pcgM.SetPreconditioner(diagM);
397 pcgM.Mult(BD, H);
398
399 hCurlMass_->RecoverFEMSolution(H, *bd_, *h_);
400
401 if (myid_ == 0) { cout << "done." << flush; }
402
403 if (myid_ == 0) { cout << " Solver done. " << endl; }
404}
405
406void
408{
409 if (myid_ == 0) { cout << "Estimating Error ... " << flush; }
410
411 // Space for the discontinuous (original) flux
412 CurlCurlIntegrator flux_integrator(*muInvCoef_);
413 RT_FECollection flux_fec(order_-1, pmesh_->SpaceDimension());
414 ParFiniteElementSpace flux_fes(pmesh_, &flux_fec);
415
416 // Space for the smoothed (conforming) flux
417 int norm_p = 1;
418 ND_FECollection smooth_flux_fec(order_, pmesh_->Dimension());
419 ParFiniteElementSpace smooth_flux_fes(pmesh_, &smooth_flux_fec);
420
421 L2ZZErrorEstimator(flux_integrator, *a_,
422 smooth_flux_fes, flux_fes, errors, norm_p);
423
424 if (myid_ == 0) { cout << "done." << endl; }
425}
426
427void
429{
430 visit_dc_ = &visit_dc;
431
432 visit_dc.RegisterField("A", a_);
433 visit_dc.RegisterField("B", b_);
434 visit_dc.RegisterField("H", h_);
435 if ( j_ ) { visit_dc.RegisterField("J", j_); }
436 if ( k_ ) { visit_dc.RegisterField("K", k_); }
437 if ( m_ ) { visit_dc.RegisterField("M", m_); }
438 if ( SurfCur_ ) { visit_dc.RegisterField("Psi", SurfCur_->GetPsi()); }
439}
440
441void
443{
444 if ( visit_dc_ )
445 {
446 if (myid_ == 0) { cout << "Writing VisIt files ..." << flush; }
447
448 HYPRE_BigInt prob_size = this->GetProblemSize();
449 visit_dc_->SetCycle(it);
450 visit_dc_->SetTime(prob_size);
451 visit_dc_->Save();
452
453 if (myid_ == 0) { cout << " done." << endl; }
454 }
455}
456
457void
459{
460 if ( myid_ == 0 ) { cout << "Opening GLVis sockets." << endl; }
461
462 socks_["A"] = new socketstream;
463 socks_["A"]->precision(8);
464
465 socks_["B"] = new socketstream;
466 socks_["B"]->precision(8);
467
468 socks_["H"] = new socketstream;
469 socks_["H"]->precision(8);
470
471 if ( j_ )
472 {
473 socks_["J"] = new socketstream;
474 socks_["J"]->precision(8);
475 }
476 if ( k_ )
477 {
478 socks_["K"] = new socketstream;
479 socks_["K"]->precision(8);
480
481 socks_["Psi"] = new socketstream;
482 socks_["Psi"]->precision(8);
483 }
484 if ( m_ )
485 {
486 socks_["M"] = new socketstream;
487 socks_["M"]->precision(8);
488 }
489 if ( myid_ == 0 ) { cout << "GLVis sockets open." << endl; }
490}
491
492void
494{
495 if (myid_ == 0) { cout << "Sending data to GLVis ..." << flush; }
496
497 char vishost[] = "localhost";
498
499 int Wx = 0, Wy = 0; // window position
500 int Ww = 350, Wh = 350; // window size
501 int offx = Ww+10, offy = Wh+45; // window offsets
502
503 VisualizeField(*socks_["A"], vishost, visport,
504 *a_, "Vector Potential (A)", Wx, Wy, Ww, Wh);
505 Wx += offx;
506
507 VisualizeField(*socks_["B"], vishost, visport,
508 *b_, "Magnetic Flux Density (B)", Wx, Wy, Ww, Wh);
509 Wx += offx;
510
511 VisualizeField(*socks_["H"], vishost, visport,
512 *h_, "Magnetic Field (H)", Wx, Wy, Ww, Wh);
513 Wx += offx;
514
515 if ( j_ )
516 {
517 VisualizeField(*socks_["J"], vishost, visport,
518 *j_, "Current Density (J)", Wx, Wy, Ww, Wh);
519 }
520
521 Wx = 0; Wy += offy; // next line
522
523 if ( k_ )
524 {
525 VisualizeField(*socks_["K"], vishost, visport,
526 *k_, "Surface Current Density (K)", Wx, Wy, Ww, Wh);
527 Wx += offx;
528
529 VisualizeField(*socks_["Psi"], vishost, visport,
530 *SurfCur_->GetPsi(),
531 "Surface Current Potential (Psi)", Wx, Wy, Ww, Wh);
532 Wx += offx;
533 }
534 if ( m_ )
535 {
536 VisualizeField(*socks_["M"], vishost, visport,
537 *m_, "Magnetization (M)", Wx, Wy, Ww, Wh);
538 // Wx += offx; // not used
539 }
540 if (myid_ == 0) { cout << " done." << endl; }
541}
542
545 Array<int> & kbcs,
546 Array<int> & vbcs, Vector & vbcv)
547 : H1FESpace_(&H1FESpace),
548 grad_(&grad),
549 kbcs_(&kbcs),
550 vbcs_(&vbcs),
551 vbcv_(&vbcv),
552 s0_(NULL),
553 psi_(NULL),
554 rhs_(NULL)
555{
556 // Initialize MPI variables
557 MPI_Comm_rank(H1FESpace_->GetParMesh()->GetComm(), &myid_);
558
559 s0_ = new ParBilinearForm(H1FESpace_);
561 s0_->Assemble();
562 s0_->Finalize();
563 S0_ = new HypreParMatrix;
564
565 AttrToMarker(H1FESpace_->GetParMesh()->bdr_attributes.Max(),
566 *vbcs_, ess_bdr_);
567 H1FESpace_->GetEssentialTrueDofs(ess_bdr_, ess_bdr_tdofs_);
568
569 non_k_bdr_.SetSize(H1FESpace_->GetParMesh()->bdr_attributes.Max());
570 non_k_bdr_ = 1;
571 for (int i=0; i<kbcs_->Size(); i++)
572 {
573 non_k_bdr_[(*kbcs_)[i]-1] = 0;
574 }
575
576 psi_ = new ParGridFunction(H1FESpace_);
577 rhs_ = new ParGridFunction(H1FESpace_);
578
579 pcg_ = NULL;
580 amg_ = NULL;
581}
582
584{
585 delete psi_;
586 delete rhs_;
587
588 delete pcg_;
589 delete amg_;
590
591 delete S0_;
592
593 delete s0_;
594}
595
596void
598{
599 delete pcg_;
600 delete amg_;
601
602 amg_ = new HypreBoomerAMG(*S0_);
603 amg_->SetPrintLevel(0);
604 pcg_ = new HyprePCG(*S0_);
605 pcg_->SetTol(1e-14);
606 pcg_->SetMaxIter(200);
607 pcg_->SetPrintLevel(0);
608 pcg_->SetPreconditioner(*amg_);
609}
610
611void
613{
614 if (myid_ == 0) { cout << "Computing K ... " << flush; }
615
616 // Apply piecewise constant voltage boundary condition
617 *psi_ = 0.0;
618 *rhs_ = 0.0;
619 Array<int> vbc_bdr_attr(H1FESpace_->GetParMesh()->bdr_attributes.Max());
620 for (int i=0; i<vbcs_->Size(); i++)
621 {
622 ConstantCoefficient voltage((*vbcv_)[i]);
623 vbc_bdr_attr = 0;
624 vbc_bdr_attr[(*vbcs_)[i]-1] = 1;
625 psi_->ProjectBdrCoefficient(voltage, vbc_bdr_attr);
626 }
627
628 // Apply essential BC and form linear system
629 s0_->FormLinearSystem(ess_bdr_tdofs_, *psi_, *rhs_, *S0_, Psi_, RHS_);
630
631 // Solve the linear system for Psi
632 if ( pcg_ == NULL ) { this->InitSolver(); }
633 pcg_->Mult(RHS_, Psi_);
634
635 // Compute the parallel grid function corresponding to Psi
636 s0_->RecoverFEMSolution(Psi_, *rhs_, *psi_);
637
638 // Compute the surface current from psi
639 grad_->Mult(*psi_, k);
640
641 // Force the tangential part of k to be zero away from the intended surfaces
642 Vector vZero(3); vZero = 0.0;
643 VectorConstantCoefficient Zero(vZero);
644 k.ProjectBdrCoefficientTangent(Zero, non_k_bdr_);
645
646 if (myid_ == 0) { cout << "done." << endl; }
647}
648
649void
651{
652 delete pcg_; pcg_ = NULL;
653 delete amg_; amg_ = NULL;
654 delete S0_; S0_ = new HypreParMatrix;
655
656 psi_->Update();
657 rhs_->Update();
658
659 s0_->Update();
660 s0_->Assemble();
661 s0_->Finalize();
662
663 H1FESpace_->GetEssentialTrueDofs(ess_bdr_, ess_bdr_tdofs_);
664}
665
666} // namespace electromagnetics
667
668} // namespace mfem
669
670#endif // MFEM_USE_MPI
T Max() const
Find the maximal element in the array, using the comparison operator < for class T.
Definition array.cpp:68
void SetSize(int nsize)
Change the logical size of the array, keep existing entries.
Definition array.hpp:758
int Size() const
Return the logical size of the array.
Definition array.hpp:147
Abstract base class BilinearFormIntegrator.
void AddDomainIntegrator(BilinearFormIntegrator *bfi)
Adds new Domain Integrator. Assumes ownership of bfi.
void Finalize(int skip_zeros=1) override
Finalizes the matrix initialization if the AssemblyLevel is AssemblyLevel::LEGACY....
void AddBoundaryIntegrator(BilinearFormIntegrator *bfi)
Adds new Boundary Integrator. Assumes ownership of bfi.
void AddMult(const Vector &x, Vector &y, const real_t a=1.0) const override
Add the matrix vector multiple to a vector: .
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
A coefficient that is constant across space and time.
Integrator for for Nedelec elements.
void SetCycle(int c)
Set time cycle (for time-dependent simulations)
void SetTime(real_t t)
Set physical time (for time-dependent simulations)
virtual void Assemble(int skip_zeros=1)
Construct the internal matrix representation of the discrete linear operator.
virtual int OrderW() const =0
Return the order of the determinant of the Jacobian (weight) of the transformation.
The Auxiliary-space Maxwell Solver in hypre.
Definition hypre.hpp:1871
void SetSingularProblem()
Set this option when solving a curl-curl problem with zero mass term.
Definition hypre.hpp:1929
The BoomerAMG solver in hypre.
Definition hypre.hpp:1717
void SetPrintLevel(int print_level)
Definition hypre.hpp:1797
Jacobi preconditioner in hypre.
Definition hypre.hpp:1507
PCG solver in hypre.
Definition hypre.hpp:1301
void Mult(const HypreParVector &b, HypreParVector &x) const override
Solve Ax=b with hypre's PCG.
Definition hypre.cpp:4242
void SetPrintLevel(int print_lvl)
Definition hypre.cpp:4214
void SetPreconditioner(HypreSolver &precond)
Set the hypre solver to be used as a preconditioner.
Definition hypre.cpp:4219
void SetMaxIter(int max_iter)
Definition hypre.cpp:4204
void SetTol(real_t tol)
Definition hypre.cpp:4194
Wrapper for hypre's ParCSR matrix class.
Definition hypre.hpp:408
Wrapper for hypre's parallel vector class.
Definition hypre.hpp:219
Class for an integration rule - an Array of IntegrationPoint.
Definition intrules.hpp:100
const IntegrationRule & Get(int GeomType, int Order)
Returns an integration rule for given GeomType and Order.
virtual void SetIntRule(const IntegrationRule *ir)
Prescribe a fixed IntegrationRule to use, or set to null to let the integrator choose an appropriate ...
Array< int > bdr_attributes
A list of all unique boundary attributes used by the Mesh.
Definition mesh.hpp:290
Geometry::Type GetTypicalElementGeometry() const
If the local mesh is not empty, return GetElementGeometry(0); otherwise, return a typical Geometry pr...
Definition mesh.cpp:1528
int Dimension() const
Dimension of the reference space used within the elements.
Definition mesh.hpp:1216
ElementTransformation * GetTypicalElementTransformation()
If the local mesh is not empty return GetElementTransformation(0); otherwise, return the identity tra...
Definition mesh.cpp:390
int SpaceDimension() const
Dimension of the physical space containing the mesh.
Definition mesh.hpp:1219
void AddMult(const Vector &x, Vector &y, const real_t a=1.0) const override
Add the matrix vector multiple to a vector: .
void Assemble(int skip_zeros=1)
void Finalize(int skip_zeros=1) override
Finalizes the matrix initialization if the AssemblyLevel is AssemblyLevel::LEGACY.
void Update()
Must be called after making changes to trial_fes or test_fes.
void Mult(const Vector &x, Vector &y) const override
Matrix multiplication: .
void AddDomainIntegrator(BilinearFormIntegrator *bfi)
Adds a domain integrator. Assumes ownership of bfi.
Arbitrary order H(curl)-conforming Nedelec finite elements.
Definition fe_coll.hpp:482
Class for parallel bilinear form.
void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x) override
void Assemble(int skip_zeros=1)
Assemble the local matrix.
void Update(FiniteElementSpace *nfes=NULL) override
Update the FiniteElementSpace and delete all data associated with the old one.
void FormLinearSystem(const Array< int > &ess_tdof_list, Vector &x, Vector &b, OperatorHandle &A, Vector &X, Vector &B, int copy_interior=0) override
Form the linear system A X = B, corresponding to this bilinear form and the linear form b(....
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:346
ParMesh * GetParMesh() const
Definition pfespace.hpp:338
void Update(bool want_transform=true) override
Class for parallel grid function.
Definition pgridfunc.hpp:50
void ProjectBdrCoefficient(Coefficient *coeff[], VectorCoefficient *vcoeff, const Array< int > &attr)
void ProjectCoefficient(Coefficient &coeff) override
Project coeff Coefficient to this GridFunction. The projection computation depends on the choice of t...
void Update() override
Transform by the Space UpdateMatrix (e.g., on Mesh change).
Definition pgridfunc.cpp:91
void ProjectBdrCoefficientTangent(VectorCoefficient &vcoeff, const Array< int > &bdr_attr) override
Project the tangential components of the given VectorCoefficient on the boundary. Only boundary attri...
Class for parallel meshes.
Definition pmesh.hpp:34
MPI_Comm GetComm() const
Definition pmesh.hpp:402
Class for parallel bilinear form using different test and trial FE spaces.
Arbitrary order H(div)-conforming Raviart-Thomas finite elements.
Definition fe_coll.hpp:403
Vector coefficient that is constant in space and time.
A general vector function coefficient.
Vector data type.
Definition vector.hpp:82
Data collection with VisIt I/O routines.
void Save() override
Save the collection and a VisIt root file.
void RegisterField(const std::string &field_name, GridFunction *gf) override
Add a grid function to the collection and update the root file.
void Mult(const Vector &x, Vector &y) const override
Operator application: y=A(x).
void ComputeSurfaceCurrent(ParGridFunction &k)
SurfaceCurrent(ParFiniteElementSpace &H1FESpace, ParDiscreteGradOperator &Grad, Array< int > &kbcs, Array< int > &vbcs, Vector &vbcv)
void DisplayToGLVis(int visport=19916)
TeslaSolver(ParMesh &pmesh, int order, Array< int > &kbcs, Array< int > &vbcs, Vector &vbcv, Coefficient &muInvCoef, void(*a_bc)(const Vector &, Vector &), void(*j_src)(const Vector &, Vector &), void(*m_src)(const Vector &, Vector &))
void RegisterVisItFields(VisItDataCollection &visit_dc)
void j_src(const Vector &x, real_t t, Vector &j)
Definition maxwell.cpp:108
HYPRE_Int HYPRE_BigInt
int Ww
int Wy
int Wx
int offx
int Wh
int offy
void VisualizeField(socketstream &sock, const char *vishost, int visport, GridFunction &gf, const char *title, int x, int y, int w, int h, const char *keys, bool vec)
void AttrToMarker(int max_attr, const Array< int > &attrs, Array< int > &marker)
Convert a set of attribute numbers to a marker array.
real_t L2ZZErrorEstimator(BilinearFormIntegrator &flux_integrator, const ParGridFunction &x, ParFiniteElementSpace &smooth_flux_fes, ParFiniteElementSpace &flux_fes, Vector &errors, int norm_p, real_t solver_tol, int solver_max_it)
IntegrationRules IntRules(0, Quadrature1D::GaussLegendre)
A global object with all integration rules (defined in intrules.cpp)
Definition intrules.hpp:492
const char vishost[]