10 #include <dolfinx/function/Function.h>
11 #include <dolfinx/la/utils.h>
60 Eigen::Array<std::int32_t, Eigen::Dynamic, Eigen::Dynamic>
62 const std::vector<std::reference_wrapper<function::FunctionSpace>>& V,
63 const int dim,
const Eigen::Ref<const Eigen::ArrayXi>& entities,
82 Eigen::Array<std::int32_t, Eigen::Dynamic, Eigen::Dynamic>
84 const std::vector<std::reference_wrapper<function::FunctionSpace>>& V,
85 const std::function<Eigen::Array<bool, Eigen::Dynamic, 1>(
86 const Eigen::Ref<
const Eigen::Array<
double, 3, Eigen::Dynamic,
87 Eigen::RowMajor>>&)>& marker);
100 template <
typename T>
113 const Eigen::Ref<
const Eigen::Array<std::int32_t, Eigen::Dynamic, 1>>&
120 const int owned_size = _function_space->dofmap()->index_map->block_size()
121 * _function_space->dofmap()->index_map->size_local();
122 auto* it = std::lower_bound(_dofs.col(0).data(),
123 _dofs.col(0).data() + _dofs.rows(), owned_size);
124 _owned_indices = std::distance(_dofs.col(0).data(), it);
138 const Eigen::Ref<
const Eigen::Array<std::int32_t, Eigen::Dynamic, 2>>&
140 std::shared_ptr<const function::FunctionSpace> V)
141 : _function_space(V), _g(g), _dofs(V_g_dofs)
143 const int owned_size = _function_space->dofmap()->index_map->block_size()
144 * _function_space->dofmap()->index_map->size_local();
145 auto* it = std::lower_bound(_dofs.col(0).data(),
146 _dofs.col(0).data() + _dofs.rows(), owned_size);
147 _owned_indices = std::distance(_dofs.col(0).data(), it);
172 return _function_space;
177 std::shared_ptr<const function::Function<T>>
value()
const {
return _g; }
181 const Eigen::Array<std::int32_t, Eigen::Dynamic, 2>&
dofs()
const
189 const Eigen::Ref<const Eigen::Array<std::int32_t, Eigen::Dynamic, 2>>
192 return _dofs.block<Eigen::Dynamic, 2>(0, 0, _owned_indices, 2);
197 void set(Eigen::Ref<Eigen::Matrix<T, Eigen::Dynamic, 1>> x,
198 double scale = 1.0)
const
202 auto& g = _g->x()->array();
203 for (Eigen::Index i = 0; i < _dofs.rows(); ++i)
205 if (_dofs(i, 0) < x.rows())
206 x[_dofs(i, 0)] = scale * g[_dofs(i, 1)];
212 void set(Eigen::Ref<Eigen::Matrix<T, Eigen::Dynamic, 1>> x,
213 const Eigen::Ref<
const Eigen::Matrix<T, Eigen::Dynamic, 1>>& x0,
214 double scale = 1.0)
const
218 auto& g = _g->x()->array();
219 assert(x.rows() <= x0.rows());
220 for (Eigen::Index i = 0; i < _dofs.rows(); ++i)
222 if (_dofs(i, 0) < x.rows())
223 x[_dofs(i, 0)] = scale * (g[_dofs(i, 1)] - x0[_dofs(i, 0)]);
230 void dof_values(Eigen::Ref<Eigen::Matrix<T, Eigen::Dynamic, 1>> values)
const
233 auto& g = _g->x()->array();
234 for (Eigen::Index i = 0; i < _dofs.rows(); ++i)
235 values[_dofs(i, 0)] = g[_dofs(i, 1)];
243 for (Eigen::Index i = 0; i < _dofs.rows(); ++i)
245 assert(_dofs(i, 0) < (std::int32_t)markers.size());
246 markers[_dofs(i, 0)] =
true;
252 std::shared_ptr<const function::FunctionSpace> _function_space;
255 std::shared_ptr<const function::Function<T>> _g;
259 Eigen::Array<std::int32_t, Eigen::Dynamic, 2> _dofs;
262 int _owned_indices = -1;