10 #include <dolfinx/common/span.hpp>
11 #include <dolfinx/fem/Function.h>
12 #include <dolfinx/fem/FunctionSpace.h>
13 #include <dolfinx/la/utils.h>
53 const std::array<std::reference_wrapper<const fem::FunctionSpace>, 2>& V,
54 const int dim,
const tcb::span<const std::int32_t>& entities,
78 std::vector<std::int32_t>
80 const tcb::span<const std::int32_t>& entities,
97 const std::array<std::reference_wrapper<const fem::FunctionSpace>, 2>& V,
98 const std::function<std::vector<bool>(
const array2d<double>&)>& marker_fn);
113 const std::function<std::vector<bool>(
const array2d<double>&)>& marker_fn);
126 template <
typename T>
145 template <
typename U>
148 _dofs0(std::forward<U>(dofs))
150 const int owned_size0 = _function_space->dofmap()->index_map->size_local();
151 auto it = std::lower_bound(_dofs0.begin(), _dofs0.end(), owned_size0);
152 const int map0_bs = _function_space->dofmap()->index_map_bs();
153 _owned_indices0 = map0_bs * std::distance(_dofs0.begin(), it);
155 const int bs = _function_space->dofmap()->bs();
159 const std::vector<std::int32_t> dof_tmp = _dofs0;
160 _dofs0.resize(bs * dof_tmp.size());
161 for (std::size_t i = 0; i < dof_tmp.size(); ++i)
163 for (
int k = 0; k < bs; ++k)
164 _dofs0[bs * i + k] = bs * dof_tmp[i] + k;
190 const std::array<std::vector<std::int32_t>, 2>& V_g_dofs,
191 std::shared_ptr<const fem::FunctionSpace> V)
192 : _function_space(V), _g(g), _dofs0(V_g_dofs[0]), _dofs1_g(V_g_dofs[1])
194 assert(_dofs0.size() == _dofs1_g.size());
195 assert(_function_space);
198 const int map0_bs = _function_space->dofmap()->index_map_bs();
199 const int map0_size = _function_space->dofmap()->index_map->size_local();
200 const int owned_size0 = map0_bs * map0_size;
201 auto it0 = std::lower_bound(_dofs0.begin(), _dofs0.end(), owned_size0);
202 _owned_indices0 = std::distance(_dofs0.begin(), it0);
227 return _function_space;
232 std::shared_ptr<const fem::Function<T>>
value()
const {
return _g; }
240 std::pair<tcb::span<const std::int32_t>, std::int32_t>
dof_indices()
const
242 return {tcb::make_span(_dofs0), _owned_indices0};
256 void set(tcb::span<T> x,
double scale = 1.0)
const
259 const std::vector<T>& g = _g->x()->array();
260 for (std::size_t i = 0; i < _dofs0.size(); ++i)
262 if (_dofs0[i] < (std::int32_t)x.size())
264 assert(_dofs1_g[i] < (std::int32_t)g.size());
265 x[_dofs0[i]] = scale * g[_dofs1_g[i]];
274 void set(tcb::span<T> x,
const tcb::span<const T>& x0,
275 double scale = 1.0)
const
278 const std::vector<T>& g = _g->x()->array();
279 assert(x.size() <= x0.size());
280 for (std::size_t i = 0; i < _dofs0.size(); ++i)
282 if (_dofs0[i] < (std::int32_t)x.size())
284 assert(_dofs1_g[i] < (std::int32_t)g.size());
285 x[_dofs0[i]] = scale * (g[_dofs1_g[i]] - x0[_dofs0[i]]);
301 const std::vector<T>& g = _g->x()->array();
302 for (std::size_t i = 0; i < _dofs1_g.size(); ++i)
303 values[_dofs0[i]] = g[_dofs1_g[i]];
314 for (std::size_t i = 0; i < _dofs0.size(); ++i)
316 assert(_dofs0[i] < (std::int32_t)markers.size());
317 markers[_dofs0[i]] =
true;
323 std::shared_ptr<const fem::FunctionSpace> _function_space;
326 std::shared_ptr<const fem::Function<T>> _g;
330 std::vector<std::int32_t> _dofs0, _dofs1_g;
333 int _owned_indices0 = -1;
334 int _owned_indices1 = -1;
This class provides a dynamic 2-dimensional row-wise array data structure.
Definition: array2d.h:21
Interface for setting (strong) Dirichlet boundary conditions.
Definition: DirichletBC.h:128
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:240
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:274
~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:256
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:189
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:225
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:146
void dof_values(tcb::span< T > values) const
Definition: DirichletBC.h:298
std::shared_ptr< const fem::Function< T > > value() const
Return boundary value function g.
Definition: DirichletBC.h:232
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:312
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:44
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< std::vector< bool >(const array2d< double > &)> &marker_fn)
Finds degrees of freedom whose geometric coordinate is true for the provided marking function.
Definition: DirichletBC.cpp:450
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:244