MFEM v4.7.0
Finite element discretization library
Loading...
Searching...
No Matches
mfem::ParLinearForm Class Reference

Class for parallel linear form. More...

#include <plinearform.hpp>

Inheritance diagram for mfem::ParLinearForm:
[legend]
Collaboration diagram for mfem::ParLinearForm:
[legend]

Public Member Functions

 ParLinearForm ()
 Create an empty ParLinearForm without an associated ParFiniteElementSpace.
 
 ParLinearForm (ParFiniteElementSpace *pf)
 Create a ParLinearForm on the FE space *pf.
 
 ParLinearForm (ParFiniteElementSpace *pf, real_t *data)
 Construct a ParLinearForm using previously allocated array data.
 
 ParLinearForm (ParFiniteElementSpace *pf, ParLinearForm *plf)
 Create a ParLinearForm on the ParFiniteElementSpace *pf, using the same integrators as the ParLinearForm *plf.
 
ParLinearFormoperator= (const ParLinearForm &rhs)
 Copy assignment. Only the data of the base class Vector is copied.
 
ParFiniteElementSpaceParFESpace () const
 
void Update (ParFiniteElementSpace *pf=NULL)
 Update the object according to the given new FE space *pf.
 
void Update (ParFiniteElementSpace *pf, Vector &v, int v_offset)
 Associate a new FE space, *pf, with this object and use the data of v, offset by v_offset, to initialize this object's Vector::data.
 
virtual void MakeRef (FiniteElementSpace *f, Vector &v, int v_offset)
 Make the ParLinearForm reference external data on a new FiniteElementSpace.
 
void MakeRef (ParFiniteElementSpace *pf, Vector &v, int v_offset)
 Make the ParLinearForm reference external data on a new ParFiniteElementSpace.
 
void Assemble ()
 Assembles the ParLinearForm i.e. sums over all domain/bdr integrators.
 
virtual bool SupportsDevice () const
 Return true if assembly on device is supported, false otherwise.
 
void AssembleSharedFaces ()
 
void ParallelAssemble (Vector &tv)
 Assemble the vector on the true dofs, i.e. P^t v.
 
HypreParVectorParallelAssemble ()
 Returns the vector assembled on the true dofs, i.e. P^t v.
 
real_t operator() (const ParGridFunction &gf) const
 Return the action of the ParLinearForm as a linear mapping.
 
ParLinearFormoperator= (real_t value)
 Assign constant values to the ParLinearForm data.
 
ParLinearFormoperator= (const Vector &v)
 Copy the data from a Vector to the ParLinearForm data.
 
- Public Member Functions inherited from mfem::LinearForm
 LinearForm (FiniteElementSpace *f)
 Creates linear form associated with FE space *f.
 
 LinearForm (FiniteElementSpace *f, LinearForm *lf)
 Create a LinearForm on the FiniteElementSpace f, using the same integrators as the LinearForm lf.
 
 LinearForm ()
 Create an empty LinearForm without an associated FiniteElementSpace.
 
 LinearForm (FiniteElementSpace *f, real_t *data)
 Construct a LinearForm using previously allocated array data.
 
LinearFormoperator= (const LinearForm &rhs)
 Copy assignment. Only the data of the base class Vector is copied.
 
MFEM_DEPRECATED FiniteElementSpaceGetFES ()
 (DEPRECATED) Return the FE space associated with the LinearForm.
 
FiniteElementSpaceFESpace ()
 Read+write access to the associated FiniteElementSpace.
 
const FiniteElementSpaceFESpace () const
 Read-only access to the associated FiniteElementSpace.
 
void AddDomainIntegrator (LinearFormIntegrator *lfi)
 Adds new Domain Integrator. Assumes ownership of lfi.
 
void AddDomainIntegrator (LinearFormIntegrator *lfi, Array< int > &elem_marker)
 
void AddBoundaryIntegrator (LinearFormIntegrator *lfi)
 Adds new Boundary Integrator. Assumes ownership of lfi.
 
void AddBoundaryIntegrator (LinearFormIntegrator *lfi, Array< int > &bdr_attr_marker)
 Add new Boundary Integrator, restricted to the given boundary attributes.
 
void AddBdrFaceIntegrator (LinearFormIntegrator *lfi)
 Adds new Boundary Face Integrator. Assumes ownership of lfi.
 
void AddBdrFaceIntegrator (LinearFormIntegrator *lfi, Array< int > &bdr_attr_marker)
 Add new Boundary Face Integrator, restricted to the given boundary attributes.
 
void AddInteriorFaceIntegrator (LinearFormIntegrator *lfi)
 Adds new Interior Face Integrator. Assumes ownership of lfi.
 
Array< LinearFormIntegrator * > * GetDLFI ()
 Access all integrators added with AddDomainIntegrator() which are not DeltaLFIntegrators or they are DeltaLFIntegrators with non-delta coefficients.
 
Array< Array< int > * > * GetDLFI_Marker ()
 Access all domain markers added with AddDomainIntegrator(). If no marker was specified when the integrator was added, the corresponding pointer (to Array<int>) will be NULL.
 
Array< DeltaLFIntegrator * > * GetDLFI_Delta ()
 Access all integrators added with AddDomainIntegrator() which are DeltaLFIntegrators with delta coefficients.
 
Array< LinearFormIntegrator * > * GetBLFI ()
 Access all integrators added with AddBoundaryIntegrator().
 
Array< LinearFormIntegrator * > * GetFLFI ()
 Access all integrators added with AddBdrFaceIntegrator().
 
Array< LinearFormIntegrator * > * GetIFLFI ()
 Access all integrators added with AddInteriorFaceIntegrator().
 
Array< Array< int > * > * GetFLFI_Marker ()
 Access all boundary markers added with AddBdrFaceIntegrator(). If no marker was specified when the integrator was added, the corresponding pointer (to Array<int>) will be NULL.
 
void UseFastAssembly (bool use_fa)
 Which assembly algorithm to use: the new device-compatible fast assembly (true), or the legacy CPU-only algorithm (false).
 
void Assemble ()
 Assembles the linear form i.e. sums over all domain/bdr integrators.
 
void AssembleDelta ()
 Assembles delta functions of the linear form.
 
void Update ()
 Update the object according to the associated FE space fes.
 
void Update (FiniteElementSpace *f)
 Associate a new FE space, *f, with this object and Update() it. *‍/.
 
void Update (FiniteElementSpace *f, Vector &v, int v_offset)
 Associate a new FE space, *f, with this object and use the data of v, offset by v_offset, to initialize this object's Vector::data.
 
real_t operator() (const GridFunction &gf) const
 Return the action of the LinearForm as a linear mapping.
 
LinearFormoperator= (real_t value)
 Redefine '=' for LinearForm = constant.
 
LinearFormoperator= (const Vector &v)
 Copy the data from v.
 
 ~LinearForm ()
 Destroys linear form.
 
- Public Member Functions inherited from mfem::Vector
 Vector ()
 
 Vector (const Vector &)
 Copy constructor. Allocates a new data array and copies the data.
 
 Vector (Vector &&v)
 Move constructor. "Steals" data from its argument.
 
 Vector (int s)
 Creates vector of size s.
 
 Vector (real_t *data_, int size_)
 Creates a vector referencing an array of doubles, owned by someone else.
 
 Vector (Vector &base, int base_offset, int size_)
 Create a Vector referencing a sub-vector of the Vector base starting at the given offset, base_offset, and size size_.
 
 Vector (int size_, MemoryType mt)
 Create a Vector of size size_ using MemoryType mt.
 
 Vector (int size_, MemoryType h_mt, MemoryType d_mt)
 Create a Vector of size size_ using host MemoryType h_mt and device MemoryType d_mt.
 
template<int N, typename T = real_t>
 Vector (const T(&values)[N])
 Create a vector using a braced initializer list.
 
virtual void UseDevice (bool use_dev) const
 Enable execution of Vector operations using the mfem::Device.
 
virtual bool UseDevice () const
 Return the device flag of the Memory object used by the Vector.
 
void Load (std::istream **in, int np, int *dim)
 Reads a vector from multiple files.
 
void Load (std::istream &in, int Size)
 Load a vector from an input stream.
 
void Load (std::istream &in)
 Load a vector from an input stream, reading the size from the stream.
 
void SetSize (int s)
 Resize the vector to size s.
 
void SetSize (int s, MemoryType mt)
 Resize the vector to size s using MemoryType mt.
 
void SetSize (int s, const Vector &v)
 Resize the vector to size s using the MemoryType of v.
 
void SetData (real_t *d)
 
void SetDataAndSize (real_t *d, int s)
 Set the Vector data and size.
 
void NewDataAndSize (real_t *d, int s)
 Set the Vector data and size, deleting the old data, if owned.
 
void NewMemoryAndSize (const Memory< real_t > &mem, int s, bool own_mem)
 Reset the Vector to use the given external Memory mem and size s.
 
void MakeRef (Vector &base, int offset, int size)
 Reset the Vector to be a reference to a sub-vector of base.
 
void MakeRef (Vector &base, int offset)
 Reset the Vector to be a reference to a sub-vector of base without changing its current size.
 
void MakeDataOwner () const
 Set the Vector data (host pointer) ownership flag.
 
void Destroy ()
 Destroy a vector.
 
void DeleteDevice (bool copy_to_host=true)
 Delete the device pointer, if owned. If copy_to_host is true and the data is valid only on device, move it to host before deleting. Invalidates the device memory.
 
int Size () const
 Returns the size of the vector.
 
int Capacity () const
 Return the size of the currently allocated data array.
 
real_tGetData () const
 Return a pointer to the beginning of the Vector data.
 
MFEM_DEPRECATED operator real_t * ()
 Conversion to double *. Deprecated.
 
MFEM_DEPRECATED operator const real_t * () const
 Conversion to const double *. Deprecated.
 
real_tbegin ()
 STL-like begin.
 
real_tend ()
 STL-like end.
 
const real_tbegin () const
 STL-like begin (const version).
 
const real_tend () const
 STL-like end (const version).
 
Memory< real_t > & GetMemory ()
 Return a reference to the Memory object used by the Vector.
 
const Memory< real_t > & GetMemory () const
 Return a reference to the Memory object used by the Vector, const version.
 
void SyncMemory (const Vector &v) const
 Update the memory location of the vector to match v.
 
void SyncAliasMemory (const Vector &v) const
 Update the alias memory location of the vector to match v.
 
bool OwnsData () const
 Read the Vector data (host pointer) ownership flag.
 
void StealData (real_t **p)
 Changes the ownership of the data; after the call the Vector is empty.
 
real_tStealData ()
 Changes the ownership of the data; after the call the Vector is empty.
 
real_tElem (int i)
 Access Vector entries. Index i = 0 .. size-1.
 
const real_tElem (int i) const
 Read only access to Vector entries. Index i = 0 .. size-1.
 
real_toperator() (int i)
 Access Vector entries using () for 0-based indexing.
 
const real_toperator() (int i) const
 Read only access to Vector entries using () for 0-based indexing.
 
real_toperator[] (int i)
 Access Vector entries using [] for 0-based indexing.
 
const real_toperator[] (int i) const
 Read only access to Vector entries using [] for 0-based indexing.
 
real_t operator* (const real_t *) const
 Dot product with a double * array.
 
real_t operator* (const Vector &v) const
 Return the inner-product.
 
Vectoroperator= (const real_t *v)
 Copy Size() entries from v.
 
Vectoroperator= (const Vector &v)
 Copy assignment.
 
Vectoroperator= (Vector &&v)
 Move assignment.
 
Vectoroperator= (real_t value)
 Redefine '=' for vector = constant.
 
Vectoroperator*= (real_t c)
 
Vectoroperator*= (const Vector &v)
 Component-wise scaling: (*this)(i) *= v(i)
 
Vectoroperator/= (real_t c)
 
Vectoroperator/= (const Vector &v)
 Component-wise division: (*this)(i) /= v(i)
 
Vectoroperator-= (real_t c)
 
Vectoroperator-= (const Vector &v)
 
Vectoroperator+= (real_t c)
 
Vectoroperator+= (const Vector &v)
 
VectorAdd (const real_t a, const Vector &Va)
 (*this) += a * Va
 
VectorSet (const real_t a, const Vector &x)
 (*this) = a * x
 
void SetVector (const Vector &v, int offset)
 
void AddSubVector (const Vector &v, int offset)
 
void Neg ()
 (*this) = -(*this)
 
void Reciprocal ()
 (*this)(i) = 1.0 / (*this)(i)
 
void Swap (Vector &other)
 Swap the contents of two Vectors.
 
void cross3D (const Vector &vin, Vector &vout) const
 
void median (const Vector &lo, const Vector &hi)
 v = median(v,lo,hi) entrywise. Implementation assumes lo <= hi.
 
void GetSubVector (const Array< int > &dofs, Vector &elemvect) const
 Extract entries listed in dofs to the output Vector elemvect.
 
void GetSubVector (const Array< int > &dofs, real_t *elem_data) const
 Extract entries listed in dofs to the output array elem_data.
 
void SetSubVector (const Array< int > &dofs, const real_t value)
 Set the entries listed in dofs to the given value.
 
void SetSubVector (const Array< int > &dofs, const Vector &elemvect)
 Set the entries listed in dofs to the values given in the elemvect Vector. Negative dof values cause the -dof-1 position in this Vector to receive the -val from elemvect.
 
void SetSubVector (const Array< int > &dofs, real_t *elem_data)
 Set the entries listed in dofs to the values given the , elem_data array. Negative dof values cause the -dof-1 position in this Vector to receive the -val from elem_data.
 
void AddElementVector (const Array< int > &dofs, const Vector &elemvect)
 Add elements of the elemvect Vector to the entries listed in dofs. Negative dof values cause the -dof-1 position in this Vector to add the -val from elemvect.
 
void AddElementVector (const Array< int > &dofs, real_t *elem_data)
 Add elements of the elem_data array to the entries listed in dofs. Negative dof values cause the -dof-1 position in this Vector to add the -val from elem_data.
 
void AddElementVector (const Array< int > &dofs, const real_t a, const Vector &elemvect)
 Add times the elements of the elemvect Vector to the entries listed in dofs. Negative dof values cause the -dof-1 position in this Vector to add the -a*val from elemvect.
 
void SetSubVectorComplement (const Array< int > &dofs, const real_t val)
 Set all vector entries NOT in the dofs Array to the given val.
 
void Print (std::ostream &out=mfem::out, int width=8) const
 Prints vector to stream out.
 
void Print (adios2stream &out, const std::string &variable_name) const
 
void Print_HYPRE (std::ostream &out) const
 Prints vector to stream out in HYPRE_Vector format.
 
void PrintHash (std::ostream &out) const
 Print the Vector size and hash of its data.
 
void Randomize (int seed=0)
 Set random values in the vector.
 
real_t Norml2 () const
 Returns the l2 norm of the vector.
 
real_t Normlinf () const
 Returns the l_infinity norm of the vector.
 
real_t Norml1 () const
 Returns the l_1 norm of the vector.
 
real_t Normlp (real_t p) const
 Returns the l_p norm of the vector.
 
real_t Max () const
 Returns the maximal element of the vector.
 
real_t Min () const
 Returns the minimal element of the vector.
 
real_t Sum () const
 Return the sum of the vector entries.
 
real_t DistanceSquaredTo (const real_t *p) const
 Compute the square of the Euclidean distance to another vector.
 
real_t DistanceSquaredTo (const Vector &p) const
 Compute the square of the Euclidean distance to another vector.
 
real_t DistanceTo (const real_t *p) const
 Compute the Euclidean distance to another vector.
 
real_t DistanceTo (const Vector &p) const
 Compute the Euclidean distance to another vector.
 
int CheckFinite () const
 Count the number of entries in the Vector for which isfinite is false, i.e. the entry is a NaN or +/-Inf.
 
virtual ~Vector ()
 Destroys vector.
 
virtual const real_tRead (bool on_dev=true) const
 Shortcut for mfem::Read(vec.GetMemory(), vec.Size(), on_dev).
 
virtual const real_tHostRead () const
 Shortcut for mfem::Read(vec.GetMemory(), vec.Size(), false).
 
virtual real_tWrite (bool on_dev=true)
 Shortcut for mfem::Write(vec.GetMemory(), vec.Size(), on_dev).
 
virtual real_tHostWrite ()
 Shortcut for mfem::Write(vec.GetMemory(), vec.Size(), false).
 
virtual real_tReadWrite (bool on_dev=true)
 Shortcut for mfem::ReadWrite(vec.GetMemory(), vec.Size(), on_dev).
 
virtual real_tHostReadWrite ()
 Shortcut for mfem::ReadWrite(vec.GetMemory(), vec.Size(), false).
 

Protected Attributes

ParFiniteElementSpacepfes
 Points to the same object as fes.
 
- Protected Attributes inherited from mfem::LinearForm
FiniteElementSpacefes
 FE space on which the LinearForm lives. Not owned.
 
LinearFormExtensionext = nullptr
 Extension for supporting different assembly levels.
 
bool fast_assembly = false
 Should we use the device-compatible fast assembly algorithm (false by default)
 
int extern_lfs = 0
 Indicates the LinearFormIntegrators stored in domain_integs, domain_delta_integs, boundary_integs, and boundary_face_integs are owned by another LinearForm.
 
Array< LinearFormIntegrator * > domain_integs
 Set of Domain Integrators to be applied.
 
Array< Array< int > * > domain_integs_marker
 
Array< DeltaLFIntegrator * > domain_delta_integs
 Separate array for integrators with delta function coefficients.
 
Array< LinearFormIntegrator * > boundary_integs
 Set of Boundary Integrators to be applied.
 
Array< Array< int > * > boundary_integs_marker
 Entries are not owned.
 
Array< LinearFormIntegrator * > boundary_face_integs
 Set of Boundary Face Integrators to be applied.
 
Array< Array< int > * > boundary_face_integs_marker
 Entries not owned.
 
Array< LinearFormIntegrator * > interior_face_integs
 Set of Internal Face Integrators to be applied.
 
Array< int > domain_delta_integs_elem_id
 The element ids where the centers of the delta functions lie.
 
Array< IntegrationPointdomain_delta_integs_ip
 The reference coordinates where the centers of the delta functions lie.
 
- Protected Attributes inherited from mfem::Vector
Memory< real_tdata
 
int size
 

Additional Inherited Members

- Protected Member Functions inherited from mfem::LinearForm
bool HaveDeltaLocations ()
 If true, the delta locations are not (re)computed during assembly.
 
void ResetDeltaLocations ()
 Force (re)computation of delta locations.
 

Detailed Description

Class for parallel linear form.

Definition at line 26 of file plinearform.hpp.

Constructor & Destructor Documentation

◆ ParLinearForm() [1/4]

mfem::ParLinearForm::ParLinearForm ( )
inline

Create an empty ParLinearForm without an associated ParFiniteElementSpace.

The associated ParFiniteElementSpace can be set later using one of the methods: Update(ParFiniteElementSpace *) or Update(ParFiniteElementSpace *, Vector &, int).

Definition at line 42 of file plinearform.hpp.

◆ ParLinearForm() [2/4]

mfem::ParLinearForm::ParLinearForm ( ParFiniteElementSpace * pf)
inline

Create a ParLinearForm on the FE space *pf.

The pointer pf is not owned by the newly constructed object.

Definition at line 46 of file plinearform.hpp.

◆ ParLinearForm() [3/4]

mfem::ParLinearForm::ParLinearForm ( ParFiniteElementSpace * pf,
real_t * data )
inline

Construct a ParLinearForm using previously allocated array data.

The ParLinearForm does not assume ownership of data which is assumed to be of size at least pf->GetVSize(). Similar to the LinearForm and Vector constructors for externally allocated array, the pointer data can be NULL. The data array can be replaced later using the method SetData().

Definition at line 54 of file plinearform.hpp.

◆ ParLinearForm() [4/4]

mfem::ParLinearForm::ParLinearForm ( ParFiniteElementSpace * pf,
ParLinearForm * plf )
inline

Create a ParLinearForm on the ParFiniteElementSpace *pf, using the same integrators as the ParLinearForm *plf.

The pointer pf is not owned by the newly constructed object.

The integrators in plf are copied as pointers and they are not owned by the newly constructed LinearForm.

Definition at line 64 of file plinearform.hpp.

Member Function Documentation

◆ Assemble()

void mfem::ParLinearForm::Assemble ( )

Assembles the ParLinearForm i.e. sums over all domain/bdr integrators.

When UseFastAssembly(true) has been called and the linear form assembly is compatible with device execution, the assembly will be executed on the device.

Definition at line 46 of file plinearform.cpp.

◆ AssembleSharedFaces()

void mfem::ParLinearForm::AssembleSharedFaces ( )

Definition at line 65 of file plinearform.cpp.

◆ MakeRef() [1/2]

void mfem::ParLinearForm::MakeRef ( FiniteElementSpace * f,
Vector & v,
int v_offset )
virtual

Make the ParLinearForm reference external data on a new FiniteElementSpace.

This method changes the FiniteElementSpace associated with the ParLinearForm to *f and sets the data of the Vector v (plus the v_offset) as external data in the ParLinearForm.

Note
This version of the method will also perform bounds checks when the build option MFEM_DEBUG is enabled.

Reimplemented from mfem::LinearForm.

Definition at line 33 of file plinearform.cpp.

◆ MakeRef() [2/2]

void mfem::ParLinearForm::MakeRef ( ParFiniteElementSpace * pf,
Vector & v,
int v_offset )

Make the ParLinearForm reference external data on a new ParFiniteElementSpace.

This method changes the ParFiniteElementSpace associated with the ParLinearForm to *pf and sets the data of the Vector v (plus the v_offset) as external data in the ParLinearForm.

Note
This version of the method will also perform bounds checks when the build option MFEM_DEBUG is enabled.

Definition at line 40 of file plinearform.cpp.

◆ operator()()

real_t mfem::ParLinearForm::operator() ( const ParGridFunction & gf) const
inline

Return the action of the ParLinearForm as a linear mapping.

Linear forms are linear functionals which map ParGridFunctions to the real numbers. This method performs this mapping which in this case is equivalent as an inner product of the ParLinearForm and ParGridFunction.

Definition at line 138 of file plinearform.hpp.

◆ operator=() [1/3]

ParLinearForm & mfem::ParLinearForm::operator= ( const ParLinearForm & rhs)
inline

Copy assignment. Only the data of the base class Vector is copied.

It is assumed that this object and rhs use ParFiniteElementSpaces that have the same size.

Note
Defining this method overwrites the implicitly defined copy assignment operator.

Definition at line 73 of file plinearform.hpp.

◆ operator=() [2/3]

ParLinearForm & mfem::ParLinearForm::operator= ( const Vector & v)
inline

Copy the data from a Vector to the ParLinearForm data.

Definition at line 148 of file plinearform.hpp.

◆ operator=() [3/3]

ParLinearForm & mfem::ParLinearForm::operator= ( real_t value)
inline

Assign constant values to the ParLinearForm data.

Definition at line 144 of file plinearform.hpp.

◆ ParallelAssemble() [1/2]

HypreParVector * mfem::ParLinearForm::ParallelAssemble ( )

Returns the vector assembled on the true dofs, i.e. P^t v.

Definition at line 101 of file plinearform.cpp.

◆ ParallelAssemble() [2/2]

void mfem::ParLinearForm::ParallelAssemble ( Vector & tv)

Assemble the vector on the true dofs, i.e. P^t v.

Definition at line 95 of file plinearform.cpp.

◆ ParFESpace()

ParFiniteElementSpace * mfem::ParLinearForm::ParFESpace ( ) const
inline

Definition at line 76 of file plinearform.hpp.

◆ SupportsDevice()

bool mfem::ParLinearForm::SupportsDevice ( ) const
virtual

Return true if assembly on device is supported, false otherwise.

Reimplemented from mfem::LinearForm.

Definition at line 57 of file plinearform.cpp.

◆ Update() [1/2]

void mfem::ParLinearForm::Update ( ParFiniteElementSpace * pf,
Vector & v,
int v_offset )

Associate a new FE space, *pf, with this object and use the data of v, offset by v_offset, to initialize this object's Vector::data.

Note
This method does not perform assembly.

Definition at line 27 of file plinearform.cpp.

◆ Update() [2/2]

void mfem::ParLinearForm::Update ( ParFiniteElementSpace * pf = NULL)

Update the object according to the given new FE space *pf.

If the pointer pf is NULL (this is the default value), the FE space already associated with this object is used.

This method should be called when the associated FE space fes has been updated, e.g. after its associated Mesh object has been refined.

Note
This method does not perform assembly.

Definition at line 21 of file plinearform.cpp.

Member Data Documentation

◆ pfes

ParFiniteElementSpace* mfem::ParLinearForm::pfes
protected

Points to the same object as fes.

Definition at line 29 of file plinearform.hpp.


The documentation for this class was generated from the following files: