Reference documentation for deal.II version GIT relicensing-446-ge11b936273 2024-04-20 12:30:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Public Attributes | List of all members
MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData Struct Reference

#include <deal.II/matrix_free/matrix_free.h>

Public Types

enum  TasksParallelScheme { none = internal::MatrixFreeFunctions::TaskInfo::none , partition_partition , partition_color , color = internal::MatrixFreeFunctions::TaskInfo::color }
 
using MatrixFreeType = MatrixFree< dim, Number, VectorizedArrayType >
 

Public Member Functions

 AdditionalData (const TasksParallelScheme tasks_parallel_scheme=partition_partition, const unsigned int tasks_block_size=0, const UpdateFlags mapping_update_flags=update_gradients|update_JxW_values, const UpdateFlags mapping_update_flags_boundary_faces=update_default, const UpdateFlags mapping_update_flags_inner_faces=update_default, const UpdateFlags mapping_update_flags_faces_by_cells=update_default, const unsigned int mg_level=numbers::invalid_unsigned_int, const bool store_plain_indices=true, const bool initialize_indices=true, const bool initialize_mapping=true, const bool overlap_communication_computation=true, const bool hold_all_faces_to_owned_cells=false, const bool cell_vectorization_categories_strict=false, const bool allow_ghosted_vectors_in_loops=true)
 
 AdditionalData (const AdditionalData &other)
 
AdditionalDataoperator= (const AdditionalData &other)
 

Public Attributes

TasksParallelScheme tasks_parallel_scheme
 
unsigned int tasks_block_size
 
UpdateFlags mapping_update_flags
 
UpdateFlags mapping_update_flags_boundary_faces
 
UpdateFlags mapping_update_flags_inner_faces
 
UpdateFlags mapping_update_flags_faces_by_cells
 
unsigned int mg_level
 
bool store_plain_indices
 
bool initialize_indices
 
bool initialize_mapping
 
bool overlap_communication_computation
 
bool hold_all_faces_to_owned_cells
 
std::vector< unsigned intcell_vectorization_category
 
bool cell_vectorization_categories_strict
 
bool allow_ghosted_vectors_in_loops
 
MPI_Comm communicator_sm
 

Detailed Description

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
struct MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData

Collects the options for initialization of the MatrixFree class. The parameter tasks_parallel_scheme specifies the parallelization options in shared memory (task-based parallelism, where one can choose between no parallelism and three schemes that avoid that cells with access to the same vector entries are accessed simultaneously), and the parameter tasks_block_size the block size for task parallel scheduling. The parameters mapping_update_flags, mapping_update_flags_boundary_faces, mapping_update_flags_inner_faces, and mapping_update_flags_faces_by_cells specify the update flags that should be stored by this class.

The parameter mg_level specifies the level in the triangulation from which the indices are to be used. If the level is set to numbers::invalid_unsigned_int, the active cells are traversed, and otherwise the cells in the given level. This option has no effect in case a DoFHandler is given.

The parameter store_plain_indices indicates whether the DoFInfo class should also allow for access to vectors without resolving constraints.

The two parameters initialize_indices and initialize_mapping allow the user to disable some of the initialization processes. For example, if only the scheduling that avoids touching the same vector/matrix indices simultaneously is to be found, the mapping needs not be initialized. Likewise, if the mapping has changed from one iteration to the next but the topology has not (like when using a deforming mesh with MappingQEulerian), it suffices to initialize the mapping only.

The two parameters cell_vectorization_categories and cell_vectorization_categories_strict control the formation of batches for vectorization over several cells. It is used implicitly when working with hp-adaptivity but can also be useful in other contexts, such as in local time stepping where one would like to control which elements together form a batch of cells. The array cell_vectorization_categories is accessed by the number given by cell->active_cell_index() when working on the active cells with mg_level set to numbers::invalid_unsigned_int and by cell->index() for the level cells. By default, the different categories in cell_vectorization_category can be mixed and the algorithm is allowed to merge lower category numbers with the next higher categories if it is necessary inside the algorithm, in order to avoid partially filled SIMD lanes as much as possible. This gives a better utilization of the vectorization but might need special treatment, in particular for face integrals. If set to ‘true’, the algorithm will instead keep different categories separate and not mix them in a single vectorized array.

Finally, allow_ghosted_vectors_in_loops allows to enable and disable checks and communicator_sm gives the MPI communicator to be used if MPI-3.0 shared-memory features should be used.

Definition at line 183 of file matrix_free.h.

Member Typedef Documentation

◆ MatrixFreeType

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
using MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::MatrixFreeType = MatrixFree<dim, Number, VectorizedArrayType>

Provide the type of the surrounding MatrixFree class.

Definition at line 188 of file matrix_free.h.

Member Enumeration Documentation

◆ TasksParallelScheme

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
enum MatrixFree::AdditionalData::TasksParallelScheme

Collects options for task parallelism. See the documentation of the member variable MatrixFree::AdditionalData::tasks_parallel_scheme for a thorough description.

Enumerator
none 

Perform application in serial.

partition_partition 

Partition the cells into two levels and afterwards form chunks.

partition_color 

Partition on the global level and color cells within the partitions.

color 

Use the traditional coloring algorithm: this is like TasksParallelScheme::partition_color, but only uses one partition.

Definition at line 195 of file matrix_free.h.

Constructor & Destructor Documentation

◆ AdditionalData() [1/2]

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::AdditionalData ( const TasksParallelScheme  tasks_parallel_scheme = partition_partition,
const unsigned int  tasks_block_size = 0,
const UpdateFlags  mapping_update_flags = update_gradients | update_JxW_values,
const UpdateFlags  mapping_update_flags_boundary_faces = update_default,
const UpdateFlags  mapping_update_flags_inner_faces = update_default,
const UpdateFlags  mapping_update_flags_faces_by_cells = update_default,
const unsigned int  mg_level = numbers::invalid_unsigned_int,
const bool  store_plain_indices = true,
const bool  initialize_indices = true,
const bool  initialize_mapping = true,
const bool  overlap_communication_computation = true,
const bool  hold_all_faces_to_owned_cells = false,
const bool  cell_vectorization_categories_strict = false,
const bool  allow_ghosted_vectors_in_loops = true 
)
inline

Constructor for AdditionalData.

Definition at line 221 of file matrix_free.h.

◆ AdditionalData() [2/2]

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::AdditionalData ( const AdditionalData other)
inline

Copy constructor.

Definition at line 258 of file matrix_free.h.

Member Function Documentation

◆ operator=()

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
AdditionalData & MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::operator= ( const AdditionalData other)
inline

Copy assignment.

Definition at line 285 of file matrix_free.h.

Member Data Documentation

◆ tasks_parallel_scheme

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
TasksParallelScheme MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::tasks_parallel_scheme

Set the scheme for task parallelism. There are four options available. If set to none, the operator application is done in serial without shared memory parallelism. If this class is used together with MPI and MPI is also used for parallelism within the nodes, this flag should be set to none. The default value is partition_partition, i.e. we actually use multithreading with the first option described below.

The first option partition_partition is to partition the cells on two levels in onion-skin-like partitions and forming chunks of tasks_block_size after the partitioning. The partitioning finds sets of independent cells that enable working in parallel without accessing the same vector entries at the same time.

The second option partition_color is to use a partition on the global level and color cells within the partitions (where all chunks within a color are independent). Here, the subdivision into chunks of cells is done before the partitioning, which might give worse partitions but better cache performance if degrees of freedom are not renumbered.

The third option color is to use a traditional algorithm of coloring on the global level. This scheme is a special case of the second option where only one partition is present. Note that for problems with hanging nodes, there are quite many colors (50 or more in 3d), which might degrade parallel performance (bad cache behavior, many synchronization points).

Note
Threading support is currently experimental for the case inner face integrals are performed and it is recommended to use MPI parallelism if possible. While the scheme has been verified to work with the partition_partition option in case of usual DG elements, no comprehensive tests have been performed for systems of more general elements, like combinations of continuous and discontinuous elements that add face integrals to all terms.

Definition at line 347 of file matrix_free.h.

◆ tasks_block_size

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
unsigned int MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::tasks_block_size

Set the number of so-called cell batches that should form one partition. If zero size is given, the class tries to find a good size for the blocks based on MultithreadInfo::n_threads() and the number of cells present. Otherwise, the given number is used. If the given number is larger than one third of the number of total cells, this means no parallelism. Note that in the case vectorization is used, a cell batch consists of more than one physical cell.

Definition at line 358 of file matrix_free.h.

◆ mapping_update_flags

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
UpdateFlags MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::mapping_update_flags

This flag determines what data needs to be computed and cached on cells.

If your computations require operations like quadrature point locations or Hessians, these need to specified here (update_quadrature_points or update_hessians, respectively). Note that values, gradients, and Jacobian determinants (JxW values) are always computed regardless of the flags specified here.

Note that some additional flags might be set automatically (for example second derivatives might be evaluated on Cartesian cells since there the Jacobian describes the mapping completely).

Definition at line 373 of file matrix_free.h.

◆ mapping_update_flags_boundary_faces

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
UpdateFlags MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::mapping_update_flags_boundary_faces

This flag determines the mapping data on boundary faces to be cached. Note that MatrixFree uses a separate loop layout for face integrals in order to effectively vectorize also in the case of hanging nodes (which require different subface settings on the two sides) or some cells in the batch of a VectorizedArray of cells that are adjacent to the boundary and others that are not.

If set to a value different from update_general (default), the face information is explicitly built. Currently, MatrixFree supports to cache the following data on faces: inverse Jacobians, Jacobian determinants (JxW), quadrature points, data for Hessians (derivative of Jacobians), and normal vectors.

Note
In order to be able to perform a boundary_face_operation in the MatrixFree::loop(), this field must be set to a value different from UpdateFlags::update_default.

Definition at line 393 of file matrix_free.h.

◆ mapping_update_flags_inner_faces

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
UpdateFlags MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::mapping_update_flags_inner_faces

This flag determines the mapping data on interior faces to be cached. Note that MatrixFree uses a separate loop layout for face integrals in order to effectively vectorize also in the case of hanging nodes (which require different subface settings on the two sides) or some cells in the batch of a VectorizedArray of cells that are adjacent to the boundary and others that are not.

If set to a value different from update_general (default), the face information is explicitly built. Currently, MatrixFree supports to cache the following data on faces: inverse Jacobians, Jacobian determinants (JxW), quadrature points, data for Hessians (derivative of Jacobians), and normal vectors.

Note
In order to be able to perform a inner_face_operation in the MatrixFree::loop(), this field must be set to a value different from UpdateFlags::update_default.

Definition at line 413 of file matrix_free.h.

◆ mapping_update_flags_faces_by_cells

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
UpdateFlags MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::mapping_update_flags_faces_by_cells

This flag determines the mapping data for faces in a different layout with respect to vectorizations. Whereas mapping_update_flags_inner_faces and mapping_update_flags_boundary_faces trigger building the data in a face-centric way with proper vectorization, the current data field attaches the face information to the cells and their way of vectorization. This is only needed in special situations, as for example for block-Jacobi methods where the full operator to a cell including its faces are evaluated. This data is accessed by FEFaceEvaluation::reinit(cell_batch_index, face_number). However, currently no coupling terms to neighbors can be computed with this approach because the neighbors are not laid out by the VectorizedArray data layout with an array-of-struct-of-array-type data structures.

Note that you should only compute this data field in case you really need it as it more than doubles the memory required by the mapping data on faces.

If set to a value different from update_general (default), the face information is explicitly built. Currently, MatrixFree supports to cache the following data on faces: inverse Jacobians, Jacobian determinants (JxW), quadrature points, data for Hessians (derivative of Jacobians), and normal vectors.

Definition at line 441 of file matrix_free.h.

◆ mg_level

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
unsigned int MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::mg_level

This option can be used to define whether we work on a certain level of the mesh, and not the active cells. If set to invalid_unsigned_int (which is the default value), the active cells are gone through, otherwise the level given by this parameter. Note that if you specify to work on a level, its dofs must be distributed by using dof_handler.distribute_mg_dofs(fe);.

Definition at line 451 of file matrix_free.h.

◆ store_plain_indices

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
bool MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::store_plain_indices

Controls whether to enable reading from vectors without resolving constraints, i.e., just read the local values of the vector. By default, this option is enabled. In case you want to use FEEvaluationBase::read_dof_values_plain, this flag needs to be set.

Definition at line 459 of file matrix_free.h.

◆ initialize_indices

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
bool MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::initialize_indices

Option to control whether the indices stored in the DoFHandler should be read and the pattern for task parallelism should be set up in the initialize method of MatrixFree. The default value is true. Can be disabled in case the mapping should be recomputed (e.g. when using a deforming mesh described through MappingEulerian) but the topology of cells has remained the same.

Definition at line 470 of file matrix_free.h.

◆ initialize_mapping

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
bool MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::initialize_mapping

Option to control whether the mapping information should be computed in the initialize method of MatrixFree. The default value is true. Can be disabled when only some indices should be set up (e.g. when only a set of independent cells should be computed).

Definition at line 479 of file matrix_free.h.

◆ overlap_communication_computation

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
bool MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::overlap_communication_computation

Option to control whether the loops should overlap communications and computations as far as possible in case the vectors passed to the loops support non-blocking data exchange. In most situations, overlapping is faster in case the amount of data to be sent is more than a few kilobytes. If less data is sent, the communication is latency bound on most clusters (point-to-point latency is around 1 microsecond on good clusters by 2016 standards). Depending on the MPI implementation and the fabric, it may be faster to not overlap and wait for the data to arrive. The default is true, i.e., communication and computation are overlapped.

Definition at line 493 of file matrix_free.h.

◆ hold_all_faces_to_owned_cells

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
bool MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::hold_all_faces_to_owned_cells

By default, the face part will only hold those faces (and ghost elements behind faces) that are going to be processed locally. In case MatrixFree should have access to all neighbors on locally owned cells, this option enables adding the respective faces at the end of the face range.

Definition at line 502 of file matrix_free.h.

◆ cell_vectorization_category

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
std::vector<unsigned int> MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::cell_vectorization_category

This data structure allows to assign a fraction of cells to different categories when building the information for vectorization. It is used implicitly when working with hp-adaptivity (with each active index being a category) but can also be useful in other contexts where one would like to control which cells together can form a batch of cells. Such an example is "local time stepping", where cells of different categories progress with different time-step sizes and, as a consequence, can only processed together with cells with the same category.

This array is accessed by the number given by cell->active_cell_index() when working on the active cells (with mg_level set to numbers::invalid_unsigned_int) and by cell->index() for the level cells.

Note
This field is empty upon construction of AdditionalData. It is the responsibility of the user to resize this field to triangulation.n_active_cells() or triangulation.n_cells(level) when filling data.

Definition at line 525 of file matrix_free.h.

◆ cell_vectorization_categories_strict

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
bool MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::cell_vectorization_categories_strict

By default, the different categories in cell_vectorization_category can be mixed and the algorithm is allowed to merge lower categories with the next higher categories if it is necessary inside the algorithm. This gives a better utilization of the vectorization but might need special treatment, in particular for face integrals. If set to true, the algorithm will instead keep different categories separate and not mix them in a single vectorized array.

Definition at line 536 of file matrix_free.h.

◆ allow_ghosted_vectors_in_loops

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
bool MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::allow_ghosted_vectors_in_loops

Assert that vectors passed to the MatrixFree loops are not ghosted. This variable is primarily intended to reveal bugs or performance problems caused by vectors that are involuntarily in ghosted mode, by adding a check that this is not the case. In terms of correctness, the MatrixFree::loop() and MatrixFree::cell_loop() methods support both cases and perform similar operations. In particular, ghost values are always updated on the source vector within the loop, and the difference is only in whether the initial non-ghosted state is restored.

Definition at line 548 of file matrix_free.h.

◆ communicator_sm

template<int dim, typename Number = double, typename VectorizedArrayType = VectorizedArray<Number>>
MPI_Comm MatrixFree< dim, Number, VectorizedArrayType >::AdditionalData::communicator_sm

Shared-memory MPI communicator. Default: MPI_COMM_SELF.

Definition at line 553 of file matrix_free.h.


The documentation for this struct was generated from the following file: