DOLFIN-X
DOLFIN-X C++ interface
Functions
dolfinx::fem::impl Namespace Reference

Implementation of assembly. More...

Functions

template<typename ScalarType >
void assemble_matrix (const std::function< int(std::int32_t, const std::int32_t *, std::int32_t, const std::int32_t *, const ScalarType *)> &mat_set_values_local, const Form &a, const std::vector< bool > &bc0, const std::vector< bool > &bc1)
 The matrix A must already be initialised. The matrix may be a proxy, i.e. a view into a larger matrix, and assembly is performed using local indices. Rows (bc0) and columns (bc1) with Dirichlet conditions are zeroed. Markers (bc0 and bc1) can be empty if not bcs are applied. Matrix is not finalised.
 
template<typename ScalarType >
void assemble_cells (const std::function< int(std::int32_t, const std::int32_t *, std::int32_t, const std::int32_t *, const ScalarType *)> &mat_set_values_local, const mesh::Mesh &mesh, const std::vector< std::int32_t > &active_cells, const graph::AdjacencyList< std::int32_t > &dofmap0, int num_dofs_per_cell0, const graph::AdjacencyList< std::int32_t > &dofmap1, int num_dofs_per_cell1, const std::vector< bool > &bc0, const std::vector< bool > &bc1, const std::function< void(ScalarType *, const ScalarType *, const ScalarType *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> &kernel, const Eigen::Array< ScalarType, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &coeffs, const Eigen::Array< ScalarType, Eigen::Dynamic, 1 > &constant_values)
 Execute kernel over cells and accumulate result in Mat.
 
template<typename ScalarType >
void assemble_exterior_facets (const std::function< int(std::int32_t, const std::int32_t *, std::int32_t, const std::int32_t *, const ScalarType *)> &mat_set_values_local, const mesh::Mesh &mesh, const std::vector< std::int32_t > &active_facets, const DofMap &dofmap0, const DofMap &dofmap1, const std::vector< bool > &bc0, const std::vector< bool > &bc1, const std::function< void(ScalarType *, const ScalarType *, const ScalarType *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> &fn, const Eigen::Array< ScalarType, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &coeffs, const Eigen::Array< ScalarType, Eigen::Dynamic, 1 > constant_values)
 Execute kernel over exterior facets and accumulate result in Mat.
 
template<typename ScalarType >
void assemble_interior_facets (const std::function< int(std::int32_t, const std::int32_t *, std::int32_t, const std::int32_t *, const ScalarType *)> &mat_set_values_local, const mesh::Mesh &mesh, const std::vector< std::int32_t > &active_facets, const DofMap &dofmap0, const DofMap &dofmap1, const std::vector< bool > &bc0, const std::vector< bool > &bc1, const std::function< void(ScalarType *, const ScalarType *, const ScalarType *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> &kernel, const Eigen::Array< ScalarType, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &coeffs, const std::vector< int > &offsets, const Eigen::Array< ScalarType, Eigen::Dynamic, 1 > &constant_values)
 Execute kernel over interior facets and accumulate result in Mat.
 
PetscScalar assemble_scalar (const fem::Form &M)
 Assemble functional into an scalar.
 
PetscScalar assemble_cells (const mesh::Mesh &mesh, const std::vector< std::int32_t > &active_cells, const std::function< void(PetscScalar *, const PetscScalar *, const PetscScalar *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> &fn, const Eigen::Array< PetscScalar, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &coeffs, const std::vector< PetscScalar > &constant_values)
 Assemble functional over cells.
 
PetscScalar assemble_exterior_facets (const mesh::Mesh &mesh, const std::vector< std::int32_t > &active_cells, const std::function< void(PetscScalar *, const PetscScalar *, const PetscScalar *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> &fn, const Eigen::Array< PetscScalar, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &coeffs, const std::vector< PetscScalar > &constant_values)
 Execute kernel over exterior facets and accumulate result.
 
PetscScalar assemble_interior_facets (const mesh::Mesh &mesh, const std::vector< std::int32_t > &active_cells, const std::function< void(PetscScalar *, const PetscScalar *, const PetscScalar *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> &fn, const Eigen::Array< PetscScalar, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &coeffs, const std::vector< int > &offsets, const std::vector< PetscScalar > &constant_values)
 Assemble functional over interior facets.
 
void assemble_vector (Eigen::Ref< Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> b, const Form &L)
 Assemble linear form into an Eigen vector. More...
 
void assemble_cells (Eigen::Ref< Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> b, const mesh::Mesh &mesh, const std::vector< std::int32_t > &active_cells, const graph::AdjacencyList< std::int32_t > &dofmap, int num_dofs_per_cell, const std::function< void(PetscScalar *, const PetscScalar *, const PetscScalar *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> &kernel, const Eigen::Array< PetscScalar, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &coeffs, const Eigen::Array< PetscScalar, Eigen::Dynamic, 1 > &constant_values)
 Execute kernel over cells and accumulate result in vector.
 
void assemble_exterior_facets (Eigen::Ref< Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> b, const mesh::Mesh &mesh, const std::vector< std::int32_t > &active_facets, const fem::DofMap &dofmap, const std::function< void(PetscScalar *, const PetscScalar *, const PetscScalar *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> &fn, const Eigen::Array< PetscScalar, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &coeffs, const Eigen::Array< PetscScalar, Eigen::Dynamic, 1 > &constant_values)
 Execute kernel over cells and accumulate result in vector.
 
void assemble_interior_facets (Eigen::Ref< Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> b, const mesh::Mesh &mesh, const std::vector< std::int32_t > &active_facets, const fem::DofMap &dofmap, const std::function< void(PetscScalar *, const PetscScalar *, const PetscScalar *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> &fn, const Eigen::Array< PetscScalar, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > &coeffs, const std::vector< int > &offsets, const Eigen::Array< PetscScalar, Eigen::Dynamic, 1 > &constant_values)
 Assemble linear form interior facet integrals into an Eigen vector.
 
void apply_lifting (Eigen::Ref< Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> b, const std::vector< std::shared_ptr< const Form >> a, const std::vector< std::vector< std::shared_ptr< const DirichletBC >>> &bcs1, const std::vector< Eigen::Ref< const Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >>> &x0, double scale)
 Modify b such that: More...
 
void lift_bc (Eigen::Ref< Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> b, const Form &a, const Eigen::Ref< const Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> &bc_values1, const std::vector< bool > &bc_markers1, double scale)
 Modify RHS vector to account for boundary condition. More...
 
void lift_bc (Eigen::Ref< Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> b, const Form &a, const Eigen::Ref< const Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> &bc_values1, const std::vector< bool > &bc_markers1, const Eigen::Ref< const Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> &x0, double scale)
 Modify RHS vector to account for boundary condition such that: b <- / b - scale * A (x_bc - x0) /. More...
 

Detailed Description

Implementation of assembly.

Function Documentation

◆ apply_lifting()

void dolfinx::fem::impl::apply_lifting ( Eigen::Ref< Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >>  b,
const std::vector< std::shared_ptr< const Form >>  a,
const std::vector< std::vector< std::shared_ptr< const DirichletBC >>> &  bcs1,
const std::vector< Eigen::Ref< const Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >>> &  x0,
double  scale 
)

Modify b such that:

b <- b - scale * A_j (g_j - x0_j)

where j is a block (nest) row index. For a non-blocked problem j = 0. The boundary conditions bc1 are on the trial spaces V_j. The forms in [a] must have the same test space as L (from which b was built), but the trial space may differ. If x0 is not supplied, then it is treated as zero.

Parameters
[in,out]bThe vector to be modified
[in]aThe bilinear forms, where a[j] is the form that generates A_j
[in]bcs1List of boundary conditions for each block, i.e. bcs1[2] are the boundary conditions applied to the columns of a[2] / x0[2] block
[in]x0The vectors used in the lifting
[in]scaleScaling to apply

◆ assemble_vector()

void dolfinx::fem::impl::assemble_vector ( Eigen::Ref< Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >>  b,
const Form L 
)

Assemble linear form into an Eigen vector.

Parameters
[in,out]bThe vector to be assembled. It will not be zeroed before assembly.
[in]LThe linear forms to assemble into b

◆ lift_bc() [1/2]

void dolfinx::fem::impl::lift_bc ( Eigen::Ref< Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >>  b,
const Form a,
const Eigen::Ref< const Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> &  bc_values1,
const std::vector< bool > &  bc_markers1,
const Eigen::Ref< const Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> &  x0,
double  scale 
)

Modify RHS vector to account for boundary condition such that: b <- / b - scale * A (x_bc - x0) /.

Parameters
[in,out]bThe vector to be modified
[in]aThe bilinear form that generates A
[in]bc_values1The boundary condition 'values'
[in]bc_markers1The indices (columns of A, rows of x) to which bcs belong
[in]x0The array used in the lifting, typically a 'current solution' in a Newton method
[in]scaleScaling to apply

◆ lift_bc() [2/2]

void dolfinx::fem::impl::lift_bc ( Eigen::Ref< Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >>  b,
const Form a,
const Eigen::Ref< const Eigen::Matrix< PetscScalar, Eigen::Dynamic, 1 >> &  bc_values1,
const std::vector< bool > &  bc_markers1,
double  scale 
)

Modify RHS vector to account for boundary condition.

b <- b - scale * A x_bc /

Parameters
[in,out]bThe vector to be modified
[in]aThe bilinear form that generates A
[in]bc_values1The boundary condition 'values'
[in]bc_markers1The indices (columns of A, rows of x) to which bcs belong
[in]scaleScaling to apply