11 #include <dolfinx/common/span.hpp>
12 #include <dolfinx/fem/Function.h>
13 #include <dolfinx/fem/FunctionSpace.h>
14 #include <dolfinx/la/utils.h>
54 const std::array<std::reference_wrapper<const fem::FunctionSpace>, 2>& V,
55 const int dim,
const tcb::span<const std::int32_t>& entities,
79 std::vector<std::int32_t>
81 const tcb::span<const std::int32_t>& entities,
98 const std::array<std::reference_wrapper<const fem::FunctionSpace>, 2>& V,
99 const std::function<Eigen::Array<bool, Eigen::Dynamic, 1>(
100 const Eigen::Ref<
const Eigen::Array<
double, 3, Eigen::Dynamic,
101 Eigen::RowMajor>>&)>& marker_fn);
116 const std::function<Eigen::Array<bool, Eigen::Dynamic, 1>(
117 const Eigen::Ref<
const Eigen::Array<
double, 3, Eigen::Dynamic,
118 Eigen::RowMajor>>&)>& marker_fn);
131 template <
typename T>
150 template <
typename U>
153 _dofs0(std::forward<U>(dofs))
155 const int owned_size0 = _function_space->dofmap()->index_map->size_local();
156 auto it = std::lower_bound(_dofs0.begin(), _dofs0.end(), owned_size0);
157 const int map0_bs = _function_space->dofmap()->index_map_bs();
158 _owned_indices0 = map0_bs * std::distance(_dofs0.begin(), it);
160 const int bs = _function_space->dofmap()->bs();
164 const std::vector<std::int32_t> dof_tmp = _dofs0;
165 _dofs0.resize(bs * dof_tmp.size());
166 for (std::size_t i = 0; i < dof_tmp.size(); ++i)
168 for (
int k = 0; k < bs; ++k)
169 _dofs0[bs * i + k] = bs * dof_tmp[i] + k;
195 const std::array<std::vector<std::int32_t>, 2>& V_g_dofs,
196 std::shared_ptr<const fem::FunctionSpace> V)
197 : _function_space(V), _g(g), _dofs0(V_g_dofs[0]), _dofs1_g(V_g_dofs[1])
199 assert(_dofs0.size() == _dofs1_g.size());
200 assert(_function_space);
203 const int map0_bs = _function_space->dofmap()->index_map_bs();
204 const int map0_size = _function_space->dofmap()->index_map->size_local();
205 const int owned_size0 = map0_bs * map0_size;
206 auto it0 = std::lower_bound(_dofs0.begin(), _dofs0.end(), owned_size0);
207 _owned_indices0 = std::distance(_dofs0.begin(), it0);
232 return _function_space;
237 std::shared_ptr<const fem::Function<T>>
value()
const {
return _g; }
245 std::pair<tcb::span<const std::int32_t>, std::int32_t>
dof_indices()
const
247 return {tcb::make_span(_dofs0), _owned_indices0};
261 void set(tcb::span<T> x,
double scale = 1.0)
const
264 const std::vector<T>& g = _g->x()->array();
265 for (std::size_t i = 0; i < _dofs0.size(); ++i)
267 if (_dofs0[i] < (std::int32_t)x.size())
269 assert(_dofs1_g[i] < (std::int32_t)g.size());
270 x[_dofs0[i]] = scale * g[_dofs1_g[i]];
279 void set(tcb::span<T> x,
const tcb::span<const T>& x0,
280 double scale = 1.0)
const
283 const std::vector<T>& g = _g->x()->array();
284 assert(x.size() <= x0.size());
285 for (std::size_t i = 0; i < _dofs0.size(); ++i)
287 if (_dofs0[i] < (std::int32_t)x.size())
289 assert(_dofs1_g[i] < (std::int32_t)g.size());
290 x[_dofs0[i]] = scale * (g[_dofs1_g[i]] - x0[_dofs0[i]]);
306 const std::vector<T>& g = _g->x()->array();
307 for (std::size_t i = 0; i < _dofs1_g.size(); ++i)
308 values[_dofs0[i]] = g[_dofs1_g[i]];
319 for (std::size_t i = 0; i < _dofs0.size(); ++i)
321 assert(_dofs0[i] < (std::int32_t)markers.size());
322 markers[_dofs0[i]] =
true;
328 std::shared_ptr<const fem::FunctionSpace> _function_space;
331 std::shared_ptr<const fem::Function<T>> _g;
335 std::vector<std::int32_t> _dofs0, _dofs1_g;
338 int _owned_indices0 = -1;
339 int _owned_indices1 = -1;
Interface for setting (strong) Dirichlet boundary conditions.
Definition: DirichletBC.h:133
DirichletBC & operator=(const DirichletBC &bc)=default
Assignment operator.
std::pair< tcb::span< const std::int32_t >, std::int32_t > dof_indices() const
Access dof indices (local indices, unrolled), including ghosts, to which a Dirichlet condition is app...
Definition: DirichletBC.h:245
DirichletBC & operator=(DirichletBC &&bc)=default
Move assignment operator.
DirichletBC(DirichletBC &&bc)=default
Move constructor.
void set(tcb::span< T > x, const tcb::span< const T > &x0, double scale=1.0) const
Set bc entries in x to scale * (x0 - x_bc)
Definition: DirichletBC.h:279
~DirichletBC()=default
Destructor.
void set(tcb::span< T > x, double scale=1.0) const
Set bc entries in x to scale * x_bc
Definition: DirichletBC.h:261
DirichletBC(const std::shared_ptr< const fem::Function< T >> &g, const std::array< std::vector< std::int32_t >, 2 > &V_g_dofs, std::shared_ptr< const fem::FunctionSpace > V)
Create a representation of a Dirichlet boundary condition where the space being constrained and the f...
Definition: DirichletBC.h:194
DirichletBC(const DirichletBC &bc)=default
Copy constructor.
std::shared_ptr< const fem::FunctionSpace > function_space() const
The function space to which boundary conditions are applied.
Definition: DirichletBC.h:230
DirichletBC(const std::shared_ptr< const fem::Function< T >> &g, U &&dofs)
Create a representation of a Dirichlet boundary condition where the space being constrained is the sa...
Definition: DirichletBC.h:151
void dof_values(tcb::span< T > values) const
Definition: DirichletBC.h:303
std::shared_ptr< const fem::Function< T > > value() const
Return boundary value function g.
Definition: DirichletBC.h:237
void mark_dofs(std::vector< bool > &markers) const
Set markers[i] = true if dof i has a boundary condition applied. Value of markers[i] is not changed o...
Definition: DirichletBC.h:317
This class represents a finite element function space defined by a mesh, a finite element,...
Definition: FunctionSpace.h:34
This class represents a function in a finite element function space , given by.
Definition: Function.h:43
A Mesh consists of a set of connected and numbered mesh topological entities, and geometry data.
Definition: Mesh.h:57
std::array< std::vector< std::int32_t >, 2 > locate_dofs_geometrical(const std::array< std::reference_wrapper< const fem::FunctionSpace >, 2 > &V, const std::function< Eigen::Array< bool, Eigen::Dynamic, 1 >(const Eigen::Ref< const Eigen::Array< double, 3, Eigen::Dynamic, Eigen::RowMajor >> &)> &marker_fn)
Finds degrees of freedom whose geometric coordinate is true for the provided marking function.
Definition: DirichletBC.cpp:455
std::array< std::vector< std::int32_t >, 2 > locate_dofs_topological(const std::array< std::reference_wrapper< const fem::FunctionSpace >, 2 > &V, const int dim, const tcb::span< const std::int32_t > &entities, bool remote=true)
Find degrees-of-freedom which belong to the provided mesh entities (topological). Note that degrees-o...
Definition: DirichletBC.cpp:249