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

Efficient batched assembly of LOR discretizations on device. More...

#include <lor_batched.hpp>

Collaboration diagram for mfem::BatchedLORAssembly:
[legend]

Public Member Functions

 BatchedLORAssembly (FiniteElementSpace &fes_ho_)
 Construct the batched assembly object corresponding to fes_ho_.
 
void Assemble (BilinearForm &a, const Array< int > ess_dofs, OperatorHandle &A)
 Assemble the given form as a matrix and place the result in A.
 
const VectorGetLORVertexCoordinates ()
 Return the vertices of the LOR mesh in E-vector format.
 
GPU kernel functions

These functions should be considered protected, but they contain mfem::forall kernels, and so they must be public (this is a compiler limitation).

int FillI (SparseMatrix &A) const
 Fill the I array of the sparse matrix A.
 
void FillJAndData (SparseMatrix &A) const
 Fill the J and data arrays of the sparse matrix A.
 
void ParAssemble (BilinearForm &a, const Array< int > &ess_dofs, OperatorHandle &A)
 Assemble the system in parallel and place the result in A.
 

Static Public Member Functions

static bool FormIsSupported (BilinearForm &a)
 Returns true if the form a supports batched assembly, false otherwise.
 
static void FormLORVertexCoordinates (FiniteElementSpace &fes_ho, Vector &X_vert)
 Compute the vertices of the LOR mesh and place the result in X_vert.
 

Protected Member Functions

void SparseIJToCSR (OperatorHandle &A) const
 After assembling the "sparse IJ" format, convert it to CSR.
 
void AssembleWithoutBC (BilinearForm &a, OperatorHandle &A)
 Assemble the system without eliminating essential DOFs.
 
template<typename LOR_KERNEL >
void AssemblyKernel (BilinearForm &a)
 Fill in sparse_ij and sparse_mapping using one of the specialized LOR assembly kernel classes.
 

Protected Attributes

FiniteElementSpacefes_ho
 The high-order space.
 
Vector X_vert
 LOR vertex coordinates.
 
Vector sparse_ij
 The elementwise LOR matrices in a sparse "ij" format.
 
Array< int > sparse_mapping
 The sparsity pattern of the element matrices.
 

Detailed Description

Efficient batched assembly of LOR discretizations on device.

This class should typically be used by the user-facing classes LORDiscretization and ParLORDiscretization. Only certain bilinear forms are supported, currently:

  • H1 diffusion + mass
  • ND curl-curl + mass
  • RT div-div + mass

Whether a form is supported can be checked with the static member function BatchedLORAssembly::FormIsSupported.

Definition at line 33 of file lor_batched.hpp.

Constructor & Destructor Documentation

◆ BatchedLORAssembly()

mfem::BatchedLORAssembly::BatchedLORAssembly ( FiniteElementSpace & fes_ho_)

Construct the batched assembly object corresponding to fes_ho_.

Definition at line 494 of file lor_batched.cpp.

Member Function Documentation

◆ Assemble()

void mfem::BatchedLORAssembly::Assemble ( BilinearForm & a,
const Array< int > ess_dofs,
OperatorHandle & A )

Assemble the given form as a matrix and place the result in A.

In serial, the result will be a SparseMatrix. In parallel, the result will be a HypreParMatrix.

Definition at line 477 of file lor_batched.cpp.

◆ AssembleWithoutBC()

void mfem::BatchedLORAssembly::AssembleWithoutBC ( BilinearForm & a,
OperatorHandle & A )
protected

Assemble the system without eliminating essential DOFs.

Definition at line 429 of file lor_batched.cpp.

◆ AssemblyKernel()

template<typename LOR_KERNEL >
void mfem::BatchedLORAssembly::AssemblyKernel ( BilinearForm & a)
protected

Fill in sparse_ij and sparse_mapping using one of the specialized LOR assembly kernel classes.

See also
Specialization classes: BatchedLOR_H1, BatchedLOR_ND, BatchedLOR_RT

Definition at line 418 of file lor_batched.cpp.

◆ FillI()

int mfem::BatchedLORAssembly::FillI ( SparseMatrix & A) const

Fill the I array of the sparse matrix A.

Note
AssemblyKernel must be called first to populate sparse_mapping.

Definition at line 137 of file lor_batched.cpp.

◆ FillJAndData()

void mfem::BatchedLORAssembly::FillJAndData ( SparseMatrix & A) const

Fill the J and data arrays of the sparse matrix A.

Note
AssemblyKernel must be called first to populate sparse_mapping and sparse_ij.

Definition at line 230 of file lor_batched.cpp.

◆ FormIsSupported()

bool mfem::BatchedLORAssembly::FormIsSupported ( BilinearForm & a)
static

Returns true if the form a supports batched assembly, false otherwise.

Definition at line 49 of file lor_batched.cpp.

◆ FormLORVertexCoordinates()

void mfem::BatchedLORAssembly::FormLORVertexCoordinates ( FiniteElementSpace & fes_ho,
Vector & X_vert )
static

Compute the vertices of the LOR mesh and place the result in X_vert.

Definition at line 72 of file lor_batched.cpp.

◆ GetLORVertexCoordinates()

const Vector & mfem::BatchedLORAssembly::GetLORVertexCoordinates ( )
inline

Return the vertices of the LOR mesh in E-vector format.

Definition at line 74 of file lor_batched.hpp.

◆ ParAssemble()

void mfem::BatchedLORAssembly::ParAssemble ( BilinearForm & a,
const Array< int > & ess_dofs,
OperatorHandle & A )

Assemble the system in parallel and place the result in A.

Definition at line 460 of file lor_batched.cpp.

◆ SparseIJToCSR()

void mfem::BatchedLORAssembly::SparseIJToCSR ( OperatorHandle & A) const
protected

After assembling the "sparse IJ" format, convert it to CSR.

Definition at line 361 of file lor_batched.cpp.

Member Data Documentation

◆ fes_ho

FiniteElementSpace& mfem::BatchedLORAssembly::fes_ho
protected

The high-order space.

Definition at line 36 of file lor_batched.hpp.

◆ sparse_ij

Vector mfem::BatchedLORAssembly::sparse_ij
protected

The elementwise LOR matrices in a sparse "ij" format.

This is interpreted to have shape (nnz_per_row, ndof_per_el, nel_ho). For index (i, j, k), this represents row j of the kth element matrix. The column index is given by sparse_mapping(i, j).

Definition at line 45 of file lor_batched.hpp.

◆ sparse_mapping

Array<int> mfem::BatchedLORAssembly::sparse_mapping
protected

The sparsity pattern of the element matrices.

This array should be interpreted as having shape (nnz_per_row, ndof_per_el). For local DOF index j, sparse_mapping(i, j) is the column index of the ith nonzero in the jth row. If the index is negative, that entry should be skipped (there is no corresponding nonzero).

Definition at line 54 of file lor_batched.hpp.

◆ X_vert

Vector mfem::BatchedLORAssembly::X_vert
protected

LOR vertex coordinates.

Definition at line 38 of file lor_batched.hpp.


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