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 <dolfinx/common/MPI.h>
12 #include <dolfinx/fem/DofMap.h>
13 #include <functional>
14 #include <map>
15 #include <memory>
16 #include <set>
17 #include <string>
18 #include <vector>
19 
20 namespace dolfinx
21 {
22 
23 namespace function
24 {
25 template <typename T>
26 class Constant;
27 template <typename T>
28 class Function;
29 class FunctionSpace;
30 } // namespace function
31 
32 namespace mesh
33 {
34 class Mesh;
35 template <typename T>
36 class MeshTags;
37 } // namespace mesh
38 
39 namespace fem
40 {
41 
65 
66 template <typename T>
67 class Form
68 {
69 public:
78  Form(const std::vector<std::shared_ptr<const function::FunctionSpace>>&
82  const std::vector<
83  std::pair<std::string, std::shared_ptr<const function::Constant<T>>>>
84  constants)
85  : _integrals(integrals), _coefficients(coefficients),
86  _constants(constants), _function_spaces(function_spaces)
87  {
88  // Set _mesh from function::FunctionSpace, and check they are the same
89  if (!function_spaces.empty())
90  _mesh = function_spaces[0]->mesh();
91  for (const auto& V : function_spaces)
92  if (_mesh != V->mesh())
93  throw std::runtime_error("Incompatible mesh");
94 
95  // Set markers for default integrals
96  if (_mesh)
97  _integrals.set_default_domains(*_mesh);
98  }
99 
108  Form(const std::vector<std::shared_ptr<const function::FunctionSpace>>&
110  bool need_mesh_permutation_data)
111  : Form(function_spaces, FormIntegrals<T>({}, need_mesh_permutation_data),
112  FormCoefficients<T>({}), {})
113  {
114  // Do nothing
115  }
116 
118  Form(const Form& form) = delete;
119 
121  Form(Form&& form) = default;
122 
124  virtual ~Form() = default;
125 
129  int rank() const { return _function_spaces.size(); }
130 
135  const std::map<std::string, std::shared_ptr<const function::Function<T>>>&
136  coefficients)
137  {
138  std::for_each(coefficients.begin(), coefficients.end(),
139  [this](auto& c) { _coefficients.set(c.first, c.second); });
140  }
141 
145  const std::map<std::string, std::shared_ptr<const function::Constant<T>>>&
146  constants)
147  {
148  for (auto const& constant : constants)
149  {
150  // Find matching string in existing constants
151  const std::string& name = constant.first;
152  auto it
153  = std::find_if(_constants.begin(), _constants.end(),
154  [&name](const auto& q) { return q.first == name; });
155  if (it != _constants.end())
156  it->second = constant.second;
157  else
158  throw std::runtime_error("Constant '" + name + "' not found in form");
159  }
160  }
161 
164  bool all_constants_set() const
165  {
166  for (const auto& constant : _constants)
167  if (!constant.second)
168  return false;
169  return true;
170  }
171 
174  std::set<std::string> get_unset_constants() const
175  {
176  std::set<std::string> unset;
177  std::for_each(_constants.begin(), _constants.end(), [&unset](auto& c) {
178  if (!c.second)
179  unset.insert(c.first);
180  });
181  return unset;
182  }
183 
190  void set_mesh(const std::shared_ptr<const mesh::Mesh>& mesh)
191  {
192  _mesh = mesh;
193  // Set markers for default integrals
194  _integrals.set_default_domains(*_mesh);
195  }
196 
199  std::shared_ptr<const mesh::Mesh> mesh() const { return _mesh; }
200 
204  std::shared_ptr<const function::FunctionSpace> function_space(int i) const
205  {
206  return _function_spaces.at(i);
207  }
208 
211  std::vector<std::shared_ptr<const function::FunctionSpace>>
213  {
214  return _function_spaces;
215  }
216 
219  IntegralType type, int i,
220  const std::function<void(T*, const T*, const T*, const double*,
221  const int*, const std::uint8_t*,
222  const std::uint32_t)>& fn)
223  {
224  _integrals.set_tabulate_tensor(type, i, fn);
225  if (i == -1 and _mesh)
226  _integrals.set_default_domains(*_mesh);
227  }
228 
230  FormCoefficients<T>& coefficients() { return _coefficients; }
231 
233  const FormCoefficients<T>& coefficients() const { return _coefficients; }
234 
236  const FormIntegrals<T>& integrals() const { return _integrals; }
237 
242  std::vector<
243  std::pair<std::string, std::shared_ptr<const function::Constant<T>>>>&
245  {
246  return _constants;
247  }
248 
253  const std::vector<
254  std::pair<std::string, std::shared_ptr<const function::Constant<T>>>>&
255  constants() const
256  {
257  return _constants;
258  }
259 
260 private:
261  // Integrals associated with the Form
262  FormIntegrals<T> _integrals;
263 
264  // Coefficients associated with the Form
265  FormCoefficients<T> _coefficients;
266 
267  // Constants associated with the Form
268  std::vector<
269  std::pair<std::string, std::shared_ptr<const function::Constant<T>>>>
270  _constants;
271 
272  // Function spaces (one for each argument)
273  std::vector<std::shared_ptr<const function::FunctionSpace>> _function_spaces;
274 
275  // The mesh (needed for functionals when we don't have any spaces)
276  std::shared_ptr<const mesh::Mesh> _mesh;
277 };
278 
279 } // namespace fem
280 } // namespace dolfinx
dolfinx::fem::Form::Form
Form(const std::vector< std::shared_ptr< const function::FunctionSpace >> &function_spaces, bool need_mesh_permutation_data)
Definition: Form.h:108
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:164
dolfinx::fem::IntegralType
IntegralType
Type of integral.
Definition: FormIntegrals.h:27
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:174
dolfinx::fem::Form::set_mesh
void set_mesh(const std::shared_ptr< const mesh::Mesh > &mesh)
Definition: Form.h:190
dolfinx::fem::Form::coefficients
const FormCoefficients< T > & coefficients() const
Access coefficients.
Definition: Form.h:233
dolfinx::fem::Form::constants
std::vector< std::pair< std::string, std::shared_ptr< const function::Constant< T > > > > & constants()
Access constants.
Definition: Form.h:244
dolfinx::fem::Form::rank
int rank() const
Rank of the form (bilinear form = 2, linear form = 1, functional = 0, etc)
Definition: Form.h:129
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. Names of the constants must agree with their names in UFL file.
Definition: Form.h:144
dolfinx::fem::Form::function_spaces
std::vector< std::shared_ptr< const function::FunctionSpace > > function_spaces() const
Return function spaces for all arguments.
Definition: Form.h:212
dolfinx::fem::Form::~Form
virtual ~Form()=default
Destructor.
dolfinx::fem::Form
Class for variational forms.
Definition: assembler.h:28
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:78
dolfinx::fem::Form::Form
Form(Form &&form)=default
Move constructor.
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:218
dolfinx::mesh::MeshTags
A MeshTags are used to associate mesh entities with values. The entity index (local to process) ident...
Definition: Form.h:36
dolfinx::fem::Form::Form
Form(const Form &form)=delete
Copy constructor.
dolfinx::function::Constant
A constant value which can be attached to a Form. Constants may be scalar (rank 0),...
Definition: Form.h:26
dolfinx::fem::Form::coefficients
FormCoefficients< T > & coefficients()
Access coefficients.
Definition: Form.h:230
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.
Definition: Form.h:134
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:255
dolfinx::fem::Form::integrals
const FormIntegrals< T > & integrals() const
Access form integrals.
Definition: Form.h:236
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:204
dolfinx::fem::Form::mesh
std::shared_ptr< const mesh::Mesh > mesh() const
Extract common mesh from form.
Definition: Form.h:199
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:39
dolfinx::function::FunctionSpace
This class represents a finite element function space defined by a mesh, a finite element,...
Definition: FunctionSpace.h:38
dolfinx::fem::FormCoefficients
Storage for the coefficients of a Form consisting of Function and the Element objects they are define...
Definition: FormCoefficients.h:34