DOLFIN-X
DOLFIN-X C++ interface
utils.h
1 // Copyright (C) 2013-2020 Johan Hake, Jan Blechta and Garth N. Wells
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 "CoordinateElement.h"
10 #include "DofMap.h"
11 #include "ElementDofLayout.h"
12 #include <dolfinx/common/types.h>
13 #include <dolfinx/la/PETScMatrix.h>
14 #include <dolfinx/la/PETScVector.h>
15 #include <dolfinx/mesh/cell_types.h>
16 #include <memory>
17 #include <string>
18 #include <vector>
19 
20 struct ufc_dofmap;
21 struct ufc_form;
22 struct ufc_coordinate_mapping;
23 struct ufc_function_space;
24 
25 namespace dolfinx
26 {
27 namespace common
28 {
29 class IndexMap;
30 }
31 
32 namespace function
33 {
34 class Constant;
35 class Function;
36 class FunctionSpace;
37 } // namespace function
38 
39 namespace mesh
40 {
41 class Mesh;
42 class Topology;
43 } // namespace mesh
44 
45 namespace fem
46 {
47 class Form;
48 
57 std::array<std::vector<std::shared_ptr<const function::FunctionSpace>>, 2>
59  const Eigen::Ref<const Eigen::Array<const fem::Form*, Eigen::Dynamic,
60  Eigen::Dynamic, Eigen::RowMajor>>& a);
61 
65 la::PETScMatrix create_matrix(const Form& a);
66 
72 la::SparsityPattern create_sparsity_pattern(const Form& a);
73 
76 la::PETScMatrix create_matrix_block(
77  const Eigen::Ref<const Eigen::Array<const fem::Form*, Eigen::Dynamic,
78  Eigen::Dynamic, Eigen::RowMajor>>& a);
79 
81 la::PETScMatrix create_matrix_nest(
82  const Eigen::Ref<const Eigen::Array<const fem::Form*, Eigen::Dynamic,
83  Eigen::Dynamic, Eigen::RowMajor>>& a);
84 
86 la::PETScVector
87 create_vector_block(const std::vector<const common::IndexMap*>& maps);
88 
90 la::PETScVector
91 create_vector_nest(const std::vector<const common::IndexMap*>& maps);
92 
96 std::int64_t get_global_offset(const std::vector<const common::IndexMap*>& maps,
97  const int field, const std::int64_t index);
98 
100 ElementDofLayout create_element_dof_layout(const ufc_dofmap& dofmap,
101  const mesh::CellType cell_type,
102  const std::vector<int>& parent_map
103  = {});
104 
109 DofMap create_dofmap(MPI_Comm comm, const ufc_dofmap& dofmap,
110  mesh::Topology& topology);
111 
118 std::shared_ptr<Form> create_form(
119  ufc_form* (*fptr)(),
120  const std::vector<std::shared_ptr<const function::FunctionSpace>>& spaces);
121 
125 Form create_form(
126  const ufc_form& ufc_form,
127  const std::vector<std::shared_ptr<const function::FunctionSpace>>& spaces);
128 
130 std::vector<std::tuple<int, std::string, std::shared_ptr<function::Function>>>
131 get_coeffs_from_ufc_form(const ufc_form& ufc_form);
132 
134 std::vector<std::pair<std::string, std::shared_ptr<const function::Constant>>>
135 get_constants_from_ufc_form(const ufc_form& ufc_form);
136 
140 fem::CoordinateElement
141 create_coordinate_map(const ufc_coordinate_mapping& ufc_cmap);
142 
147 fem::CoordinateElement
148 create_coordinate_map(ufc_coordinate_mapping* (*fptr)());
149 
159 std::shared_ptr<function::FunctionSpace>
160 create_functionspace(ufc_function_space* (*fptr)(const char*),
161  const std::string function_name,
162  std::shared_ptr<mesh::Mesh> mesh);
163 
164 // NOTE: This is subject to change
166 Eigen::Array<PetscScalar, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor>
167 pack_coefficients(const fem::Form& form);
168 
169 // NOTE: This is subject to change
171 Eigen::Array<PetscScalar, Eigen::Dynamic, 1>
172 pack_constants(const fem::Form& form);
173 
174 } // namespace fem
175 } // namespace dolfinx
dolfinx::fem::create_form
std::shared_ptr< Form > create_form(ufc_form *(*fptr)(), const std::vector< std::shared_ptr< const function::FunctionSpace >> &spaces)
Create a form from a form_create function returning a pointer to a ufc_form, taking care of memory al...
Definition: utils.cpp:569
dolfinx::fem::get_coeffs_from_ufc_form
std::vector< std::tuple< int, std::string, std::shared_ptr< function::Function > > > get_coeffs_from_ufc_form(const ufc_form &ufc_form)
Extract coefficients from a UFC form.
Definition: utils.cpp:545
dolfinx::fem::create_vector_nest
la::PETScVector create_vector_nest(const std::vector< const common::IndexMap * > &maps)
Create nested (VecNest) vector. Vector is not zeroed.
Definition: utils.cpp:387
dolfinx::mesh::CellType
CellType
Cell type identifier.
Definition: cell_types.h:22
dolfinx::fem::create_sparsity_pattern
la::SparsityPattern create_sparsity_pattern(const Form &a)
Create a sparsity pattern for a given form. The pattern is not finalised, i.e. the caller is responsi...
Definition: utils.cpp:141
dolfinx::fem::create_coordinate_map
fem::CoordinateElement create_coordinate_map(const ufc_coordinate_mapping &ufc_cmap)
Create a CoordinateElement from ufc.
Definition: utils.cpp:671
dolfinx::fem::block_function_spaces
std::array< std::vector< std::shared_ptr< const function::FunctionSpace > >, 2 > block_function_spaces(const Eigen::Ref< const Eigen::Array< const fem::Form *, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> &a)
Extract FunctionSpaces for (0) rows blocks and (1) columns blocks from a rectangular array of bilinea...
Definition: utils.cpp:93
dolfinx::fem::create_functionspace
std::shared_ptr< function::FunctionSpace > create_functionspace(ufc_function_space *(*fptr)(const char *), const std::string function_name, std::shared_ptr< mesh::Mesh > mesh)
Create FunctionSpace from UFC.
Definition: utils.cpp:707
dolfinx::fem::create_matrix_nest
la::PETScMatrix create_matrix_nest(const Eigen::Ref< const Eigen::Array< const fem::Form *, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> &a)
Create nested (MatNest) matrix. Matrix is not zeroed.
Definition: utils.cpp:315
dolfinx::fem::create_vector_block
la::PETScVector create_vector_block(const std::vector< const common::IndexMap * > &maps)
Initialise monolithic vector. Vector is not zeroed.
Definition: utils.cpp:355
dolfinx::fem::get_global_offset
std::int64_t get_global_offset(const std::vector< const common::IndexMap * > &maps, const int field, const std::int64_t index)
Definition: utils.cpp:412
dolfinx::fem::create_element_dof_layout
ElementDofLayout create_element_dof_layout(const ufc_dofmap &dofmap, const mesh::CellType cell_type, const std::vector< int > &parent_map={})
Create an ElementDofLayout from a ufc_dofmap.
Definition: utils.cpp:444
dolfinx::fem::create_dofmap
DofMap create_dofmap(MPI_Comm comm, const ufc_dofmap &dofmap, mesh::Topology &topology)
Create dof map on mesh from a ufc_dofmap.
Definition: utils.cpp:514
dolfinx::fem::pack_coefficients
Eigen::Array< PetscScalar, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor > pack_coefficients(const fem::Form &form)
Pack form coefficients ready for assembly.
Definition: utils.cpp:726
dolfinx::fem::create_matrix
la::PETScMatrix create_matrix(const Form &a)
Create a matrix.
Definition: utils.cpp:194
dolfinx::fem::pack_constants
Eigen::Array< PetscScalar, Eigen::Dynamic, 1 > pack_constants(const fem::Form &form)
Pack form constants ready for assembly.
Definition: utils.cpp:782
dolfinx::fem::create_matrix_block
la::PETScMatrix create_matrix_block(const Eigen::Ref< const Eigen::Array< const fem::Form *, Eigen::Dynamic, Eigen::Dynamic, Eigen::RowMajor >> &a)
Initialise monolithic matrix for an array for bilinear forms. Matrix is not zeroed.
Definition: utils.cpp:210
dolfinx::fem::get_constants_from_ufc_form
std::vector< std::pair< std::string, std::shared_ptr< const function::Constant > > > get_constants_from_ufc_form(const ufc_form &ufc_form)
Extract coefficients from a UFC form.
Definition: utils.cpp:559