DOLFIN-X
DOLFIN-X C++ interface
FormIntegrals.h
1 // Copyright (C) 2018 Chris Richardson
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 <array>
10 #include <functional>
11 #include <memory>
12 #include <petscsys.h>
13 #include <vector>
14 
15 namespace dolfinx
16 {
17 namespace mesh
18 {
19 class Mesh;
20 template <typename T>
21 class MeshTags;
22 } // namespace mesh
23 
24 namespace fem
25 {
26 
29 
31 {
32 public:
34  enum class Type : std::int8_t
35  {
36  cell = 0,
37  exterior_facet = 1,
38  interior_facet = 2,
39  vertex = 3
40  };
41 
43  FormIntegrals();
44 
50  const std::function<void(PetscScalar*, const PetscScalar*, const PetscScalar*,
51  const double*, const int*, const std::uint8_t*,
52  const std::uint32_t)>&
53  get_tabulate_tensor(FormIntegrals::Type type, int i) const;
54 
61  FormIntegrals::Type type, int i,
62  std::function<void(PetscScalar*, const PetscScalar*, const PetscScalar*,
63  const double*, const int*, const std::uint8_t*,
64  const std::uint32_t)>
65  fn);
66 
71 
77  std::vector<int> integral_ids(FormIntegrals::Type t) const;
78 
87  const std::vector<std::int32_t>& integral_domains(FormIntegrals::Type t,
88  int i) const;
89 
96  const mesh::MeshTags<int>& marker);
97 
103  void set_default_domains(const mesh::Mesh& mesh);
104 
105 private:
106  // Collect together the function, id, and indices of entities to
107  // integrate on
108  struct Integral
109  {
110  std::function<void(PetscScalar*, const PetscScalar*, const PetscScalar*,
111  const double*, const int*, const std::uint8_t*,
112  const std::uint32_t)>
113  tabulate;
114  int id;
115  std::vector<std::int32_t> active_entities;
116  };
117 
118  // Array of vectors of integrals, arranged by type (see Type enum, and
119  // struct Integral above)
120  std::array<std::vector<struct Integral>, 4> _integrals;
121 };
122 } // namespace fem
123 } // namespace dolfinx
dolfinx::fem::FormIntegrals::set_domains
void set_domains(FormIntegrals::Type type, const mesh::MeshTags< int > &marker)
Set the valid domains for the integrals of a given type from a MeshTags "marker". Note the MeshTags i...
Definition: FormIntegrals.cpp:89
dolfinx::fem::FormIntegrals::set_tabulate_tensor
void set_tabulate_tensor(FormIntegrals::Type type, int i, std::function< void(PetscScalar *, const PetscScalar *, const PetscScalar *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> fn)
Set the function for 'tabulate_tensor' for integral i of given type.
Definition: FormIntegrals.cpp:34
dolfinx::mesh::MeshTags
A MeshTags are used to associate mesh entities with values. The entity index (local to process) ident...
Definition: Form.h:35
dolfinx::mesh::Mesh
A Mesh consists of a set of connected and numbered mesh topological entities, and geometry data.
Definition: Mesh.h:46
dolfinx::fem::FormIntegrals::Type
Type
Type of integral.
Definition: FormIntegrals.h:34
dolfinx::fem::FormIntegrals::integral_ids
std::vector< int > integral_ids(FormIntegrals::Type t) const
Get the integer IDs of integrals of type t. The IDs correspond to the domains which the integrals are...
Definition: FormIntegrals.cpp:72
dolfinx::fem::FormIntegrals::integral_domains
const std::vector< std::int32_t > & integral_domains(FormIntegrals::Type t, int i) const
Get the list of active entities for the ith integral of type t. Note, these are not retrieved by ID,...
Definition: FormIntegrals.cpp:83
dolfinx::fem::FormIntegrals::FormIntegrals
FormIntegrals()
Construct empty object.
Definition: FormIntegrals.cpp:18
dolfinx::fem::FormIntegrals::num_integrals
int num_integrals(FormIntegrals::Type t) const
Number of integrals of given type.
Definition: FormIntegrals.cpp:67
dolfinx::fem::FormIntegrals::get_tabulate_tensor
const std::function< void(PetscScalar *, const PetscScalar *, const PetscScalar *, const double *, const int *, const std::uint8_t *, const std::uint32_t)> & get_tabulate_tensor(FormIntegrals::Type type, int i) const
Get the function for 'tabulate_tensor' for integral i of given type.
Definition: FormIntegrals.cpp:26
dolfinx::fem::FormIntegrals
Integrals of a Form, including those defined over cells, interior and exterior facets,...
Definition: FormIntegrals.h:30
dolfinx::fem::FormIntegrals::set_default_domains
void set_default_domains(const mesh::Mesh &mesh)
If there exists a default integral of any type, set the list of entities for those integrals from the...
Definition: FormIntegrals.cpp:179