9 #include "assemble_matrix_impl.h"
10 #include "assemble_scalar_impl.h"
11 #include "assemble_vector_impl.h"
12 #include <Eigen/Dense>
40 return fem::impl::assemble_scalar(M);
53 fem::impl::assemble_vector(b, L);
76 Eigen::Ref<Eigen::Matrix<T, Eigen::Dynamic, 1>> b,
77 const std::vector<std::shared_ptr<
const Form<T>>>& a,
78 const std::vector<std::vector<std::shared_ptr<
const DirichletBC<T>>>>& bcs1,
79 const std::vector<Eigen::Ref<
const Eigen::Matrix<T, Eigen::Dynamic, 1>>>&
83 fem::impl::apply_lifting(b, a, bcs1, x0, scale);
96 const std::function<
int(std::int32_t,
const std::int32_t*, std::int32_t,
97 const std::int32_t*,
const T*)>& mat_add,
107 std::vector<bool> dof_marker0, dof_marker1;
110 = map0->block_size() * (map0->size_local() + map0->num_ghosts());
113 = map1->block_size() * (map1->size_local() + map1->num_ghosts());
114 for (std::size_t k = 0; k < bcs.size(); ++k)
117 assert(bcs[k]->function_space());
120 dof_marker0.resize(dim0,
false);
121 bcs[k]->mark_dofs(dof_marker0);
125 dof_marker1.resize(dim1,
false);
126 bcs[k]->mark_dofs(dof_marker1);
131 impl::assemble_matrix(mat_add, a, dof_marker0, dof_marker1);
144 template <
typename T>
146 const std::function<
int(std::int32_t,
const std::int32_t*, std::int32_t,
147 const std::int32_t*,
const T*)>& mat_add,
148 const Form<T>& a,
const std::vector<bool>& dof_marker0,
149 const std::vector<bool>& dof_marker1)
152 impl::assemble_matrix(mat_add, a, dof_marker0, dof_marker1);
165 template <
typename T>
167 const std::function<
int(std::int32_t,
const std::int32_t*, std::int32_t,
168 const std::int32_t*,
const T*)>& mat_add,
169 const Eigen::Ref<
const Eigen::Array<std::int32_t, Eigen::Dynamic, 1>>& rows,
172 for (Eigen::Index i = 0; i < rows.size(); ++i)
174 const std::int32_t row = rows(i);
175 mat_add(1, &row, 1, &row, &diagonal);
194 template <
typename T>
196 const std::function<
int(std::int32_t,
const std::int32_t*, std::int32_t,
197 const std::int32_t*,
const T*)>& mat_add,
202 for (
const auto& bc : bcs)
205 if (V.
contains(*bc->function_space()))
206 add_diagonal<T>(mat_add, bc->dofs_owned().col(0), diagonal);
220 template <
typename T>
221 void set_bc(Eigen::Ref<Eigen::Matrix<T, Eigen::Dynamic, 1>> b,
223 const Eigen::Ref<
const Eigen::Matrix<T, Eigen::Dynamic, 1>>& x0,
226 if (b.rows() > x0.rows())
227 throw std::runtime_error(
"Size mismatch between b and x0 vectors.");
228 for (
const auto& bc : bcs)
231 bc->set(b, x0, scale);
238 template <
typename T>
239 void set_bc(Eigen::Ref<Eigen::Matrix<T, Eigen::Dynamic, 1>> b,
243 for (
const auto& bc : bcs)
260 template <
typename T>
261 std::vector<std::vector<std::shared_ptr<const fem::DirichletBC<T>>>>
266 std::vector<std::vector<std::shared_ptr<const fem::DirichletBC<T>>>> bcs0(
268 for (std::size_t i = 0; i < L.size(); ++i)
270 if (L[i]->function_spaces()[0]->contains(*bc->function_space()))
271 bcs0[i].push_back(bc);
285 template <
typename T>
287 std::vector<std::vector<std::shared_ptr<const fem::DirichletBC<T>>>>>
293 std::vector<std::vector<std::shared_ptr<const fem::DirichletBC<T>>>>>
295 for (std::size_t i = 0; i < a.size(); ++i)
297 for (std::size_t j = 0; j < a[i].size(); ++j)
299 bcs1[i].resize(a[j].size());
305 if (a[i][j]->function_spaces()[1]->contains(*bc->function_space()))
306 bcs1[i][j].push_back(bc);
Interface for setting (strong) Dirichlet boundary conditions.
Definition: DirichletBC.h:102
This class represents a finite element function space defined by a mesh, a finite element,...
Definition: FunctionSpace.h:38
bool contains(const FunctionSpace &V) const
Check whether V is subspace of this, or this itself.
Definition: FunctionSpace.cpp:255
void assemble_matrix(const std::function< int(std::int32_t, const std::int32_t *, std::int32_t, const std::int32_t *, const T *)> &mat_add, const Form< T > &a, const std::vector< std::shared_ptr< const DirichletBC< T >>> &bcs)
Assemble bilinear form into a matrix.
Definition: assembler.h:95
void add_diagonal(const std::function< int(std::int32_t, const std::int32_t *, std::int32_t, const std::int32_t *, const T *)> &mat_add, const Eigen::Ref< const Eigen::Array< std::int32_t, Eigen::Dynamic, 1 >> &rows, T diagonal=1.0)
Adds a value to the diagonal of a matrix for specified rows. It is typically called after assembly....
Definition: assembler.h:166
void apply_lifting(Eigen::Ref< Eigen::Matrix< T, Eigen::Dynamic, 1 >> b, const std::vector< std::shared_ptr< const Form< T >>> &a, const std::vector< std::vector< std::shared_ptr< const DirichletBC< T >>>> &bcs1, const std::vector< Eigen::Ref< const Eigen::Matrix< T, Eigen::Dynamic, 1 >>> &x0, double scale)
Modify b such that:
Definition: assembler.h:75
void set_bc(Eigen::Ref< Eigen::Matrix< T, Eigen::Dynamic, 1 >> b, const std::vector< std::shared_ptr< const DirichletBC< T >>> &bcs, const Eigen::Ref< const Eigen::Matrix< T, Eigen::Dynamic, 1 >> &x0, double scale=1.0)
Set bc values in owned (local) part of the PETScVector, multiplied by 'scale'. The vectors b and x0 m...
Definition: assembler.h:221
T assemble_scalar(const Form< T > &M)
Assemble functional into scalar. Caller is responsible for accumulation across processes.
Definition: assembler.h:38
std::vector< std::vector< std::shared_ptr< const fem::DirichletBC< T > > > > bcs_rows(const std::vector< const Form< T > * > &L, const std::vector< std::shared_ptr< const fem::DirichletBC< T >>> &bcs)
Arrange boundary conditions by block.
Definition: assembler.h:262
std::vector< std::vector< std::vector< std::shared_ptr< const fem::DirichletBC< T > > > > > bcs_cols(const std::vector< std::vector< std::shared_ptr< const Form< T >>>> &a, const std::vector< std::shared_ptr< const DirichletBC< T >>> &bcs)
Arrange boundary conditions by block.
Definition: assembler.h:288
void assemble_vector(Eigen::Ref< Eigen::Matrix< T, Eigen::Dynamic, 1 >> b, const Form< T > &L)
Assemble linear form into an Eigen vector.
Definition: assembler.h:50