DOLFIN-X
DOLFIN-X C++ interface
Form.h
1 // Copyright (C) 2007-2014 Anders Logg
2 //
3 // This file is part of DOLFINX (https://www.fenicsproject.org)
4 //
5 // SPDX-License-Identifier: LGPL-3.0-or-later
6 
7 #pragma once
8 
9 #include "FormCoefficients.h"
10 #include "FormIntegrals.h"
11 #include <functional>
12 #include <map>
13 #include <memory>
14 #include <set>
15 #include <string>
16 #include <vector>
17 
18 // Forward declaration
19 struct ufc_form;
20 
21 namespace dolfinx
22 {
23 
24 namespace function
25 {
26 template <typename T>
27 class Constant;
28 template <typename T>
29 class Function;
30 class FunctionSpace;
31 } // namespace function
32 
33 namespace mesh
34 {
35 class Mesh;
36 template <typename T>
37 class MeshTags;
38 } // namespace mesh
39 
40 namespace fem
41 {
42 
66 
67 template <typename T>
68 class Form
69 {
70 public:
79  Form(const std::vector<std::shared_ptr<const function::FunctionSpace>>&
83  const std::vector<
84  std::pair<std::string, std::shared_ptr<const function::Constant<T>>>>
85  constants)
86  : _integrals(integrals), _coefficients(coefficients),
87  _constants(constants), _function_spaces(function_spaces)
88  {
89  // Set _mesh from function::FunctionSpace, and check they are the same
90  if (!function_spaces.empty())
91  _mesh = function_spaces[0]->mesh();
92  for (const auto& V : function_spaces)
93  if (_mesh != V->mesh())
94  throw std::runtime_error("Incompatible mesh");
95 
96  // Set markers for default integrals
97  if (_mesh)
98  _integrals.set_default_domains(*_mesh);
99  }
100 
106  explicit Form(
107  const std::vector<std::shared_ptr<const function::FunctionSpace>>&
110  {
111  // Do nothing
112  }
113 
115  Form(Form&& form) = default;
116 
118  virtual ~Form() = default;
119 
123  int rank() const { return _function_spaces.size(); }
124 
129  const std::map<int, std::shared_ptr<const function::Function<T>>>&
130  coefficients)
131  {
132  for (const auto& c : coefficients)
133  _coefficients.set(c.first, c.second);
134  }
135 
140  const std::map<std::string, std::shared_ptr<const function::Function<T>>>&
141  coefficients)
142  {
143  for (const auto& c : coefficients)
144  _coefficients.set(c.first, c.second);
145  }
146 
154  const std::map<std::string, std::shared_ptr<const function::Constant<T>>>&
155  constants)
156  {
157  for (auto const& constant : constants)
158  {
159  // Find matching string in existing constants
160  const std::string name = constant.first;
161  const auto it = std::find_if(
162  _constants.begin(), _constants.end(),
163  [&](const std::pair<
164  std::string, std::shared_ptr<const function::Constant<T>>>& q) {
165  return (q.first == name);
166  });
167  if (it != _constants.end())
168  it->second = constant.second;
169  else
170  throw std::runtime_error("Constant '" + name + "' not found in form");
171  }
172  }
173 
182  void
183  set_constants(const std::vector<std::shared_ptr<const function::Constant<T>>>&
184  constants)
185  {
186  if (constants.size() != _constants.size())
187  throw std::runtime_error("Incorrect number of constants.");
188 
189  // Loop over each constant that user wants to attach
190  for (std::size_t i = 0; i < constants.size(); ++i)
191  {
192  // In this case, the constants don't have names
193  _constants[i] = std::pair("", constants[i]);
194  }
195  }
196 
199  bool all_constants_set() const
200  {
201  for (const auto& constant : _constants)
202  if (!constant.second)
203  return false;
204  return true;
205  }
206 
209  std::set<std::string> get_unset_constants() const
210  {
211  std::set<std::string> unset;
212  for (const auto& constant : _constants)
213  if (!constant.second)
214  unset.insert(constant.first);
215  return unset;
216  }
217 
221  void set_mesh(const std::shared_ptr<const mesh::Mesh>& mesh)
222  {
223  _mesh = mesh;
224  // Set markers for default integrals
225  _integrals.set_default_domains(*_mesh);
226  }
227 
230  std::shared_ptr<const mesh::Mesh> mesh() const { return _mesh; }
231 
235  std::shared_ptr<const function::FunctionSpace> function_space(int i) const
236  {
237  return _function_spaces.at(i);
238  }
239 
242  std::vector<std::shared_ptr<const function::FunctionSpace>>
244  {
245  return _function_spaces;
246  }
247 
250  IntegralType type, int i,
251  const std::function<void(T*, const T*, const T*, const double*,
252  const int*, const std::uint8_t*,
253  const std::uint32_t)>& fn)
254  {
255  _integrals.set_tabulate_tensor(type, i, fn);
256  if (i == -1 and _mesh)
257  _integrals.set_default_domains(*_mesh);
258  }
259 
261  FormCoefficients<T>& coefficients() { return _coefficients; }
262 
264  const FormCoefficients<T>& coefficients() const { return _coefficients; }
265 
267  const FormIntegrals<T>& integrals() const { return _integrals; }
268 
273  const std::vector<
274  std::pair<std::string, std::shared_ptr<const function::Constant<T>>>>&
275  constants() const
276  {
277  return _constants;
278  }
279 
280 private:
281  // Integrals associated with the Form
282  FormIntegrals<T> _integrals;
283 
284  // Coefficients associated with the Form
285  FormCoefficients<T> _coefficients;
286 
287  // Constants associated with the Form
288  std::vector<
289  std::pair<std::string, std::shared_ptr<const function::Constant<T>>>>
290  _constants;
291 
292  // Function spaces (one for each argument)
293  std::vector<std::shared_ptr<const function::FunctionSpace>> _function_spaces;
294 
295  // The mesh (needed for functionals when we don't have any spaces)
296  std::shared_ptr<const mesh::Mesh> _mesh;
297 };
298 } // namespace fem
299 } // namespace dolfinx
dolfinx::fem::Form::all_constants_set
bool all_constants_set() const
Check if all constants associated with the form have been set.
Definition: Form.h:199
dolfinx::fem::IntegralType
IntegralType
Type of integral.
Definition: FormIntegrals.h:26
dolfinx::fem::Form::get_unset_constants
std::set< std::string > get_unset_constants() const
Return names of any constants that have not been set.
Definition: Form.h:209
dolfinx::fem::Form::set_mesh
void set_mesh(const std::shared_ptr< const mesh::Mesh > &mesh)
Set mesh, necessary for functionals when there are no function spaces.
Definition: Form.h:221
dolfinx::fem::Form::Form
Form(const std::vector< std::shared_ptr< const function::FunctionSpace >> &function_spaces)
Create form (no UFC integrals). Integrals can be attached later using FormIntegrals::set_cell_tabulat...
Definition: Form.h:106
dolfinx::fem::Form::coefficients
const FormCoefficients< T > & coefficients() const
Access coefficients.
Definition: Form.h:264
dolfinx::fem::Form::rank
int rank() const
Rank of form (bilinear form = 2, linear form = 1, functional = 0, etc)
Definition: Form.h:123
dolfinx::fem::Form::set_constants
void set_constants(const std::map< std::string, std::shared_ptr< const function::Constant< T >>> &constants)
Set constants based on their names.
Definition: Form.h:153
dolfinx::fem::Form::function_spaces
std::vector< std::shared_ptr< const function::FunctionSpace > > function_spaces() const
Return function spaces for each argument.
Definition: Form.h:243
dolfinx::fem::Form::~Form
virtual ~Form()=default
Destructor.
dolfinx::fem::Form
Class for variational forms.
Definition: assembler.h:29
dolfinx::fem::Form::Form
Form(const std::vector< std::shared_ptr< const function::FunctionSpace >> &function_spaces, const FormIntegrals< T > &integrals, const FormCoefficients< T > &coefficients, const std::vector< std::pair< std::string, std::shared_ptr< const function::Constant< T >>>> constants)
Create form.
Definition: Form.h:79
dolfinx::fem::Form::set_tabulate_tensor
void set_tabulate_tensor(IntegralType type, int i, const std::function< void(T *, const T *, const T *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> &fn)
Register the function for 'tabulate_tensor' for cell integral i.
Definition: Form.h:249
dolfinx::mesh::MeshTags
A MeshTags are used to associate mesh entities with values. The entity index (local to process) ident...
Definition: Form.h:37
dolfinx::function::Constant
A constant value which can be attached to a Form. Constants may be scalar (rank 0),...
Definition: Form.h:27
dolfinx::fem::Form::coefficients
FormCoefficients< T > & coefficients()
Access coefficients.
Definition: Form.h:261
dolfinx::fem::Form::set_coefficients
void set_coefficients(const std::map< std::string, std::shared_ptr< const function::Function< T >>> &coefficients)
Set coefficient with given name (shared pointer version)
Definition: Form.h:139
dolfinx::fem::Form::constants
const std::vector< std::pair< std::string, std::shared_ptr< const function::Constant< T > > > > & constants() const
Access constants.
Definition: Form.h:275
dolfinx::fem::Form::integrals
const FormIntegrals< T > & integrals() const
Access form integrals.
Definition: Form.h:267
dolfinx::fem::Form::set_constants
void set_constants(const std::vector< std::shared_ptr< const function::Constant< T >>> &constants)
Set constants based on their order (without names)
Definition: Form.h:183
dolfinx::fem::Form::function_space
std::shared_ptr< const function::FunctionSpace > function_space(int i) const
Return function space for given argument.
Definition: Form.h:235
dolfinx::fem::Form::mesh
std::shared_ptr< const mesh::Mesh > mesh() const
Extract common mesh from form.
Definition: Form.h:230
dolfinx::function::Function
This class represents a function in a finite element function space , given by.
Definition: DirichletBC.h:22
dolfinx::fem::FormIntegrals
Integrals of a Form, including those defined over cells, interior and exterior facets,...
Definition: FormIntegrals.h:38
dolfinx::fem::Form::set_coefficients
void set_coefficients(const std::map< int, std::shared_ptr< const function::Function< T >>> &coefficients)
Set coefficient with given number (shared pointer version)
Definition: Form.h:128
dolfinx::function::FunctionSpace
This class represents a finite element function space defined by a mesh, a finite element,...
Definition: FunctionSpace.h:37
dolfinx::fem::FormCoefficients
Storage for the coefficients of a Form consisting of Function and the Element objects they are define...
Definition: FormCoefficients.h:33