MFEM v4.7.0
Finite element discretization library
Loading...
Searching...
No Matches
pcomplexweakform.hpp
Go to the documentation of this file.
1// Copyright (c) 2010-2024, 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_PCOMPLEX_DPGWEAKFORM
13#define MFEM_PCOMPLEX_DPGWEAKFORM
14
15#include "mfem.hpp"
16#include "complexweakform.hpp"
17
18#ifdef MFEM_USE_MPI
19
20namespace mfem
21{
22
23/** @brief Class representing the parallel weak formulation.
24 (Convenient for DPG Equations) */
26{
27
28protected:
29 /// Trial FE spaces
31
32 /// ess_tdof list for each space
34
35 /** split ess_tdof_list given in global tdof (for all spaces)
36 to individual lists for each space */
37 void FillEssTdofLists(const Array<int> & ess_tdof_list);
38
39 // Block Prolongation
40 BlockOperator * P = nullptr;
41 // Block Restriction
42 BlockMatrix * R = nullptr;
43
49
50 void BuildProlongation();
51
52public:
53
55
56 /// Creates bilinear form associated with FE spaces @a trial_pfes_.
63
66 {
67 trial_pfes = trial_pfes_;
69
71 for (int i = 0; i<trial_sfes.Size(); i++)
72 {
73 trial_sfes[i] = (FiniteElementSpace *)trial_pfes[i];
74 ess_tdofs[i] = new Array<int>();
75 }
76 SetSpaces(trial_sfes,fecol_);
77 }
78
79
80 /// Assembles the form i.e. sums over all domain integrators.
81 void Assemble(int skip_zeros = 1);
82
83 /// Returns the matrix assembled on the true dofs, i.e. P^t A P.
84 /** The returned matrix has to be deleted by the caller. */
86
87 void FormLinearSystem(const Array<int> &ess_tdof_list, Vector &x,
89 Vector &X, Vector &B,
90 int copy_interior = 0);
91
92 void FormSystemMatrix(const Array<int> &ess_tdof_list, OperatorHandle &A);
93
94 /** Call this method after solving a linear system constructed using the
95 FormLinearSystem method to recover the solution as a ParGridFunction-size
96 vector in x. Use the same arguments as in the FormLinearSystem call. */
97 virtual void RecoverFEMSolution(const Vector &X, Vector &x);
98
99 virtual void Update();
100
101 /// Destroys bilinear form.
102 virtual ~ParComplexDPGWeakForm();
103
104
105};
106
107} // namespace mfem
108
109
110#endif // MFEM_USE_MPI
111
112#endif
void SetSize(int nsize)
Change the logical size of the array, keep existing entries.
Definition array.hpp:697
int Size() const
Return the logical size of the array.
Definition array.hpp:144
A class to handle Block systems in a matrix-free implementation.
Class representing the DPG weak formulation for complex valued systems (see the class DPGWeakForm).
void SetSpaces(Array< FiniteElementSpace * > &fes_, Array< FiniteElementCollection * > &fecol_)
BlockMatrix * mat_r
Block matrix to be associated with the real/imag Block bilinear form. Owned.
Mimic the action of a complex operator using two real operators.
Class FiniteElementSpace - responsible for providing FEM view of the mesh, mainly managing the set of...
Definition fespace.hpp:220
Pointer to an Operator of a specified type.
Definition handle.hpp:34
Class representing the parallel weak formulation. (Convenient for DPG Equations)
void FormSystemMatrix(const Array< int > &ess_tdof_list, OperatorHandle &A)
Array< Array< int > * > ess_tdofs
ess_tdof list for each space
virtual ~ParComplexDPGWeakForm()
Destroys bilinear form.
Array< ParFiniteElementSpace * > trial_pfes
Trial FE spaces.
void FillEssTdofLists(const Array< int > &ess_tdof_list)
ParComplexDPGWeakForm(Array< ParFiniteElementSpace * > &trial_pfes_, Array< FiniteElementCollection * > &fecol_)
Creates bilinear form associated with FE spaces trial_pfes_.
void FormLinearSystem(const Array< int > &ess_tdof_list, Vector &x, OperatorHandle &A, Vector &X, Vector &B, int copy_interior=0)
void Assemble(int skip_zeros=1)
Assembles the form i.e. sums over all domain integrators.
virtual void RecoverFEMSolution(const Vector &X, Vector &x)
void SetParSpaces(Array< ParFiniteElementSpace * > &trial_pfes_, Array< FiniteElementCollection * > &fecol_)
void ParallelAssemble(BlockMatrix *mat_r, BlockMatrix *mat_i)
Returns the matrix assembled on the true dofs, i.e. P^t A P.
Vector data type.
Definition vector.hpp:80