MFEM  v4.6.0
Finite element discretization library
complex_fem.hpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2023, 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 #ifndef MFEM_COMPLEX_FEM
13 #define MFEM_COMPLEX_FEM
14 
15 #include "../linalg/complex_operator.hpp"
16 #include "gridfunc.hpp"
17 #include "linearform.hpp"
18 #include "bilinearform.hpp"
19 #ifdef MFEM_USE_MPI
20 #include "pgridfunc.hpp"
21 #include "plinearform.hpp"
22 #include "pbilinearform.hpp"
23 #endif
24 #include <complex>
25 
26 namespace mfem
27 {
28 
29 /// Class for complex-valued grid function - real + imaginary part Vector with
30 /// associated FE space.
32 {
33 private:
34  GridFunction * gfr;
35  GridFunction * gfi;
36 
37 protected:
38  void Destroy() { delete gfr; delete gfi; }
39 
40 public:
41  /** @brief Construct a ComplexGridFunction associated with the
42  FiniteElementSpace @a *f. */
44 
45  void Update();
46 
47  /// Assign constant values to the ComplexGridFunction data.
48  ComplexGridFunction &operator=(const std::complex<double> & value)
49  { *gfr = value.real(); *gfi = value.imag(); return *this; }
50 
51  virtual void ProjectCoefficient(Coefficient &real_coeff,
52  Coefficient &imag_coeff);
53  virtual void ProjectCoefficient(VectorCoefficient &real_vcoeff,
54  VectorCoefficient &imag_vcoeff);
55 
56  virtual void ProjectBdrCoefficient(Coefficient &real_coeff,
57  Coefficient &imag_coeff,
58  Array<int> &attr);
59  virtual void ProjectBdrCoefficientNormal(VectorCoefficient &real_coeff,
60  VectorCoefficient &imag_coeff,
61  Array<int> &attr);
62  virtual void ProjectBdrCoefficientTangent(VectorCoefficient &real_coeff,
63  VectorCoefficient &imag_coeff,
64  Array<int> &attr);
65 
66  FiniteElementSpace *FESpace() { return gfr->FESpace(); }
67  const FiniteElementSpace *FESpace() const { return gfr->FESpace(); }
68 
69  GridFunction & real() { return *gfr; }
70  GridFunction & imag() { return *gfi; }
71  const GridFunction & real() const { return *gfr; }
72  const GridFunction & imag() const { return *gfi; }
73 
74  /// Update the memory location of the real and imaginary GridFunction @a gfr
75  /// and @a gfi to match the ComplexGridFunction.
76  void Sync() { gfr->SyncMemory(*this); gfi->SyncMemory(*this); }
77 
78  /// Update the alias memory location of the real and imaginary GridFunction
79  /// @a gfr and @a gfi to match the ComplexGridFunction.
80  void SyncAlias() { gfr->SyncAliasMemory(*this); gfi->SyncAliasMemory(*this); }
81 
82  /// Destroys the grid function.
83  virtual ~ComplexGridFunction() { Destroy(); }
84 
85 };
86 
87 /** Class for a complex-valued linear form
88 
89  The @a convention argument in the class's constructor is documented in the
90  mfem::ComplexOperator class found in linalg/complex_operator.hpp.
91 
92  When supplying integrators to the ComplexLinearForm either the real or
93  imaginary integrator can be NULL. This indicates that the corresponding
94  portion of the complex-valued field is equal to zero.
95  */
96 class ComplexLinearForm : public Vector
97 {
98 private:
100 
101 protected:
104 
105 public:
108  convention = ComplexOperator::HERMITIAN);
109 
110  /** @brief Create a ComplexLinearForm on the FiniteElementSpace @a fes, using
111  the same integrators as the LinearForms @a lf_r (real) and @a lf_i (imag).
112 
113  The pointer @a fes is not owned by the newly constructed object.
114 
115  The integrators are copied as pointers and they are not owned by the
116  newly constructed ComplexLinearForm. */
119  convention = ComplexOperator::HERMITIAN);
120 
121  virtual ~ComplexLinearForm();
122 
125  convention) { conv = convention; }
126 
127  /// Adds new Domain Integrator.
129  LinearFormIntegrator *lfi_imag);
130 
131  /// Adds new Domain Integrator, restricted to the given attributes.
133  LinearFormIntegrator *lfi_imag,
134  Array<int> &elem_attr_marker);
135 
136  /// Adds new Boundary Integrator.
138  LinearFormIntegrator *lfi_imag);
139 
140  /** @brief Add new Boundary Integrator, restricted to the given boundary
141  attributes.
142 
143  Assumes ownership of @a lfi_real and @a lfi_imag.
144 
145  The array @a bdr_attr_marker is stored internally as a pointer to the
146  given Array<int> object. */
148  LinearFormIntegrator *lfi_imag,
149  Array<int> &bdr_attr_marker);
150 
151  /// Adds new Boundary Face Integrator. Assumes ownership of @a lfi.
153  LinearFormIntegrator *lfi_imag);
154 
155  /** @brief Add new Boundary Face Integrator, restricted to the given boundary
156  attributes.
157 
158  Assumes ownership of @a lfi_real and @a lfi_imag.
159 
160  The array @a bdr_attr_marker is stored internally as a pointer to the
161  given Array<int> object. */
163  LinearFormIntegrator *lfi_imag,
164  Array<int> &bdr_attr_marker);
165 
166  FiniteElementSpace *FESpace() const { return lfr->FESpace(); }
167 
168  LinearForm & real() { return *lfr; }
169  LinearForm & imag() { return *lfi; }
170  const LinearForm & real() const { return *lfr; }
171  const LinearForm & imag() const { return *lfi; }
172 
173  /// Update the memory location of the real and imaginary LinearForm @a lfr
174  /// and @a lfi to match the ComplexLinearForm.
175  void Sync() { lfr->SyncMemory(*this); lfi->SyncMemory(*this); }
176 
177  /// Update the alias memory location of the real and imaginary LinearForm @a
178  /// lfr and @a lfi to match the ComplexLinearForm.
179  void SyncAlias() { lfr->SyncAliasMemory(*this); lfi->SyncAliasMemory(*this); }
180 
181  void Update();
182  void Update(FiniteElementSpace *f);
183 
184  /// Assembles the linear form i.e. sums over all domain/bdr integrators.
185  void Assemble();
186 
187  std::complex<double> operator()(const ComplexGridFunction &gf) const;
188 };
189 
190 
191 /** Class for sesquilinear form
192 
193  A sesquilinear form is a generalization of a bilinear form to complex-valued
194  fields. Sesquilinear forms are linear in the second argument but the first
195  argument involves a complex conjugate in the sense that:
196 
197  a(alpha u, beta v) = conj(alpha) beta a(u, v)
198 
199  The @a convention argument in the class's constructor is documented in the
200  mfem::ComplexOperator class found in linalg/complex_operator.hpp.
201 
202  When supplying integrators to the SesquilinearForm either the real or
203  imaginary integrator can be NULL. This indicates that the corresponding
204  portion of the complex-valued material coefficient is equal to zero.
205 */
207 {
208 private:
210 
211  /** This data member allows one to specify what should be done to the
212  diagonal matrix entries and corresponding RHS values upon elimination of
213  the constrained DoFs. */
215 
216  BilinearForm *blfr;
217  BilinearForm *blfi;
218 
219  /* These methods check if the real/imag parts of the sesquilinear form are
220  not empty */
221  bool RealInteg();
222  bool ImagInteg();
223 
224 public:
227  convention = ComplexOperator::HERMITIAN);
228  /** @brief Create a SesquilinearForm on the FiniteElementSpace @a fes, using
229  the same integrators as the BilinearForms @a bfr and @a bfi .
230 
231  The pointer @a fes is not owned by the newly constructed object.
232 
233  The integrators are copied as pointers and they are not owned by the
234  newly constructed SesquilinearForm. */
237  convention = ComplexOperator::HERMITIAN);
238 
241  convention) { conv = convention; }
242 
243  /// Set the desired assembly level.
244  /** Valid choices are:
245 
246  - AssemblyLevel::LEGACY (default)
247  - AssemblyLevel::FULL
248  - AssemblyLevel::PARTIAL
249  - AssemblyLevel::ELEMENT
250  - AssemblyLevel::NONE
251 
252  This method must be called before assembly. */
253  void SetAssemblyLevel(AssemblyLevel assembly_level)
254  {
255  blfr->SetAssemblyLevel(assembly_level);
256  blfi->SetAssemblyLevel(assembly_level);
257  }
258 
259  BilinearForm & real() { return *blfr; }
260  BilinearForm & imag() { return *blfi; }
261  const BilinearForm & real() const { return *blfr; }
262  const BilinearForm & imag() const { return *blfi; }
263 
264  /// Adds new Domain Integrator.
266  BilinearFormIntegrator *bfi_imag);
267 
268  /// Adds new Domain Integrator, restricted to the given attributes.
270  BilinearFormIntegrator *bfi_imag,
271  Array<int> &elem_marker);
272 
273  /// Adds new Boundary Integrator.
275  BilinearFormIntegrator *bfi_imag);
276 
277  /// Adds new Boundary Integrator, restricted to specific boundary attributes.
279  BilinearFormIntegrator *bfi_imag,
280  Array<int> &bdr_marker);
281 
282  /// Adds new interior Face Integrator. Assumes ownership of @a bfi.
284  BilinearFormIntegrator *bfi_imag);
285 
286  /// Adds new boundary Face Integrator. Assumes ownership of @a bfi.
288  BilinearFormIntegrator *bfi_imag);
289 
290  /** @brief Adds new boundary Face Integrator, restricted to specific boundary
291  attributes.
292 
293  Assumes ownership of @a bfi.
294 
295  The array @a bdr_marker is stored internally as a pointer to the given
296  Array<int> object. */
298  BilinearFormIntegrator *bfi_imag,
299  Array<int> &bdr_marker);
300 
301  /// Assemble the local matrix
302  void Assemble(int skip_zeros = 1);
303 
304  /// Finalizes the matrix initialization.
305  void Finalize(int skip_zeros = 1);
306 
307  /// Returns the matrix assembled on the true dofs, i.e. P^t A P.
308  /** The returned matrix has to be deleted by the caller. */
310 
311  /// Return the parallel FE space associated with the ParBilinearForm.
312  FiniteElementSpace *FESpace() const { return blfr->FESpace(); }
313 
314  void FormLinearSystem(const Array<int> &ess_tdof_list, Vector &x, Vector &b,
315  OperatorHandle &A, Vector &X, Vector &B,
316  int copy_interior = 0);
317 
318  void FormSystemMatrix(const Array<int> &ess_tdof_list,
319  OperatorHandle &A);
320 
321  /** Call this method after solving a linear system constructed using the
322  FormLinearSystem method to recover the solution as a ParGridFunction-size
323  vector in x. Use the same arguments as in the FormLinearSystem call. */
324  virtual void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x);
325 
326  virtual void Update(FiniteElementSpace *nfes = NULL);
327 
328  /// Sets diagonal policy used upon construction of the linear system
330 
331  /// Returns the diagonal policy of the sesquilinear form
332  Matrix::DiagonalPolicy GetDiagonalPolicy() const {return diag_policy;}
333 
334  virtual ~SesquilinearForm();
335 };
336 
337 #ifdef MFEM_USE_MPI
338 
339 /// Class for parallel complex-valued grid function - real + imaginary part
340 /// Vector with associated parallel FE space.
342 {
343 private:
344 
345  ParGridFunction * pgfr;
346  ParGridFunction * pgfi;
347 
348 protected:
349  void Destroy() { delete pgfr; delete pgfi; }
350 
351 public:
352 
353  /** @brief Construct a ParComplexGridFunction associated with the
354  ParFiniteElementSpace @a *pf. */
356 
357  void Update();
358 
359  /// Assign constant values to the ParComplexGridFunction data.
360  ParComplexGridFunction &operator=(const std::complex<double> & value)
361  { *pgfr = value.real(); *pgfi = value.imag(); return *this; }
362 
363  virtual void ProjectCoefficient(Coefficient &real_coeff,
364  Coefficient &imag_coeff);
365  virtual void ProjectCoefficient(VectorCoefficient &real_vcoeff,
366  VectorCoefficient &imag_vcoeff);
367 
368  virtual void ProjectBdrCoefficient(Coefficient &real_coeff,
369  Coefficient &imag_coeff,
370  Array<int> &attr);
371  virtual void ProjectBdrCoefficientNormal(VectorCoefficient &real_coeff,
372  VectorCoefficient &imag_coeff,
373  Array<int> &attr);
374  virtual void ProjectBdrCoefficientTangent(VectorCoefficient &real_coeff,
375  VectorCoefficient &imag_coeff,
376  Array<int> &attr);
377 
378  void Distribute(const Vector *tv);
379  void Distribute(const Vector &tv) { Distribute(&tv); }
380 
381  /// Returns the vector restricted to the true dofs.
382  void ParallelProject(Vector &tv) const;
383 
384  FiniteElementSpace *FESpace() { return pgfr->FESpace(); }
385  const FiniteElementSpace *FESpace() const { return pgfr->FESpace(); }
386 
388  const ParFiniteElementSpace *ParFESpace() const { return pgfr->ParFESpace(); }
389 
390  ParGridFunction & real() { return *pgfr; }
391  ParGridFunction & imag() { return *pgfi; }
392  const ParGridFunction & real() const { return *pgfr; }
393  const ParGridFunction & imag() const { return *pgfi; }
394 
395  /// Update the memory location of the real and imaginary ParGridFunction @a
396  /// pgfr and @a pgfi to match the ParComplexGridFunction.
397  void Sync() { pgfr->SyncMemory(*this); pgfi->SyncMemory(*this); }
398 
399  /// Update the alias memory location of the real and imaginary
400  /// ParGridFunction @a pgfr and @a pgfi to match the ParComplexGridFunction.
401  void SyncAlias() { pgfr->SyncAliasMemory(*this); pgfi->SyncAliasMemory(*this); }
402 
403 
404  virtual double ComputeL2Error(Coefficient &exsolr, Coefficient &exsoli,
405  const IntegrationRule *irs[] = NULL) const
406  {
407  double err_r = pgfr->ComputeL2Error(exsolr, irs);
408  double err_i = pgfi->ComputeL2Error(exsoli, irs);
409  return sqrt(err_r * err_r + err_i * err_i);
410  }
411 
412  virtual double ComputeL2Error(VectorCoefficient &exsolr,
413  VectorCoefficient &exsoli,
414  const IntegrationRule *irs[] = NULL,
415  Array<int> *elems = NULL) const
416  {
417  double err_r = pgfr->ComputeL2Error(exsolr, irs, elems);
418  double err_i = pgfi->ComputeL2Error(exsoli, irs, elems);
419  return sqrt(err_r * err_r + err_i * err_i);
420  }
421 
422 
423  /// Destroys grid function.
425 
426 };
427 
428 /** Class for a complex-valued, parallel linear form
429 
430  The @a convention argument in the class's constructor is documented in the
431  mfem::ComplexOperator class found in linalg/complex_operator.hpp.
432 
433  When supplying integrators to the ParComplexLinearForm either the real or
434  imaginary integrator can be NULL. This indicates that the corresponding
435  portion of the complex-valued field is equal to zero.
436  */
438 {
439 private:
441 
442 protected:
445 
447 
448 public:
449 
452  convention = ComplexOperator::HERMITIAN);
453 
454  /** @brief Create a ParComplexLinearForm on the ParFiniteElementSpace @a pf,
455  using the same integrators as the LinearForms @a plf_r (real) and
456  @a plf_i (imag).
457 
458  The pointer @a fes is not owned by the newly constructed object.
459 
460  The integrators are copied as pointers and they are not owned by the newly
461  constructed ParComplexLinearForm. */
463  ParLinearForm *plf_i,
465  convention = ComplexOperator::HERMITIAN);
466 
467  virtual ~ParComplexLinearForm();
468 
471  convention) { conv = convention; }
472 
473  /// Adds new Domain Integrator.
475  LinearFormIntegrator *lfi_imag);
476 
477  /// Adds new Domain Integrator, restricted to specific attributes.
479  LinearFormIntegrator *lfi_imag,
480  Array<int> &elem_attr_marker);
481 
482  /// Adds new Boundary Integrator.
484  LinearFormIntegrator *lfi_imag);
485 
486  /** @brief Add new Boundary Integrator, restricted to the given boundary
487  attributes.
488 
489  Assumes ownership of @a lfi_real and @a lfi_imag.
490 
491  The array @a bdr_attr_marker is stored internally as a pointer to the
492  given Array<int> object. */
494  LinearFormIntegrator *lfi_imag,
495  Array<int> &bdr_attr_marker);
496 
497  /// Adds new Boundary Face Integrator. Assumes ownership of @a lfi.
499  LinearFormIntegrator *lfi_imag);
500 
501  /** @brief Add new Boundary Face Integrator, restricted to the given boundary
502  attributes.
503 
504  Assumes ownership of @a lfi_real and @a lfi_imag.
505 
506  The array @a bdr_attr_marker is stored internally as a pointer to the
507  given Array<int> object. */
509  LinearFormIntegrator *lfi_imag,
510  Array<int> &bdr_attr_marker);
511 
513 
514  ParLinearForm & real() { return *plfr; }
515  ParLinearForm & imag() { return *plfi; }
516  const ParLinearForm & real() const { return *plfr; }
517  const ParLinearForm & imag() const { return *plfi; }
518 
519  /// Update the memory location of the real and imaginary ParLinearForm @a lfr
520  /// and @a lfi to match the ParComplexLinearForm.
521  void Sync() { plfr->SyncMemory(*this); plfi->SyncMemory(*this); }
522 
523  /// Update the alias memory location of the real and imaginary ParLinearForm
524  /// @a plfr and @a plfi to match the ParComplexLinearForm.
525  void SyncAlias() { plfr->SyncAliasMemory(*this); plfi->SyncAliasMemory(*this); }
526 
527  void Update(ParFiniteElementSpace *pf = NULL);
528 
529  /// Assembles the linear form i.e. sums over all domain/bdr integrators.
530  void Assemble();
531 
532  /// Assemble the vector on the true dofs, i.e. P^t v.
533  void ParallelAssemble(Vector &tv);
534 
535  /// Returns the vector assembled on the true dofs, i.e. P^t v.
537 
538  std::complex<double> operator()(const ParComplexGridFunction &gf) const;
539 
540 };
541 
542 /** Class for a parallel sesquilinear form
543 
544  A sesquilinear form is a generalization of a bilinear form to complex-valued
545  fields. Sesquilinear forms are linear in the second argument but the
546  first argument involves a complex conjugate in the sense that:
547 
548  a(alpha u, beta v) = conj(alpha) beta a(u, v)
549 
550  The @a convention argument in the class's constructor is documented in the
551  mfem::ComplexOperator class found in linalg/complex_operator.hpp.
552 
553  When supplying integrators to the ParSesquilinearForm either the real or
554  imaginary integrator can be NULL. This indicates that the corresponding
555  portion of the complex-valued material coefficient is equal to zero.
556 */
558 {
559 private:
561 
562  ParBilinearForm *pblfr;
563  ParBilinearForm *pblfi;
564 
565  /* These methods check if the real/imag parts of the sesqulinear form are not
566  empty */
567  bool RealInteg();
568  bool ImagInteg();
569 
570 public:
573  convention = ComplexOperator::HERMITIAN);
574 
575  /** @brief Create a ParSesquilinearForm on the ParFiniteElementSpace @a pf,
576  using the same integrators as the ParBilinearForms @a pbfr and @a pbfi .
577 
578  The pointer @a pf is not owned by the newly constructed object.
579 
580  The integrators are copied as pointers and they are not owned by the
581  newly constructed ParSesquilinearForm. */
583  ParBilinearForm *pbfi,
585  convention = ComplexOperator::HERMITIAN);
586 
589  convention) { conv = convention; }
590 
591  /// Set the desired assembly level.
592  /** Valid choices are:
593 
594  - AssemblyLevel::LEGACY (default)
595  - AssemblyLevel::FULL
596  - AssemblyLevel::PARTIAL
597  - AssemblyLevel::ELEMENT
598  - AssemblyLevel::NONE
599 
600  This method must be called before assembly. */
601  void SetAssemblyLevel(AssemblyLevel assembly_level)
602  {
603  pblfr->SetAssemblyLevel(assembly_level);
604  pblfi->SetAssemblyLevel(assembly_level);
605  }
606 
607  ParBilinearForm & real() { return *pblfr; }
608  ParBilinearForm & imag() { return *pblfi; }
609  const ParBilinearForm & real() const { return *pblfr; }
610  const ParBilinearForm & imag() const { return *pblfi; }
611 
612  /// Adds new Domain Integrator.
614  BilinearFormIntegrator *bfi_imag);
615 
616  /// Adds new Domain Integrator, restricted to specific attributes.
618  BilinearFormIntegrator *bfi_imag,
619  Array<int> &elem_marker);
620 
621  /// Adds new Boundary Integrator.
623  BilinearFormIntegrator *bfi_imag);
624 
625  /** @brief Adds new boundary Integrator, restricted to specific boundary
626  attributes.
627 
628  Assumes ownership of @a bfi.
629 
630  The array @a bdr_marker is stored internally as a pointer to the given
631  Array<int> object. */
633  BilinearFormIntegrator *bfi_imag,
634  Array<int> &bdr_marker);
635 
636  /// Adds new interior Face Integrator. Assumes ownership of @a bfi.
638  BilinearFormIntegrator *bfi_imag);
639 
640  /// Adds new boundary Face Integrator. Assumes ownership of @a bfi.
642  BilinearFormIntegrator *bfi_imag);
643 
644  /** @brief Adds new boundary Face Integrator, restricted to specific boundary
645  attributes.
646 
647  Assumes ownership of @a bfi.
648 
649  The array @a bdr_marker is stored internally as a pointer to the given
650  Array<int> object. */
652  BilinearFormIntegrator *bfi_imag,
653  Array<int> &bdr_marker);
654 
655  /// Assemble the local matrix
656  void Assemble(int skip_zeros = 1);
657 
658  /// Finalizes the matrix initialization.
659  void Finalize(int skip_zeros = 1);
660 
661  /// Returns the matrix assembled on the true dofs, i.e. P^t A P.
662  /** The returned matrix has to be deleted by the caller. */
664 
665  /// Return the parallel FE space associated with the ParBilinearForm.
666  ParFiniteElementSpace *ParFESpace() const { return pblfr->ParFESpace(); }
667 
668  void FormLinearSystem(const Array<int> &ess_tdof_list, Vector &x, Vector &b,
669  OperatorHandle &A, Vector &X, Vector &B,
670  int copy_interior = 0);
671 
672  void FormSystemMatrix(const Array<int> &ess_tdof_list,
673  OperatorHandle &A);
674 
675  /** Call this method after solving a linear system constructed using the
676  FormLinearSystem method to recover the solution as a ParGridFunction-size
677  vector in x. Use the same arguments as in the FormLinearSystem call. */
678  virtual void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x);
679 
680  virtual void Update(FiniteElementSpace *nfes = NULL);
681 
682  virtual ~ParSesquilinearForm();
683 };
684 
685 #endif // MFEM_USE_MPI
686 
687 }
688 
689 #endif // MFEM_COMPLEX_FEM
void SetConvention(const ComplexOperator::Convention &convention)
ComplexOperator::Convention GetConvention() const
const LinearForm & real() const
ParFiniteElementSpace * ParFESpace() const
Return the parallel FE space associated with the ParBilinearForm.
AssemblyLevel
Enumeration defining the assembly level for bilinear and nonlinear form classes derived from Operator...
virtual void ProjectBdrCoefficientTangent(VectorCoefficient &real_coeff, VectorCoefficient &imag_coeff, Array< int > &attr)
Class for an integration rule - an Array of IntegrationPoint.
Definition: intrules.hpp:96
Class for grid function - Vector with associated FE space.
Definition: gridfunc.hpp:30
void SetDiagonalPolicy(mfem::Matrix::DiagonalPolicy dpolicy)
Sets diagonal policy used upon construction of the linear system.
void FormSystemMatrix(const Array< int > &ess_tdof_list, OperatorHandle &A)
Base class for vector Coefficients that optionally depend on time and space.
virtual void ProjectBdrCoefficient(Coefficient &real_coeff, Coefficient &imag_coeff, Array< int > &attr)
const LinearForm & imag() const
ParBilinearForm & real()
const ParLinearForm & imag() const
FiniteElementSpace * FESpace()
Return the FE space associated with the BilinearForm.
std::complex< double > operator()(const ComplexGridFunction &gf) const
void SyncMemory(const Vector &v) const
Update the memory location of the vector to match v.
Definition: vector.hpp:235
ComplexSparseMatrix * AssembleComplexSparseMatrix()
Returns the matrix assembled on the true dofs, i.e. P^t A P.
const ParLinearForm & real() const
virtual void ProjectCoefficient(Coefficient &real_coeff, Coefficient &imag_coeff)
Definition: complex_fem.cpp:88
void SetConvention(const ComplexOperator::Convention &convention)
void SetAssemblyLevel(AssemblyLevel assembly_level)
Set the desired assembly level.
Pointer to an Operator of a specified type.
Definition: handle.hpp:33
void Finalize(int skip_zeros=1)
Finalizes the matrix initialization.
FiniteElementSpace * FESpace()
void AddDomainIntegrator(LinearFormIntegrator *lfi_real, LinearFormIntegrator *lfi_imag)
Adds new Domain Integrator.
const ParGridFunction & real() const
void AddDomainIntegrator(LinearFormIntegrator *lfi_real, LinearFormIntegrator *lfi_imag)
Adds new Domain Integrator.
ParGridFunction & imag()
void SetAssemblyLevel(AssemblyLevel assembly_level)
Set the desired assembly level.
Abstract parallel finite element space.
Definition: pfespace.hpp:28
ParFiniteElementSpace * ParFESpace() const
HypreParVector * ParallelAssemble()
Returns the vector assembled on the true dofs, i.e. P^t v.
ParBilinearForm & imag()
ParComplexLinearForm(ParFiniteElementSpace *pf, ComplexOperator::Convention convention=ComplexOperator::HERMITIAN)
Abstract base class LinearFormIntegrator.
Definition: lininteg.hpp:24
std::complex< double > operator()(const ParComplexGridFunction &gf) const
ParLinearForm & imag()
virtual double ComputeL2Error(VectorCoefficient &exsolr, VectorCoefficient &exsoli, const IntegrationRule *irs[]=NULL, Array< int > *elems=NULL) const
virtual void ProjectBdrCoefficientNormal(VectorCoefficient &real_coeff, VectorCoefficient &imag_coeff, Array< int > &attr)
FiniteElementSpace * FESpace()
Read+write access to the associated FiniteElementSpace.
Definition: linearform.hpp:129
std::function< double(const Vector &)> f(double mass_coeff)
Definition: lor_mms.hpp:30
void FormLinearSystem(const Array< int > &ess_tdof_list, Vector &x, Vector &b, OperatorHandle &A, Vector &X, Vector &B, int copy_interior=0)
ComplexLinearForm(FiniteElementSpace *fes, ComplexOperator::Convention convention=ComplexOperator::HERMITIAN)
ParFiniteElementSpace * ParFESpace()
SesquilinearForm(FiniteElementSpace *fes, ComplexOperator::Convention convention=ComplexOperator::HERMITIAN)
ParSesquilinearForm(ParFiniteElementSpace *pf, ComplexOperator::Convention convention=ComplexOperator::HERMITIAN)
ParComplexGridFunction & operator=(const std::complex< double > &value)
Assign constant values to the ParComplexGridFunction data.
void ParallelProject(Vector &tv) const
Returns the vector restricted to the true dofs.
ComplexHypreParMatrix * ParallelAssemble()
Returns the matrix assembled on the true dofs, i.e. P^t A P.
const FiniteElementSpace * FESpace() const
Definition: complex_fem.hpp:67
Class for parallel linear form.
Definition: plinearform.hpp:26
void SetAssemblyLevel(AssemblyLevel assembly_level)
Set the desired assembly level.
void Distribute(const Vector &tv)
void AddBdrFaceIntegrator(BilinearFormIntegrator *bfi_real, BilinearFormIntegrator *bfi_imag)
Adds new boundary Face Integrator. Assumes ownership of bfi.
ComplexGridFunction(FiniteElementSpace *f)
Construct a ComplexGridFunction associated with the FiniteElementSpace *f.
Definition: complex_fem.cpp:20
double b
Definition: lissajous.cpp:42
const ParBilinearForm & real() const
virtual void ProjectBdrCoefficient(Coefficient &real_coeff, Coefficient &imag_coeff, Array< int > &attr)
virtual ~ComplexGridFunction()
Destroys the grid function.
Definition: complex_fem.hpp:83
void FormSystemMatrix(const Array< int > &ess_tdof_list, OperatorHandle &A)
void Update(ParFiniteElementSpace *pf=NULL)
void SetConvention(const ComplexOperator::Convention &convention)
virtual void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x)
GridFunction & real()
Definition: complex_fem.hpp:69
Set the diagonal value to one.
Definition: operator.hpp:50
FiniteElementSpace * FESpace()
Definition: gridfunc.hpp:691
virtual double ComputeL2Error(Coefficient *exsol[], const IntegrationRule *irs[]=NULL, const Array< int > *elems=NULL) const
Definition: pgridfunc.hpp:285
ComplexOperator::Convention GetConvention() const
BilinearForm & real()
Wrapper for hypre&#39;s parallel vector class.
Definition: hypre.hpp:161
void Distribute(const Vector *tv)
ParFiniteElementSpace * ParFESpace() const
Definition: plinearform.hpp:76
void SetConvention(const ComplexOperator::Convention &convention)
void AddInteriorFaceIntegrator(BilinearFormIntegrator *bfi_real, BilinearFormIntegrator *bfi_imag)
Adds new interior Face Integrator. Assumes ownership of bfi.
virtual void ProjectBdrCoefficientNormal(VectorCoefficient &real_coeff, VectorCoefficient &imag_coeff, Array< int > &attr)
Abstract base class BilinearFormIntegrator.
Definition: bilininteg.hpp:24
void Finalize(int skip_zeros=1)
Finalizes the matrix initialization.
const BilinearForm & real() const
Specialization of the ComplexOperator built from a pair of Sparse Matrices.
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:219
GridFunction & imag()
Definition: complex_fem.hpp:70
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
Definition: coefficient.hpp:41
BilinearForm & imag()
ParFiniteElementSpace * ParFESpace() const
Definition: pgridfunc.hpp:113
const GridFunction & real() const
Definition: complex_fem.hpp:71
void AddBoundaryIntegrator(BilinearFormIntegrator *bfi_real, BilinearFormIntegrator *bfi_imag)
Adds new Boundary Integrator.
HYPRE_Int HYPRE_BigInt
virtual ~ParComplexGridFunction()
Destroys grid function.
const FiniteElementSpace * FESpace() const
const ParGridFunction & imag() const
void AddBoundaryIntegrator(LinearFormIntegrator *lfi_real, LinearFormIntegrator *lfi_imag)
Adds new Boundary Integrator.
ParComplexGridFunction(ParFiniteElementSpace *pf)
Construct a ParComplexGridFunction associated with the ParFiniteElementSpace *pf. ...
ParLinearForm & real()
virtual void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x)
void AddDomainIntegrator(BilinearFormIntegrator *bfi_real, BilinearFormIntegrator *bfi_imag)
Adds new Domain Integrator.
ParFiniteElementSpace * ParFESpace() const
Return the parallel FE space associated with the ParBilinearForm.
FiniteElementSpace * FESpace()
Definition: complex_fem.hpp:66
A "square matrix" operator for the associated FE space and BLFIntegrators The sum of all the BLFInteg...
void Assemble()
Assembles the linear form i.e. sums over all domain/bdr integrators.
void AddBdrFaceIntegrator(LinearFormIntegrator *lfi_real, LinearFormIntegrator *lfi_imag)
Adds new Boundary Face Integrator. Assumes ownership of lfi.
virtual void ProjectBdrCoefficientTangent(VectorCoefficient &real_coeff, VectorCoefficient &imag_coeff, Array< int > &attr)
void Assemble(int skip_zeros=1)
Assemble the local matrix.
void FormLinearSystem(const Array< int > &ess_tdof_list, Vector &x, Vector &b, OperatorHandle &A, Vector &X, Vector &B, int copy_interior=0)
void SyncAliasMemory(const Vector &v) const
Update the alias memory location of the vector to match v.
Definition: vector.hpp:238
void Assemble(int skip_zeros=1)
Assemble the local matrix.
virtual double ComputeL2Error(Coefficient &exsolr, Coefficient &exsoli, const IntegrationRule *irs[]=NULL) const
DiagonalPolicy
Defines operator diagonal policy upon elimination of rows and/or columns.
Definition: operator.hpp:47
Specialization of the ComplexOperator built from a pair of HypreParMatrices.
Native convention for Hermitian operators.
void AddDomainIntegrator(BilinearFormIntegrator *bfi_real, BilinearFormIntegrator *bfi_imag)
Adds new Domain Integrator.
const ParFiniteElementSpace * ParFESpace() const
Matrix::DiagonalPolicy GetDiagonalPolicy() const
Returns the diagonal policy of the sesquilinear form.
Class for parallel bilinear form.
const ParBilinearForm & imag() const
Vector data type.
Definition: vector.hpp:58
ComplexOperator::Convention GetConvention() const
const BilinearForm & imag() const
virtual void Update(FiniteElementSpace *nfes=NULL)
void AddBdrFaceIntegrator(BilinearFormIntegrator *bfi_real, BilinearFormIntegrator *bfi_imag)
Adds new boundary Face Integrator. Assumes ownership of bfi.
ParGridFunction & real()
Vector with associated FE space and LinearFormIntegrators.
Definition: linearform.hpp:24
ComplexGridFunction & operator=(const std::complex< double > &value)
Assign constant values to the ComplexGridFunction data.
Definition: complex_fem.hpp:48
void AddBdrFaceIntegrator(LinearFormIntegrator *lfi_real, LinearFormIntegrator *lfi_imag)
Adds new Boundary Face Integrator. Assumes ownership of lfi.
Class for parallel grid function.
Definition: pgridfunc.hpp:32
void Assemble()
Assembles the linear form i.e. sums over all domain/bdr integrators.
void AddBoundaryIntegrator(BilinearFormIntegrator *bfi_real, BilinearFormIntegrator *bfi_imag)
Adds new Boundary Integrator.
FiniteElementSpace * FESpace() const
Return the parallel FE space associated with the ParBilinearForm.
virtual void Update(FiniteElementSpace *nfes=NULL)
FiniteElementSpace * FESpace() const
void AddInteriorFaceIntegrator(BilinearFormIntegrator *bfi_real, BilinearFormIntegrator *bfi_imag)
Adds new interior Face Integrator. Assumes ownership of bfi.
const GridFunction & imag() const
Definition: complex_fem.hpp:72
void AddBoundaryIntegrator(LinearFormIntegrator *lfi_real, LinearFormIntegrator *lfi_imag)
Adds new Boundary Integrator.
virtual void ProjectCoefficient(Coefficient &real_coeff, Coefficient &imag_coeff)
ComplexOperator::Convention GetConvention() const