MFEM  v4.2.0
Finite element discretization library
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Pages
complex_fem.cpp
Go to the documentation of this file.
1 // Copyright (c) 2010-2020, 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 "complex_fem.hpp"
13 #include "../general/forall.hpp"
14 
15 using namespace std;
16 
17 namespace mfem
18 {
19 
20 ComplexGridFunction::ComplexGridFunction(FiniteElementSpace *fes)
21  : Vector(2*(fes->GetVSize()))
22 {
23  UseDevice(true);
24  this->Vector::operator=(0.0);
25 
26  gfr = new GridFunction();
27  gfr->MakeRef(fes, *this, 0);
28 
29  gfi = new GridFunction();
30  gfi->MakeRef(fes, *this, fes->GetVSize());
31 }
32 
33 void
35 {
36  FiniteElementSpace *fes = gfr->FESpace();
37  const int vsize = fes->GetVSize();
38 
39  const Operator *T = fes->GetUpdateOperator();
40  if (T)
41  {
42  // Update the individual GridFunction objects. This will allocate new data
43  // arrays for each GridFunction.
44  gfr->Update();
45  gfi->Update();
46 
47  // Our data array now contains old data as well as being the wrong size so
48  // reallocate it.
49  UseDevice(true);
50  this->SetSize(2 * vsize);
51  this->Vector::operator=(0.0);
52 
53  // Create temporary vectors which point to the new data array
54  Vector gf_r; gf_r.MakeRef(*this, 0, vsize);
55  Vector gf_i; gf_i.MakeRef(*this, vsize, vsize);
56 
57  // Copy the updated GridFunctions into the new data array
58  gf_r = *gfr;
59  gf_i = *gfi;
60  gf_r.SyncAliasMemory(*this);
61  gf_i.SyncAliasMemory(*this);
62 
63  // Replace the individual data arrays with pointers into the new data
64  // array
65  gfr->MakeRef(*this, 0, vsize);
66  gfi->MakeRef(*this, vsize, vsize);
67  }
68  else
69  {
70  // The existing data will not be transferred to the new GridFunctions so
71  // delete it and allocate a new array
72  UseDevice(true);
73  this->SetSize(2 * vsize);
74  this->Vector::operator=(0.0);
75 
76  // Point the individual GridFunctions to the new data array
77  gfr->MakeRef(*this, 0, vsize);
78  gfi->MakeRef(*this, vsize, vsize);
79 
80  // These updates will only set the proper 'sequence' value within the
81  // individual GridFunction objects because their sizes are already correct
82  gfr->Update();
83  gfi->Update();
84  }
85 }
86 
87 void
89  Coefficient &imag_coeff)
90 {
91  gfr->SyncMemory(*this);
92  gfi->SyncMemory(*this);
93  gfr->ProjectCoefficient(real_coeff);
94  gfi->ProjectCoefficient(imag_coeff);
95  gfr->SyncAliasMemory(*this);
96  gfi->SyncAliasMemory(*this);
97 }
98 
99 void
101  VectorCoefficient &imag_vcoeff)
102 {
103  gfr->SyncMemory(*this);
104  gfi->SyncMemory(*this);
105  gfr->ProjectCoefficient(real_vcoeff);
106  gfi->ProjectCoefficient(imag_vcoeff);
107  gfr->SyncAliasMemory(*this);
108  gfi->SyncAliasMemory(*this);
109 }
110 
111 void
113  Coefficient &imag_coeff,
114  Array<int> &attr)
115 {
116  gfr->SyncMemory(*this);
117  gfi->SyncMemory(*this);
118  gfr->ProjectBdrCoefficient(real_coeff, attr);
119  gfi->ProjectBdrCoefficient(imag_coeff, attr);
120  gfr->SyncAliasMemory(*this);
121  gfi->SyncAliasMemory(*this);
122 }
123 
124 void
126  VectorCoefficient &imag_vcoeff,
127  Array<int> &attr)
128 {
129  gfr->SyncMemory(*this);
130  gfi->SyncMemory(*this);
131  gfr->ProjectBdrCoefficientNormal(real_vcoeff, attr);
132  gfi->ProjectBdrCoefficientNormal(imag_vcoeff, attr);
133  gfr->SyncAliasMemory(*this);
134  gfi->SyncAliasMemory(*this);
135 }
136 
137 void
139  &real_vcoeff,
141  &imag_vcoeff,
142  Array<int> &attr)
143 {
144  gfr->SyncMemory(*this);
145  gfi->SyncMemory(*this);
146  gfr->ProjectBdrCoefficientTangent(real_vcoeff, attr);
147  gfi->ProjectBdrCoefficientTangent(imag_vcoeff, attr);
148  gfr->SyncAliasMemory(*this);
149  gfi->SyncAliasMemory(*this);
150 }
151 
152 
154  ComplexOperator::Convention convention)
155  : Vector(2*(fes->GetVSize())),
156  conv(convention)
157 {
158  UseDevice(true);
159  this->Vector::operator=(0.0);
160 
161  lfr = new LinearForm();
162  lfr->MakeRef(fes, *this, 0);
163 
164  lfi = new LinearForm();
165  lfi->MakeRef(fes, *this, fes->GetVSize());
166 }
167 
169  LinearForm *lf_r, LinearForm *lf_i,
170  ComplexOperator::Convention convention)
171  : Vector(2*(fes->GetVSize())),
172  conv(convention)
173 {
174  UseDevice(true);
175  this->Vector::operator=(0.0);
176 
177  lfr = new LinearForm(fes, lf_r);
178  lfi = new LinearForm(fes, lf_i);
179 
180  lfr->MakeRef(fes, *this, 0);
181  lfi->MakeRef(fes, *this, fes->GetVSize());
182 }
183 
185 {
186  delete lfr;
187  delete lfi;
188 }
189 
190 void
192  LinearFormIntegrator *lfi_imag)
193 {
194  if ( lfi_real ) { lfr->AddDomainIntegrator(lfi_real); }
195  if ( lfi_imag ) { lfi->AddDomainIntegrator(lfi_imag); }
196 }
197 
198 void
200  LinearFormIntegrator *lfi_imag)
201 {
202  if ( lfi_real ) { lfr->AddBoundaryIntegrator(lfi_real); }
203  if ( lfi_imag ) { lfi->AddBoundaryIntegrator(lfi_imag); }
204 }
205 
206 void
208  LinearFormIntegrator *lfi_imag,
209  Array<int> &bdr_attr_marker)
210 {
211  if ( lfi_real ) { lfr->AddBoundaryIntegrator(lfi_real, bdr_attr_marker); }
212  if ( lfi_imag ) { lfi->AddBoundaryIntegrator(lfi_imag, bdr_attr_marker); }
213 }
214 
215 void
217  LinearFormIntegrator *lfi_imag)
218 {
219  if ( lfi_real ) { lfr->AddBdrFaceIntegrator(lfi_real); }
220  if ( lfi_imag ) { lfi->AddBdrFaceIntegrator(lfi_imag); }
221 }
222 
223 void
225  LinearFormIntegrator *lfi_imag,
226  Array<int> &bdr_attr_marker)
227 {
228  if ( lfi_real ) { lfr->AddBdrFaceIntegrator(lfi_real, bdr_attr_marker); }
229  if ( lfi_imag ) { lfi->AddBdrFaceIntegrator(lfi_imag, bdr_attr_marker); }
230 }
231 
232 void
234 {
235  FiniteElementSpace *fes = lfr->FESpace();
236  this->Update(fes);
237 }
238 
239 void
241 {
242  UseDevice(true);
243  SetSize(2 * fes->GetVSize());
244  this->Vector::operator=(0.0);
245 
246  lfr->MakeRef(fes, *this, 0);
247  lfi->MakeRef(fes, *this, fes->GetVSize());
248 }
249 
250 void
252 {
253  lfr->SyncMemory(*this);
254  lfi->SyncMemory(*this);
255  lfr->Assemble();
256  lfi->Assemble();
257  if (conv == ComplexOperator::BLOCK_SYMMETRIC) { *lfi *= -1.0; }
258  lfr->SyncAliasMemory(*this);
259  lfi->SyncAliasMemory(*this);
260 }
261 
262 complex<double>
264 {
265  double s = (conv == ComplexOperator::HERMITIAN) ? 1.0 : -1.0;
266  lfr->SyncMemory(*this);
267  lfi->SyncMemory(*this);
268  return complex<double>((*lfr)(gf.real()) - s * (*lfi)(gf.imag()),
269  (*lfr)(gf.imag()) + s * (*lfi)(gf.real()));
270 }
271 
272 
273 bool SesquilinearForm::RealInteg()
274 {
275  int nint = blfr->GetFBFI()->Size() + blfr->GetDBFI()->Size() +
276  blfr->GetBBFI()->Size() + blfr->GetBFBFI()->Size();
277  return (nint != 0);
278 }
279 
280 bool SesquilinearForm::ImagInteg()
281 {
282  int nint = blfi->GetFBFI()->Size() + blfi->GetDBFI()->Size() +
283  blfi->GetBBFI()->Size() + blfi->GetBFBFI()->Size();
284  return (nint != 0);
285 }
286 
288  ComplexOperator::Convention convention)
289  : conv(convention),
290  blfr(new BilinearForm(f)),
291  blfi(new BilinearForm(f))
292 {}
293 
295  BilinearForm *bfr, BilinearForm *bfi,
296  ComplexOperator::Convention convention)
297  : conv(convention),
298  blfr(new BilinearForm(f,bfr)),
299  blfi(new BilinearForm(f,bfi))
300 {}
301 
303 {
304  diag_policy = dpolicy;
305 }
306 
308 {
309  delete blfr;
310  delete blfi;
311 }
312 
314  BilinearFormIntegrator *bfi_imag)
315 {
316  if (bfi_real) { blfr->AddDomainIntegrator(bfi_real); }
317  if (bfi_imag) { blfi->AddDomainIntegrator(bfi_imag); }
318 }
319 
320 void
322  BilinearFormIntegrator *bfi_imag)
323 {
324  if (bfi_real) { blfr->AddBoundaryIntegrator(bfi_real); }
325  if (bfi_imag) { blfi->AddBoundaryIntegrator(bfi_imag); }
326 }
327 
328 void
330  BilinearFormIntegrator *bfi_imag,
331  Array<int> & bdr_marker)
332 {
333  if (bfi_real) { blfr->AddBoundaryIntegrator(bfi_real, bdr_marker); }
334  if (bfi_imag) { blfi->AddBoundaryIntegrator(bfi_imag, bdr_marker); }
335 }
336 
337 void
339  BilinearFormIntegrator *bfi_imag)
340 {
341  if (bfi_real) { blfr->AddInteriorFaceIntegrator(bfi_real); }
342  if (bfi_imag) { blfi->AddInteriorFaceIntegrator(bfi_imag); }
343 }
344 
346  BilinearFormIntegrator *bfi_imag)
347 {
348  if (bfi_real) { blfr->AddBdrFaceIntegrator(bfi_real); }
349  if (bfi_imag) { blfi->AddBdrFaceIntegrator(bfi_imag); }
350 }
351 
353  BilinearFormIntegrator *bfi_imag,
354  Array<int> &bdr_marker)
355 {
356  if (bfi_real) { blfr->AddBdrFaceIntegrator(bfi_real, bdr_marker); }
357  if (bfi_imag) { blfi->AddBdrFaceIntegrator(bfi_imag, bdr_marker); }
358 }
359 
360 void
362 {
363  blfr->Assemble(skip_zeros);
364  blfi->Assemble(skip_zeros);
365 }
366 
367 void
369 {
370  blfr->Finalize(skip_zeros);
371  blfi->Finalize(skip_zeros);
372 }
373 
376 {
377  return new ComplexSparseMatrix(&blfr->SpMat(),
378  &blfi->SpMat(),
379  false, false, conv);
380 }
381 
382 void
384  Vector &x, Vector &b,
385  OperatorHandle &A,
386  Vector &X, Vector &B,
387  int ci)
388 {
389  FiniteElementSpace *fes = blfr->FESpace();
390  const int vsize = fes->GetVSize();
391 
392  // Allocate temporary vector
393  Vector b_0;
394  b_0.UseDevice(true);
395  b_0.SetSize(vsize);
396  b_0 = 0.0;
397 
398  // Extract the real and imaginary parts of the input vectors
399  MFEM_ASSERT(x.Size() == 2 * vsize, "Input GridFunction of incorrect size!");
400  x.Read();
401  Vector x_r; x_r.MakeRef(x, 0, vsize);
402  Vector x_i; x_i.MakeRef(x, vsize, vsize);
403 
404  MFEM_ASSERT(b.Size() == 2 * vsize, "Input LinearForm of incorrect size!");
405  b.Read();
406  Vector b_r; b_r.MakeRef(b, 0, vsize);
407  Vector b_i; b_i.MakeRef(b, vsize, vsize);
408 
409  if (conv == ComplexOperator::BLOCK_SYMMETRIC) { b_i *= -1.0; }
410 
411  const int tvsize = fes->GetTrueVSize();
412  OperatorHandle A_r, A_i;
413 
414  X.UseDevice(true);
415  X.SetSize(2 * tvsize);
416  X = 0.0;
417 
418  B.UseDevice(true);
419  B.SetSize(2 * tvsize);
420  B = 0.0;
421 
422  Vector X_r; X_r.MakeRef(X, 0, tvsize);
423  Vector X_i; X_i.MakeRef(X, tvsize, tvsize);
424  Vector B_r; B_r.MakeRef(B, 0, tvsize);
425  Vector B_i; B_i.MakeRef(B, tvsize, tvsize);
426 
427  Vector X_0, B_0;
428 
429  if (RealInteg())
430  {
431  blfr->SetDiagonalPolicy(diag_policy);
432 
433  b_0 = b_r;
434  blfr->FormLinearSystem(ess_tdof_list, x_r, b_0, A_r, X_0, B_0, ci);
435  X_r = X_0; B_r = B_0;
436 
437  b_0 = b_i;
438  blfr->FormLinearSystem(ess_tdof_list, x_i, b_0, A_r, X_0, B_0, ci);
439  X_i = X_0; B_i = B_0;
440 
441  if (ImagInteg())
442  {
443  blfi->SetDiagonalPolicy(mfem::Matrix::DiagonalPolicy::DIAG_ZERO);
444 
445  b_0 = 0.0;
446  blfi->FormLinearSystem(ess_tdof_list, x_i, b_0, A_i, X_0, B_0, false);
447  B_r -= B_0;
448 
449  b_0 = 0.0;
450  blfi->FormLinearSystem(ess_tdof_list, x_r, b_0, A_i, X_0, B_0, false);
451  B_i += B_0;
452  }
453  }
454  else if (ImagInteg())
455  {
456  blfi->SetDiagonalPolicy(diag_policy);
457 
458  b_0 = b_i;
459  blfi->FormLinearSystem(ess_tdof_list, x_r, b_0, A_i, X_0, B_0, ci);
460  X_r = X_0; B_i = B_0;
461 
462  b_0 = b_r; b_0 *= -1.0;
463  blfi->FormLinearSystem(ess_tdof_list, x_i, b_0, A_i, X_0, B_0, ci);
464  X_i = X_0; B_r = B_0; B_r *= -1.0;
465  }
466  else
467  {
468  MFEM_ABORT("Real and Imaginary part of the Sesquilinear form are empty");
469  }
470 
471  if (RealInteg() && ImagInteg())
472  {
473  // Modify RHS and offdiagonal blocks (imaginary parts of the matrix) to
474  // conform with standard essential BC treatment
475  if (A_i.Is<ConstrainedOperator>())
476  {
477  const int n = ess_tdof_list.Size();
478  auto d_B_r = B_r.Write();
479  auto d_B_i = B_i.Write();
480  auto d_X_r = X_r.Read();
481  auto d_X_i = X_i.Read();
482  auto d_idx = ess_tdof_list.Read();
483  MFEM_FORALL(i, n,
484  {
485  const int j = d_idx[i];
486  d_B_r[j] = d_X_r[j];
487  d_B_i[j] = d_X_i[j];
488  });
490  (mfem::Operator::DiagonalPolicy::DIAG_ZERO);
491  }
492  }
493 
495  {
496  B_i *= -1.0;
497  b_i *= -1.0;
498  }
499 
500  x_r.SyncAliasMemory(x);
501  x_i.SyncAliasMemory(x);
502  b_r.SyncAliasMemory(b);
503  b_i.SyncAliasMemory(b);
504 
505  X_r.SyncAliasMemory(X);
506  X_i.SyncAliasMemory(X);
507  B_r.SyncAliasMemory(B);
508  B_i.SyncAliasMemory(B);
509 
510  // A = A_r + i A_i
511  A.Clear();
512  if ( A_r.Type() == Operator::MFEM_SPARSEMAT ||
513  A_i.Type() == Operator::MFEM_SPARSEMAT )
514  {
515  ComplexSparseMatrix * A_sp =
517  A_i.As<SparseMatrix>(),
518  A_r.OwnsOperator(),
519  A_i.OwnsOperator(),
520  conv);
521  A.Reset<ComplexSparseMatrix>(A_sp, true);
522  }
523  else
524  {
525  ComplexOperator * A_op =
526  new ComplexOperator(A_r.Ptr(),
527  A_i.Ptr(),
528  A_r.OwnsOperator(),
529  A_i.OwnsOperator(),
530  conv);
531  A.Reset<ComplexOperator>(A_op, true);
532  }
533  A_r.SetOperatorOwner(false);
534  A_i.SetOperatorOwner(false);
535 }
536 
537 void
539  OperatorHandle &A)
540 
541 {
542  OperatorHandle A_r, A_i;
543  if (RealInteg())
544  {
545  blfr->SetDiagonalPolicy(diag_policy);
546  blfr->FormSystemMatrix(ess_tdof_list, A_r);
547  }
548  if (ImagInteg())
549  {
550  blfi->SetDiagonalPolicy(RealInteg() ?
551  mfem::Matrix::DiagonalPolicy::DIAG_ZERO :
552  diag_policy);
553  blfi->FormSystemMatrix(ess_tdof_list, A_i);
554  }
555  if (!RealInteg() && !ImagInteg())
556  {
557  MFEM_ABORT("Both Real and Imaginary part of the Sesquilinear form are empty");
558  }
559 
560  if (RealInteg() && ImagInteg())
561  {
562  // Modify offdiagonal blocks (imaginary parts of the matrix) to conform
563  // with standard essential BC treatment
564  if (A_i.Is<ConstrainedOperator>())
565  {
567  (mfem::Operator::DiagonalPolicy::DIAG_ZERO);
568  }
569  }
570 
571  // A = A_r + i A_i
572  A.Clear();
573  if ( A_r.Type() == Operator::MFEM_SPARSEMAT ||
574  A_i.Type() == Operator::MFEM_SPARSEMAT )
575  {
576  ComplexSparseMatrix * A_sp =
578  A_i.As<SparseMatrix>(),
579  A_r.OwnsOperator(),
580  A_i.OwnsOperator(),
581  conv);
582  A.Reset<ComplexSparseMatrix>(A_sp, true);
583  }
584  else
585  {
586  ComplexOperator * A_op =
587  new ComplexOperator(A_r.Ptr(),
588  A_i.Ptr(),
589  A_r.OwnsOperator(),
590  A_i.OwnsOperator(),
591  conv);
592  A.Reset<ComplexOperator>(A_op, true);
593  }
594  A_r.SetOperatorOwner(false);
595  A_i.SetOperatorOwner(false);
596 }
597 
598 void
600  Vector &x)
601 {
602  FiniteElementSpace *fes = blfr->FESpace();
603 
604  const SparseMatrix *P = fes->GetConformingProlongation();
605  if (!P)
606  {
607  x = X;
608  return;
609  }
610 
611  const int vsize = fes->GetVSize();
612  const int tvsize = X.Size() / 2;
613 
614  X.Read();
615  Vector X_r; X_r.MakeRef(const_cast<Vector&>(X), 0, tvsize);
616  Vector X_i; X_i.MakeRef(const_cast<Vector&>(X), tvsize, tvsize);
617 
618  x.Write();
619  Vector x_r; x_r.MakeRef(x, 0, vsize);
620  Vector x_i; x_i.MakeRef(x, vsize, vsize);
621 
622  // Apply conforming prolongation
623  P->Mult(X_r, x_r);
624  P->Mult(X_i, x_i);
625 
626  x_r.SyncAliasMemory(x);
627  x_i.SyncAliasMemory(x);
628 }
629 
630 void
632 {
633  if ( blfr ) { blfr->Update(nfes); }
634  if ( blfi ) { blfi->Update(nfes); }
635 }
636 
637 
638 #ifdef MFEM_USE_MPI
639 
641  : Vector(2*(pfes->GetVSize()))
642 {
643  UseDevice(true);
644  this->Vector::operator=(0.0);
645 
646  pgfr = new ParGridFunction();
647  pgfr->MakeRef(pfes, *this, 0);
648 
649  pgfi = new ParGridFunction();
650  pgfi->MakeRef(pfes, *this, pfes->GetVSize());
651 }
652 
653 void
655 {
656  ParFiniteElementSpace *pfes = pgfr->ParFESpace();
657  const int vsize = pfes->GetVSize();
658 
659  const Operator *T = pfes->GetUpdateOperator();
660  if (T)
661  {
662  // Update the individual GridFunction objects. This will allocate new data
663  // arrays for each GridFunction.
664  pgfr->Update();
665  pgfi->Update();
666 
667  // Our data array now contains old data as well as being the wrong size so
668  // reallocate it.
669  UseDevice(true);
670  this->SetSize(2 * vsize);
671  this->Vector::operator=(0.0);
672 
673  // Create temporary vectors which point to the new data array
674  Vector gf_r; gf_r.MakeRef(*this, 0, vsize);
675  Vector gf_i; gf_i.MakeRef(*this, vsize, vsize);
676 
677  // Copy the updated GridFunctions into the new data array
678  gf_r = *pgfr; gf_r.SyncAliasMemory(*this);
679  gf_i = *pgfi; gf_i.SyncAliasMemory(*this);
680 
681  // Replace the individual data arrays with pointers into the new data
682  // array
683  pgfr->MakeRef(*this, 0, vsize);
684  pgfi->MakeRef(*this, vsize, vsize);
685  }
686  else
687  {
688  // The existing data will not be transferred to the new GridFunctions so
689  // delete it and allocate a new array
690  UseDevice(true);
691  this->SetSize(2 * vsize);
692  this->Vector::operator=(0.0);
693 
694  // Point the individual GridFunctions to the new data array
695  pgfr->MakeRef(*this, 0, vsize);
696  pgfi->MakeRef(*this, vsize, vsize);
697 
698  // These updates will only set the proper 'sequence' value within the
699  // individual GridFunction objects because their sizes are already correct
700  pgfr->Update();
701  pgfi->Update();
702  }
703 }
704 
705 void
707  Coefficient &imag_coeff)
708 {
709  pgfr->SyncMemory(*this);
710  pgfi->SyncMemory(*this);
711  pgfr->ProjectCoefficient(real_coeff);
712  pgfi->ProjectCoefficient(imag_coeff);
713  pgfr->SyncAliasMemory(*this);
714  pgfi->SyncAliasMemory(*this);
715 }
716 
717 void
719  VectorCoefficient &imag_vcoeff)
720 {
721  pgfr->SyncMemory(*this);
722  pgfi->SyncMemory(*this);
723  pgfr->ProjectCoefficient(real_vcoeff);
724  pgfi->ProjectCoefficient(imag_vcoeff);
725  pgfr->SyncAliasMemory(*this);
726  pgfi->SyncAliasMemory(*this);
727 }
728 
729 void
731  Coefficient &imag_coeff,
732  Array<int> &attr)
733 {
734  pgfr->SyncMemory(*this);
735  pgfi->SyncMemory(*this);
736  pgfr->ProjectBdrCoefficient(real_coeff, attr);
737  pgfi->ProjectBdrCoefficient(imag_coeff, attr);
738  pgfr->SyncAliasMemory(*this);
739  pgfi->SyncAliasMemory(*this);
740 }
741 
742 void
744  &real_vcoeff,
746  &imag_vcoeff,
747  Array<int> &attr)
748 {
749  pgfr->SyncMemory(*this);
750  pgfi->SyncMemory(*this);
751  pgfr->ProjectBdrCoefficientNormal(real_vcoeff, attr);
752  pgfi->ProjectBdrCoefficientNormal(imag_vcoeff, attr);
753  pgfr->SyncAliasMemory(*this);
754  pgfi->SyncAliasMemory(*this);
755 }
756 
757 void
759  &real_vcoeff,
761  &imag_vcoeff,
762  Array<int> &attr)
763 {
764  pgfr->SyncMemory(*this);
765  pgfi->SyncMemory(*this);
766  pgfr->ProjectBdrCoefficientTangent(real_vcoeff, attr);
767  pgfi->ProjectBdrCoefficientTangent(imag_vcoeff, attr);
768  pgfr->SyncAliasMemory(*this);
769  pgfi->SyncAliasMemory(*this);
770 }
771 
772 void
774 {
775  ParFiniteElementSpace *pfes = pgfr->ParFESpace();
776  const int tvsize = pfes->GetTrueVSize();
777 
778  tv->Read();
779  Vector tvr; tvr.MakeRef(const_cast<Vector&>(*tv), 0, tvsize);
780  Vector tvi; tvi.MakeRef(const_cast<Vector&>(*tv), tvsize, tvsize);
781 
782  pgfr->SyncMemory(*this);
783  pgfi->SyncMemory(*this);
784  pgfr->Distribute(tvr);
785  pgfi->Distribute(tvi);
786  pgfr->SyncAliasMemory(*this);
787  pgfi->SyncAliasMemory(*this);
788 }
789 
790 void
792 {
793  ParFiniteElementSpace *pfes = pgfr->ParFESpace();
794  const int tvsize = pfes->GetTrueVSize();
795 
796  tv.Write();
797  Vector tvr; tvr.MakeRef(tv, 0, tvsize);
798  Vector tvi; tvi.MakeRef(tv, tvsize, tvsize);
799 
800  pgfr->SyncMemory(*this);
801  pgfi->SyncMemory(*this);
802  pgfr->ParallelProject(tvr);
803  pgfi->ParallelProject(tvi);
804  pgfr->SyncAliasMemory(*this);
805  pgfi->SyncAliasMemory(*this);
806 
807  tvr.SyncAliasMemory(tv);
808  tvi.SyncAliasMemory(tv);
809 }
810 
811 
814  convention)
815  : Vector(2*(pfes->GetVSize())),
816  conv(convention)
817 {
818  UseDevice(true);
819  this->Vector::operator=(0.0);
820 
821  plfr = new ParLinearForm();
822  plfr->MakeRef(pfes, *this, 0);
823 
824  plfi = new ParLinearForm();
825  plfi->MakeRef(pfes, *this, pfes->GetVSize());
826 
827  HYPRE_Int *tdof_offsets_fes = pfes->GetTrueDofOffsets();
828 
829  int n = (HYPRE_AssumedPartitionCheck()) ? 2 : pfes->GetNRanks();
830  tdof_offsets = new HYPRE_Int[n+1];
831 
832  for (int i = 0; i <= n; i++)
833  {
834  tdof_offsets[i] = 2 * tdof_offsets_fes[i];
835  }
836 }
837 
838 
840  ParLinearForm *plf_r,
841  ParLinearForm *plf_i,
843  convention)
844  : Vector(2*(pfes->GetVSize())),
845  conv(convention)
846 {
847  UseDevice(true);
848  this->Vector::operator=(0.0);
849 
850  plfr = new ParLinearForm(pfes, plf_r);
851  plfi = new ParLinearForm(pfes, plf_i);
852 
853  plfr->MakeRef(pfes, *this, 0);
854  plfi->MakeRef(pfes, *this, pfes->GetVSize());
855 
856  HYPRE_Int *tdof_offsets_fes = pfes->GetTrueDofOffsets();
857 
858  int n = (HYPRE_AssumedPartitionCheck()) ? 2 : pfes->GetNRanks();
859  tdof_offsets = new HYPRE_Int[n+1];
860 
861  for (int i = 0; i <= n; i++)
862  {
863  tdof_offsets[i] = 2 * tdof_offsets_fes[i];
864  }
865 }
866 
868 {
869  delete plfr;
870  delete plfi;
871  delete [] tdof_offsets;
872 }
873 
874 void
876  LinearFormIntegrator *lfi_imag)
877 {
878  if ( lfi_real ) { plfr->AddDomainIntegrator(lfi_real); }
879  if ( lfi_imag ) { plfi->AddDomainIntegrator(lfi_imag); }
880 }
881 
882 void
884  LinearFormIntegrator *lfi_imag)
885 {
886  if ( lfi_real ) { plfr->AddBoundaryIntegrator(lfi_real); }
887  if ( lfi_imag ) { plfi->AddBoundaryIntegrator(lfi_imag); }
888 }
889 
890 void
892  LinearFormIntegrator *lfi_imag,
893  Array<int> &bdr_attr_marker)
894 {
895  if ( lfi_real ) { plfr->AddBoundaryIntegrator(lfi_real, bdr_attr_marker); }
896  if ( lfi_imag ) { plfi->AddBoundaryIntegrator(lfi_imag, bdr_attr_marker); }
897 }
898 
899 void
901  LinearFormIntegrator *lfi_imag)
902 {
903  if ( lfi_real ) { plfr->AddBdrFaceIntegrator(lfi_real); }
904  if ( lfi_imag ) { plfi->AddBdrFaceIntegrator(lfi_imag); }
905 }
906 
907 void
909  LinearFormIntegrator *lfi_imag,
910  Array<int> &bdr_attr_marker)
911 {
912  if ( lfi_real ) { plfr->AddBdrFaceIntegrator(lfi_real, bdr_attr_marker); }
913  if ( lfi_imag ) { plfi->AddBdrFaceIntegrator(lfi_imag, bdr_attr_marker); }
914 }
915 
916 void
918 {
919  ParFiniteElementSpace *pfes = (pf != NULL) ? pf : plfr->ParFESpace();
920 
921  UseDevice(true);
922  SetSize(2 * pfes->GetVSize());
923  this->Vector::operator=(0.0);
924 
925  plfr->MakeRef(pfes, *this, 0);
926  plfi->MakeRef(pfes, *this, pfes->GetVSize());
927 }
928 
929 void
931 {
932  plfr->SyncMemory(*this);
933  plfi->SyncMemory(*this);
934  plfr->Assemble();
935  plfi->Assemble();
936  if (conv == ComplexOperator::BLOCK_SYMMETRIC) { *plfi *= -1.0; }
937  plfr->SyncAliasMemory(*this);
938  plfi->SyncAliasMemory(*this);
939 }
940 
941 void
943 {
944  const int tvsize = plfr->ParFESpace()->GetTrueVSize();
945 
946  tv.Write();
947  Vector tvr; tvr.MakeRef(tv, 0, tvsize);
948  Vector tvi; tvi.MakeRef(tv, tvsize, tvsize);
949 
950  plfr->SyncMemory(*this);
951  plfi->SyncMemory(*this);
952  plfr->ParallelAssemble(tvr);
953  plfi->ParallelAssemble(tvi);
954  plfr->SyncAliasMemory(*this);
955  plfi->SyncAliasMemory(*this);
956 
957  tvr.SyncAliasMemory(tv);
958  tvi.SyncAliasMemory(tv);
959 }
960 
963 {
964  const ParFiniteElementSpace *pfes = plfr->ParFESpace();
965  const int tvsize = pfes->GetTrueVSize();
966 
967  HypreParVector *tv = new HypreParVector(pfes->GetComm(),
968  2*(pfes->GlobalTrueVSize()),
969  tdof_offsets);
970 
971  tv->Write();
972  Vector tvr; tvr.MakeRef(*tv, 0, tvsize);
973  Vector tvi; tvi.MakeRef(*tv, tvsize, tvsize);
974 
975  plfr->SyncMemory(*this);
976  plfi->SyncMemory(*this);
977  plfr->ParallelAssemble(tvr);
978  plfi->ParallelAssemble(tvi);
979  plfr->SyncAliasMemory(*this);
980  plfi->SyncAliasMemory(*this);
981 
982  tvr.SyncAliasMemory(*tv);
983  tvi.SyncAliasMemory(*tv);
984 
985  return tv;
986 }
987 
988 complex<double>
990 {
991  plfr->SyncMemory(*this);
992  plfi->SyncMemory(*this);
993  double s = (conv == ComplexOperator::HERMITIAN) ? 1.0 : -1.0;
994  return complex<double>((*plfr)(gf.real()) - s * (*plfi)(gf.imag()),
995  (*plfr)(gf.imag()) + s * (*plfi)(gf.real()));
996 }
997 
998 
999 bool ParSesquilinearForm::RealInteg()
1000 {
1001  int nint = pblfr->GetFBFI()->Size() + pblfr->GetDBFI()->Size() +
1002  pblfr->GetBBFI()->Size() + pblfr->GetBFBFI()->Size();
1003  return (nint != 0);
1004 }
1005 
1006 bool ParSesquilinearForm::ImagInteg()
1007 {
1008  int nint = pblfi->GetFBFI()->Size() + pblfi->GetDBFI()->Size() +
1009  pblfi->GetBBFI()->Size() + pblfi->GetBFBFI()->Size();
1010  return (nint != 0);
1011 }
1012 
1015  convention)
1016  : conv(convention),
1017  pblfr(new ParBilinearForm(pf)),
1018  pblfi(new ParBilinearForm(pf))
1019 {}
1020 
1022  ParBilinearForm *pbfr,
1023  ParBilinearForm *pbfi,
1024  ComplexOperator::Convention convention)
1025  : conv(convention),
1026  pblfr(new ParBilinearForm(pf,pbfr)),
1027  pblfi(new ParBilinearForm(pf,pbfi))
1028 {}
1029 
1031 {
1032  delete pblfr;
1033  delete pblfi;
1034 }
1035 
1037  BilinearFormIntegrator *bfi_imag)
1038 {
1039  if (bfi_real) { pblfr->AddDomainIntegrator(bfi_real); }
1040  if (bfi_imag) { pblfi->AddDomainIntegrator(bfi_imag); }
1041 }
1042 
1043 void
1045  BilinearFormIntegrator *bfi_imag)
1046 {
1047  if (bfi_real) { pblfr->AddBoundaryIntegrator(bfi_real); }
1048  if (bfi_imag) { pblfi->AddBoundaryIntegrator(bfi_imag); }
1049 }
1050 
1051 void
1053  BilinearFormIntegrator *bfi_imag,
1054  Array<int> & bdr_marker)
1055 {
1056  if (bfi_real) { pblfr->AddBoundaryIntegrator(bfi_real, bdr_marker); }
1057  if (bfi_imag) { pblfi->AddBoundaryIntegrator(bfi_imag, bdr_marker); }
1058 }
1059 
1060 void
1062  BilinearFormIntegrator *bfi_imag)
1063 {
1064  if (bfi_real) { pblfr->AddInteriorFaceIntegrator(bfi_real); }
1065  if (bfi_imag) { pblfi->AddInteriorFaceIntegrator(bfi_imag); }
1066 }
1067 
1068 void
1070  BilinearFormIntegrator *bfi_imag)
1071 {
1072  if (bfi_real) { pblfr->AddBdrFaceIntegrator(bfi_real); }
1073  if (bfi_imag) { pblfi->AddBdrFaceIntegrator(bfi_imag); }
1074 }
1075 
1076 void
1078  BilinearFormIntegrator *bfi_imag,
1079  Array<int> &bdr_marker)
1080 {
1081  if (bfi_real) { pblfr->AddBdrFaceIntegrator(bfi_real, bdr_marker); }
1082  if (bfi_imag) { pblfi->AddBdrFaceIntegrator(bfi_imag, bdr_marker); }
1083 }
1084 
1085 void
1087 {
1088  pblfr->Assemble(skip_zeros);
1089  pblfi->Assemble(skip_zeros);
1090 }
1091 
1092 void
1094 {
1095  pblfr->Finalize(skip_zeros);
1096  pblfi->Finalize(skip_zeros);
1097 }
1098 
1101 {
1102  return new ComplexHypreParMatrix(pblfr->ParallelAssemble(),
1103  pblfi->ParallelAssemble(),
1104  true, true, conv);
1105 }
1106 
1107 void
1109  Vector &x, Vector &b,
1110  OperatorHandle &A,
1111  Vector &X, Vector &B,
1112  int ci)
1113 {
1114  ParFiniteElementSpace *pfes = pblfr->ParFESpace();
1115  const int vsize = pfes->GetVSize();
1116 
1117  // Allocate temporary vector
1118  Vector b_0;
1119  b_0.UseDevice(true);
1120  b_0.SetSize(vsize);
1121  b_0 = 0.0;
1122 
1123  // Extract the real and imaginary parts of the input vectors
1124  MFEM_ASSERT(x.Size() == 2 * vsize, "Input GridFunction of incorrect size!");
1125  x.Read();
1126  Vector x_r; x_r.MakeRef(x, 0, vsize);
1127  Vector x_i; x_i.MakeRef(x, vsize, vsize);
1128 
1129  MFEM_ASSERT(b.Size() == 2 * vsize, "Input LinearForm of incorrect size!");
1130  b.Read();
1131  Vector b_r; b_r.MakeRef(b, 0, vsize);
1132  Vector b_i; b_i.MakeRef(b, vsize, vsize);
1133 
1134  if (conv == ComplexOperator::BLOCK_SYMMETRIC) { b_i *= -1.0; }
1135 
1136  const int tvsize = pfes->GetTrueVSize();
1137  OperatorHandle A_r, A_i;
1138 
1139  X.UseDevice(true);
1140  X.SetSize(2 * tvsize);
1141  X = 0.0;
1142 
1143  B.UseDevice(true);
1144  B.SetSize(2 * tvsize);
1145  B = 0.0;
1146 
1147  Vector X_r; X_r.MakeRef(X, 0, tvsize);
1148  Vector X_i; X_i.MakeRef(X, tvsize, tvsize);
1149  Vector B_r; B_r.MakeRef(B, 0, tvsize);
1150  Vector B_i; B_i.MakeRef(B, tvsize, tvsize);
1151 
1152  Vector X_0, B_0;
1153 
1154  if (RealInteg())
1155  {
1156  b_0 = b_r;
1157  pblfr->FormLinearSystem(ess_tdof_list, x_r, b_0, A_r, X_0, B_0, ci);
1158  X_r = X_0; B_r = B_0;
1159 
1160  b_0 = b_i;
1161  pblfr->FormLinearSystem(ess_tdof_list, x_i, b_0, A_r, X_0, B_0, ci);
1162  X_i = X_0; B_i = B_0;
1163 
1164  if (ImagInteg())
1165  {
1166  b_0 = 0.0;
1167  pblfi->FormLinearSystem(ess_tdof_list, x_i, b_0, A_i, X_0, B_0, false);
1168  B_r -= B_0;
1169 
1170  b_0 = 0.0;
1171  pblfi->FormLinearSystem(ess_tdof_list, x_r, b_0, A_i, X_0, B_0, false);
1172  B_i += B_0;
1173  }
1174  }
1175  else if (ImagInteg())
1176  {
1177  b_0 = b_i;
1178  pblfi->FormLinearSystem(ess_tdof_list, x_r, b_0, A_i, X_0, B_0, ci);
1179  X_r = X_0; B_i = B_0;
1180 
1181  b_0 = b_r; b_0 *= -1.0;
1182  pblfi->FormLinearSystem(ess_tdof_list, x_i, b_0, A_i, X_0, B_0, ci);
1183  X_i = X_0; B_r = B_0; B_r *= -1.0;
1184  }
1185  else
1186  {
1187  MFEM_ABORT("Real and Imaginary part of the Sesquilinear form are empty");
1188  }
1189 
1190  if (RealInteg() && ImagInteg())
1191  {
1192  // Modify RHS to conform with standard essential BC treatment
1193  const int n = ess_tdof_list.Size();
1194  auto d_B_r = B_r.Write();
1195  auto d_B_i = B_i.Write();
1196  auto d_X_r = X_r.Read();
1197  auto d_X_i = X_i.Read();
1198  auto d_idx = ess_tdof_list.Read();
1199  MFEM_FORALL(i, n,
1200  {
1201  const int j = d_idx[i];
1202  d_B_r[j] = d_X_r[j];
1203  d_B_i[j] = d_X_i[j];
1204  });
1205  // Modify offdiagonal blocks (imaginary parts of the matrix) to conform
1206  // with standard essential BC treatment
1207  if (A_i.Type() == Operator::Hypre_ParCSR)
1208  {
1209  HypreParMatrix * Ah;
1210  A_i.Get(Ah);
1211  hypre_ParCSRMatrix *Aih = *Ah;
1212  for (int k = 0; k < n; k++)
1213  {
1214  const int j = ess_tdof_list[k];
1215  Aih->diag->data[Aih->diag->i[j]] = 0.0;
1216  }
1217  }
1218  else
1219  {
1220  A_i.As<ConstrainedOperator>()->SetDiagonalPolicy
1221  (mfem::Operator::DiagonalPolicy::DIAG_ZERO);
1222  }
1223  }
1224 
1226  {
1227  B_i *= -1.0;
1228  b_i *= -1.0;
1229  }
1230 
1231  x_r.SyncAliasMemory(x);
1232  x_i.SyncAliasMemory(x);
1233  b_r.SyncAliasMemory(b);
1234  b_i.SyncAliasMemory(b);
1235 
1236  X_r.SyncAliasMemory(X);
1237  X_i.SyncAliasMemory(X);
1238  B_r.SyncAliasMemory(B);
1239  B_i.SyncAliasMemory(B);
1240 
1241  // A = A_r + i A_i
1242  A.Clear();
1243  if ( A_r.Type() == Operator::Hypre_ParCSR ||
1244  A_i.Type() == Operator::Hypre_ParCSR )
1245  {
1246  ComplexHypreParMatrix * A_hyp =
1248  A_i.As<HypreParMatrix>(),
1249  A_r.OwnsOperator(),
1250  A_i.OwnsOperator(),
1251  conv);
1252  A.Reset<ComplexHypreParMatrix>(A_hyp, true);
1253  }
1254  else
1255  {
1256  ComplexOperator * A_op =
1257  new ComplexOperator(A_r.As<Operator>(),
1258  A_i.As<Operator>(),
1259  A_r.OwnsOperator(),
1260  A_i.OwnsOperator(),
1261  conv);
1262  A.Reset<ComplexOperator>(A_op, true);
1263  }
1264  A_r.SetOperatorOwner(false);
1265  A_i.SetOperatorOwner(false);
1266 }
1267 
1268 void
1270  OperatorHandle &A)
1271 {
1272  OperatorHandle A_r, A_i;
1273  if (RealInteg())
1274  {
1275  pblfr->FormSystemMatrix(ess_tdof_list, A_r);
1276  }
1277  if (ImagInteg())
1278  {
1279  pblfi->FormSystemMatrix(ess_tdof_list, A_i);
1280  }
1281  if (!RealInteg() && !ImagInteg())
1282  {
1283  MFEM_ABORT("Both Real and Imaginary part of the Sesquilinear form are empty");
1284  }
1285 
1286  if (RealInteg() && ImagInteg())
1287  {
1288  // Modify offdiagonal blocks (imaginary parts of the matrix) to conform
1289  // with standard essential BC treatment
1290  if ( A_i.Type() == Operator::Hypre_ParCSR )
1291  {
1292  int n = ess_tdof_list.Size();
1293  HypreParMatrix * Ah;
1294  A_i.Get(Ah);
1295  hypre_ParCSRMatrix * Aih = *Ah;
1296  for (int k = 0; k < n; k++)
1297  {
1298  int j = ess_tdof_list[k];
1299  Aih->diag->data[Aih->diag->i[j]] = 0.0;
1300  }
1301  }
1302  else
1303  {
1304  A_i.As<ConstrainedOperator>()->SetDiagonalPolicy
1305  (mfem::Operator::DiagonalPolicy::DIAG_ZERO);
1306  }
1307  }
1308 
1309  // A = A_r + i A_i
1310  A.Clear();
1311  if ( A_r.Type() == Operator::Hypre_ParCSR ||
1312  A_i.Type() == Operator::Hypre_ParCSR )
1313  {
1314  ComplexHypreParMatrix * A_hyp =
1316  A_i.As<HypreParMatrix>(),
1317  A_r.OwnsOperator(),
1318  A_i.OwnsOperator(),
1319  conv);
1320  A.Reset<ComplexHypreParMatrix>(A_hyp, true);
1321  }
1322  else
1323  {
1324  ComplexOperator * A_op =
1325  new ComplexOperator(A_r.As<Operator>(),
1326  A_i.As<Operator>(),
1327  A_r.OwnsOperator(),
1328  A_i.OwnsOperator(),
1329  conv);
1330  A.Reset<ComplexOperator>(A_op, true);
1331  }
1332  A_r.SetOperatorOwner(false);
1333  A_i.SetOperatorOwner(false);
1334 }
1335 
1336 void
1338  Vector &x)
1339 {
1340  ParFiniteElementSpace *pfes = pblfr->ParFESpace();
1341 
1342  const Operator &P = *pfes->GetProlongationMatrix();
1343 
1344  const int vsize = pfes->GetVSize();
1345  const int tvsize = X.Size() / 2;
1346 
1347  X.Read();
1348  Vector X_r; X_r.MakeRef(const_cast<Vector&>(X), 0, tvsize);
1349  Vector X_i; X_i.MakeRef(const_cast<Vector&>(X), tvsize, tvsize);
1350 
1351  x.Write();
1352  Vector x_r; x_r.MakeRef(x, 0, vsize);
1353  Vector x_i; x_i.MakeRef(x, vsize, vsize);
1354 
1355  // Apply conforming prolongation
1356  P.Mult(X_r, x_r);
1357  P.Mult(X_i, x_i);
1358 
1359  x_r.SyncAliasMemory(x);
1360  x_i.SyncAliasMemory(x);
1361 }
1362 
1363 void
1365 {
1366  if ( pblfr ) { pblfr->Update(nfes); }
1367  if ( pblfi ) { pblfi->Update(nfes); }
1368 }
1369 
1370 #endif // MFEM_USE_MPI
1371 
1372 }
int Size() const
Return the logical size of the array.
Definition: array.hpp:124
int GetVSize() const
Return the number of vector dofs, i.e. GetNDofs() x GetVDim().
Definition: fespace.hpp:400
OpType * As() const
Return the Operator pointer statically cast to a specified OpType. Similar to the method Get()...
Definition: handle.hpp:96
Array< BilinearFormIntegrator * > * GetBBFI()
Access all the integrators added with AddBoundaryIntegrator().
virtual void ProjectBdrCoefficientTangent(VectorCoefficient &real_coeff, VectorCoefficient &imag_coeff, Array< int > &attr)
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.
const SparseMatrix * GetConformingProlongation() const
The returned SparseMatrix is owned by the FiniteElementSpace.
Definition: fespace.cpp:874
void FormSystemMatrix(const Array< int > &ess_tdof_list, OperatorHandle &A)
void Assemble(int skip_zeros=1)
Assembles the form i.e. sums over all domain/bdr integrators.
Base class for vector Coefficients that optionally depend on time and space.
virtual void ProjectBdrCoefficient(Coefficient &real_coeff, Coefficient &imag_coeff, Array< int > &attr)
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(...
FiniteElementSpace * FESpace()
Return the FE space associated with the BilinearForm.
virtual void MakeRef(FiniteElementSpace *f, Vector &v, int v_offset)
Make the ParLinearForm reference external data on a new FiniteElementSpace.
Definition: plinearform.cpp:33
void SyncMemory(const Vector &v)
Update the memory location of the vector to match v.
Definition: vector.hpp:190
virtual const Operator * GetProlongationMatrix() const
The returned Operator is owned by the FiniteElementSpace.
Definition: pfespace.cpp:895
Array< BilinearFormIntegrator * > * GetDBFI()
Access all the integrators added with AddDomainIntegrator().
ComplexSparseMatrix * AssembleComplexSparseMatrix()
Returns the matrix assembled on the true dofs, i.e. P^t A P.
void SetSize(int s)
Resize the vector to size s.
Definition: vector.hpp:459
virtual void MakeRef(FiniteElementSpace *f, Vector &v, int v_offset)
Make the LinearForm reference external data on a new FiniteElementSpace.
Definition: linearform.cpp:207
Mimic the action of a complex operator using two real operators.
virtual void ProjectCoefficient(Coefficient &real_coeff, Coefficient &imag_coeff)
Definition: complex_fem.cpp:88
Array< BilinearFormIntegrator * > * GetBFBFI()
Access all integrators added with AddBdrFaceIntegrator().
void Assemble()
Assembles the linear form i.e. sums over all domain/bdr integrators.
Definition: linearform.cpp:79
Pointer to an Operator of a specified type.
Definition: handle.hpp:33
void Finalize(int skip_zeros=1)
Finalizes the matrix initialization.
Operator::Type Type() const
Get the currently set operator type id.
Definition: handle.hpp:91
void UseDevice(bool use_dev) const
Enable execution of Vector operations using the mfem::Device.
Definition: vector.hpp:89
void AddDomainIntegrator(LinearFormIntegrator *lfi_real, LinearFormIntegrator *lfi_imag)
Adds new Domain Integrator.
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(...
int Size() const
Returns the size of the vector.
Definition: vector.hpp:160
void AddDomainIntegrator(LinearFormIntegrator *lfi_real, LinearFormIntegrator *lfi_imag)
Adds new Domain Integrator.
ParGridFunction & imag()
virtual void Mult(const Vector &x, Vector &y) const =0
Operator application: y=A(x).
Abstract parallel finite element space.
Definition: pfespace.hpp:28
virtual void ProjectCoefficient(Coefficient &coeff)
Definition: pgridfunc.cpp:474
HypreParVector * ParallelAssemble()
Returns the vector assembled on the true dofs, i.e. P^t v.
double * Write(bool on_dev=true)
Shortcut for mfem::Write(vec.GetMemory(), vec.Size(), on_dev).
Definition: vector.hpp:380
Array< BilinearFormIntegrator * > * GetFBFI()
Access all integrators added with AddInteriorFaceIntegrator().
Alternate convention for damping operators.
ParComplexLinearForm(ParFiniteElementSpace *pf, ComplexOperator::Convention convention=ComplexOperator::HERMITIAN)
Abstract base class LinearFormIntegrator.
Definition: lininteg.hpp:22
HypreParMatrix * ParallelAssemble()
Returns the matrix assembled on the true dofs, i.e. P^t A P.
void SetDiagonalPolicy(DiagonalPolicy policy)
Sets diagonal policy used upon construction of the linear system.
virtual void Update(FiniteElementSpace *nfes=NULL)
Update the FiniteElementSpace and delete all data associated with the old one.
virtual void ProjectBdrCoefficientNormal(VectorCoefficient &real_coeff, VectorCoefficient &imag_coeff, Array< int > &attr)
FiniteElementSpace * FESpace()
Read+write access to the associated FiniteElementSpace.
Definition: linearform.hpp:106
virtual int GetTrueVSize() const
Return the number of local vector true dofs.
Definition: pfespace.hpp:255
std::complex< double > operator()(const ComplexGridFunction &gf) const
Vector & operator=(const double *v)
Copy Size() entries from v.
Definition: vector.cpp:120
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)
SesquilinearForm(FiniteElementSpace *fes, ComplexOperator::Convention convention=ComplexOperator::HERMITIAN)
virtual void MakeRef(FiniteElementSpace *f, double *v)
Make the ParGridFunction reference external data on a new FiniteElementSpace.
Definition: pgridfunc.cpp:102
ID for class SparseMatrix.
Definition: operator.hpp:240
ParSesquilinearForm(ParFiniteElementSpace *pf, ComplexOperator::Convention convention=ComplexOperator::HERMITIAN)
ComplexHypreParMatrix * ParallelAssemble()
Returns the matrix assembled on the true dofs, i.e. P^t A P.
Class for parallel linear form.
Definition: plinearform.hpp:26
void ProjectBdrCoefficientNormal(VectorCoefficient &vcoeff, Array< int > &bdr_attr)
Definition: gridfunc.cpp:2479
OpType * Is() const
Return the Operator pointer dynamically cast to a specified OpType.
Definition: handle.hpp:100
MPI_Comm GetComm() const
Definition: pfespace.hpp:239
bool UseDevice() const
Return the device flag of the Memory object used by the Vector.
Definition: vector.hpp:92
void AddBdrFaceIntegrator(BilinearFormIntegrator *bfi_real, BilinearFormIntegrator *bfi_imag)
Adds new boundary Face Integrator. Assumes ownership of bfi.
Data type sparse matrix.
Definition: sparsemat.hpp:46
virtual void Update()
Transform by the Space UpdateMatrix (e.g., on Mesh change).
Definition: pgridfunc.cpp:81
ParFiniteElementSpace * ParFESpace() const
Return the parallel FE space associated with the ParBilinearForm.
double b
Definition: lissajous.cpp:42
ParFiniteElementSpace * ParFESpace() const
Definition: plinearform.hpp:76
virtual void ProjectBdrCoefficient(Coefficient &real_coeff, Coefficient &imag_coeff, Array< int > &attr)
HYPRE_Int GlobalTrueVSize() const
Definition: pfespace.hpp:251
void Assemble(int skip_zeros=1)
Assemble the local matrix.
const T * Read(bool on_dev=true) const
Shortcut for mfem::Read(a.GetMemory(), a.Size(), on_dev).
Definition: array.hpp:290
void FormSystemMatrix(const Array< int > &ess_tdof_list, OperatorHandle &A)
virtual void Update(FiniteElementSpace *nfes=NULL)
Update the FiniteElementSpace and delete all data associated with the old one.
void Update(ParFiniteElementSpace *pf=NULL)
virtual void MakeRef(FiniteElementSpace *f, double *v)
Make the GridFunction reference external data on a new FiniteElementSpace.
Definition: gridfunc.cpp:188
virtual void FormSystemMatrix(const Array< int > &ess_tdof_list, OperatorHandle &A)
Form the linear system matrix A, see FormLinearSystem() for details.
void AddBdrFaceIntegrator(LinearFormIntegrator *lfi)
Adds new Boundary Face Integrator. Assumes ownership of lfi.
Definition: linearform.cpp:66
virtual void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x)
virtual int GetTrueVSize() const
Return the number of vector true (conforming) dofs.
Definition: fespace.hpp:403
GridFunction & real()
Definition: complex_fem.hpp:69
void AddBoundaryIntegrator(LinearFormIntegrator *lfi)
Adds new Boundary Integrator. Assumes ownership of lfi.
Definition: linearform.cpp:53
FiniteElementSpace * FESpace()
Definition: gridfunc.hpp:595
Operator * Ptr() const
Access the underlying Operator pointer.
Definition: handle.hpp:82
Memory< T > data
Pointer to data.
Definition: array.hpp:48
Wrapper for hypre&#39;s parallel vector class.
Definition: hypre.hpp:70
void Distribute(const Vector *tv)
void AddDomainIntegrator(LinearFormIntegrator *lfi)
Adds new Domain Integrator. Assumes ownership of lfi.
Definition: linearform.cpp:39
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:31
void Finalize(int skip_zeros=1)
Finalizes the matrix initialization.
Specialization of the ComplexOperator built from a pair of Sparse Matrices.
void ProjectBdrCoefficient(Coefficient *coeff[], VectorCoefficient *vcoeff, Array< int > &attr)
Definition: pgridfunc.cpp:582
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition: fespace.hpp:87
GridFunction & imag()
Definition: complex_fem.hpp:70
virtual void Mult(const Vector &x, Vector &y) const
Matrix vector multiplication.
Definition: sparsemat.cpp:594
bool OwnsOperator() const
Return true if the OperatorHandle owns the held Operator.
Definition: handle.hpp:109
Base class Coefficients that optionally depend on space and time. These are used by the BilinearFormI...
Definition: coefficient.hpp:39
void AddBoundaryIntegrator(BilinearFormIntegrator *bfi_real, BilinearFormIntegrator *bfi_imag)
Adds new Boundary Integrator.
void Clear()
Clear the OperatorHandle, deleting the held Operator (if owned), while leaving the type id unchanged...
Definition: handle.hpp:116
void AddBoundaryIntegrator(LinearFormIntegrator *lfi_real, LinearFormIntegrator *lfi_imag)
Adds new Boundary Integrator.
std::complex< double > operator()(const ParComplexGridFunction &gf) const
void Distribute(const Vector *tv)
Definition: pgridfunc.cpp:132
ParComplexGridFunction(ParFiniteElementSpace *pf)
Construct a ParComplexGridFunction associated with the ParFiniteElementSpace *pf. ...
virtual void RecoverFEMSolution(const Vector &X, const Vector &b, Vector &x)
virtual void Update()
Transform by the Space UpdateMatrix (e.g., on Mesh change).
Definition: gridfunc.cpp:152
void AddDomainIntegrator(BilinearFormIntegrator *bfi_real, BilinearFormIntegrator *bfi_imag)
Adds new Domain Integrator.
virtual void FormSystemMatrix(const Array< int > &ess_tdof_list, OperatorHandle &A)
Form the linear system matrix A, see FormLinearSystem() for details.
A &quot;square matrix&quot; 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 SetOperatorOwner(bool own=true)
Set the ownership flag for the held Operator.
Definition: handle.hpp:112
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)
const double * Read(bool on_dev=true) const
Shortcut for mfem::Read(vec.GetMemory(), vec.Size(), on_dev).
Definition: vector.hpp:372
virtual void Finalize(int skip_zeros=1)
Finalizes the matrix initialization.
void Assemble(int skip_zeros=1)
Assemble the local matrix.
void AddBoundaryIntegrator(BilinearFormIntegrator *bfi)
Adds new Boundary Integrator. Assumes ownership of bfi.
DiagonalPolicy
Defines operator diagonal policy upon elimination of rows and/or columns.
Definition: operator.hpp:45
void AddInteriorFaceIntegrator(BilinearFormIntegrator *bfi)
Adds new interior Face Integrator. Assumes ownership of bfi.
virtual void ProjectBdrCoefficientTangent(VectorCoefficient &vcoeff, Array< int > &bdr_attr)
Project the tangential components of the given VectorCoefficient on the boundary. Only boundary attri...
Definition: pgridfunc.cpp:620
Specialization of the ComplexOperator built from a pair of HypreParMatrices.
void AddDomainIntegrator(BilinearFormIntegrator *bfi)
Adds new Domain Integrator. Assumes ownership of bfi.
Native convention for Hermitian operators.
const Operator * GetUpdateOperator()
Get the GridFunction update operator.
Definition: fespace.hpp:683
void AddDomainIntegrator(BilinearFormIntegrator *bfi_real, BilinearFormIntegrator *bfi_imag)
Adds new Domain Integrator.
virtual void ProjectBdrCoefficientTangent(VectorCoefficient &vcoeff, Array< int > &bdr_attr)
Project the tangential components of the given VectorCoefficient on the boundary. Only boundary attri...
Definition: gridfunc.cpp:2550
virtual void ProjectCoefficient(Coefficient &coeff)
Definition: gridfunc.cpp:2252
Class for parallel bilinear form.
HYPRE_Int * GetTrueDofOffsets() const
Definition: pfespace.hpp:248
void Get(OpType *&A) const
Return the Operator pointer statically cast to a given OpType.
Definition: handle.hpp:106
Vector data type.
Definition: vector.hpp:51
void MakeRef(Vector &base, int offset, int size)
Reset the Vector to be a reference to a sub-vector of base.
Definition: vector.hpp:517
ID for class HypreParMatrix.
Definition: operator.hpp:241
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:23
const SparseMatrix & SpMat() const
Returns a const reference to the sparse matrix.
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 ParallelProject(Vector &tv) const
Returns the vector restricted to the true dofs.
Square Operator for imposing essential boundary conditions using only the action, Mult()...
Definition: operator.hpp:800
Abstract operator.
Definition: operator.hpp:24
Wrapper for hypre&#39;s ParCSR matrix class.
Definition: hypre.hpp:181
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.
virtual void Update(FiniteElementSpace *nfes=NULL)
void ParallelAssemble(Vector &tv)
Assemble the vector on the true dofs, i.e. P^t v.
Definition: plinearform.cpp:46
void AddInteriorFaceIntegrator(BilinearFormIntegrator *bfi_real, BilinearFormIntegrator *bfi_imag)
Adds new interior Face Integrator. Assumes ownership of bfi.
void AddBoundaryIntegrator(LinearFormIntegrator *lfi_real, LinearFormIntegrator *lfi_imag)
Adds new Boundary Integrator.
void Reset(OpType *A, bool own_A=true)
Reset the OperatorHandle to the given OpType pointer, A.
Definition: handle.hpp:137
virtual void ProjectCoefficient(Coefficient &real_coeff, Coefficient &imag_coeff)
void AddBdrFaceIntegrator(BilinearFormIntegrator *bfi)
Adds new boundary Face Integrator. Assumes ownership of bfi.
void ProjectBdrCoefficient(Coefficient &coeff, Array< int > &attr)
Project a Coefficient on the GridFunction, modifying only DOFs on the boundary associated with the bo...
Definition: gridfunc.hpp:399
void ParallelProject(Vector &tv) const
Returns the vector restricted to the true dofs.
Definition: pgridfunc.cpp:171
double f(const Vector &p)
ParFiniteElementSpace * ParFESpace() const
Definition: pgridfunc.hpp:108
void SyncAliasMemory(const Vector &v)
Update the alias memory location of the vector to match v.
Definition: vector.hpp:193