MFEM v4.7.0
Finite element discretization library
Loading...
Searching...
No Matches
solvers-atpmg.hpp
Go to the documentation of this file.
1// Copyright (c) 2010-2024, Lawrence Livermore National Security, LLC. Produced
2// at the Lawrence Livermore National Laboratory. All Rights reserved. See files
3// LICENSE and NOTICE for details. LLNL-CODE-806117.
4//
5// This file is part of the MFEM library. For more information and source code
6// availability visit https://mfem.org.
7//
8// MFEM is free software; you can redistribute it and/or modify it under the
9// terms of the BSD-3 license. We welcome feedback and contributions, see file
10// CONTRIBUTING.md for details.
11
12#ifndef MFEM_CEEDSOLVERS_ATPMG_H
13#define MFEM_CEEDSOLVERS_ATPMG_H
14
15#include "../interface/ceed.hpp"
16
17#ifdef MFEM_USE_CEED
18
19namespace mfem
20{
21
22namespace ceed
23{
24
25/** @brief Take given (high-order) CeedElemRestriction and make a new
26 CeedElemRestriction, which corresponds to a lower-order problem.
27
28 Assumes a Gauss-Lobatto basis and tensor product elements, and assumes that
29 the nodes in er_in are ordered in a tensor-product way.
30
31 This is a setup routine that operates on the host.
32
33 The caller is responsible for freeing er_out and dof_map. */
34int CeedATPMGElemRestriction(int order,
35 int order_reduction,
36 CeedElemRestriction er_in,
37 CeedElemRestriction* er_out,
38 CeedInt *&dof_map);
39
40/** @brief Create coarse-to-fine basis, given number of input nodes and order
41 reduction.
42
43 Assumes Gauss-Lobatto basis. This is useful because it does not require an
44 input CeedBasis object, which depends on choice of quadrature rule, whereas
45 the coarse-to-fine operator is independent of quadrature. */
46int CeedBasisATPMGCoarseToFine(Ceed ceed, int P1d, int dim, int order_reduction,
47 CeedBasis *basisc2f);
48
49/** @brief Given basis basisin, reduces its order by order_reduction and return
50 basisout (which has the same height (Q1d) but is narrower (smaller P1d))
51
52 The algorithm takes the locations of the fine nodes as input, but this
53 particular implementation simply assumes Gauss-Lobatto, and furthermore
54 assumes the MFEM [0, 1] reference element (rather than the Ceed/Petsc [-1,
55 1] element) */
56int CeedBasisATPMGCoarsen(CeedBasis basisin, CeedBasis* basisout,
57 CeedBasis* basis_ctof,
58 int order_reduction);
59
60/** @brief Coarsen a CeedOperator using semi-algebraic p-multigrid
61
62 This implementation does not coarsen the integration points at all.
63
64 @param[in] oper the operator to coarsen
65 @param[in] order_reduction how much to coarsen (order p)
66 @param[in] coarse_er CeedElemRestriction for coarse operator
67 (see CeedATPMGElemRestriction)
68 @param[out] coarse_basis_out CeedBasis for coarser operator
69 @param[out] basis_ctof_out CeedBasis describing interpolation from coarse to fine
70 @param[out] out coarsened CeedOperator
71*/
72int CeedATPMGOperator(CeedOperator oper, int order_reduction,
73 CeedElemRestriction coarse_er,
74 CeedBasis* coarse_basis_out,
75 CeedBasis* basis_ctof_out,
76 CeedOperator* out);
77
78/** @brief Given (fine) CeedOperator, produces everything you need for a coarse
79 level (operator and interpolation).
80
81 @param[in] oper Fine CeedOperator to coarsen
82 @param[in] order_reduction Amount to reduce the order (p) of the operator
83 @param[out] coarse_basis_out CeedBasis for coarse operator
84 @param[out] basis_ctof_out CeedBasis describing interpolation from coarse to fine
85 @param[out] er_out CeedElemRestriction for coarse operator
86 @param[out] coarse_oper coarse operator itself
87 @param[out] dof_map maps high-order ldof to low-order ldof, needed for
88 further coarsening
89*/
90int CeedATPMGBundle(CeedOperator oper, int order_reduction,
91 CeedBasis* coarse_basis_out,
92 CeedBasis* basis_ctof_out,
93 CeedElemRestriction* er_out,
94 CeedOperator* coarse_oper,
95 CeedInt *&dof_map);
96
97} // namespace ceed
98
99} // namespace mfem
100
101#endif // MFEM_USE_CEED
102
103#endif // include guard
int dim
Definition ex24.cpp:53
int CeedBasisATPMGCoarseToFine(Ceed ceed, int P1d, int dim, int order_reduction, CeedBasis *basisc2f)
Create coarse-to-fine basis, given number of input nodes and order reduction.
int CeedATPMGOperator(CeedOperator oper, int order_reduction, CeedElemRestriction coarse_er, CeedBasis coarse_basis_in, CeedBasis basis_ctof_in, CeedOperator *out)
int CeedBasisATPMGCoarsen(CeedBasis basisin, CeedBasis basisc2f, CeedBasis *basisout, int order_reduction)
int CeedATPMGBundle(CeedOperator oper, int order_reduction, CeedBasis *coarse_basis_out, CeedBasis *basis_ctof_out, CeedElemRestriction *er_out, CeedOperator *coarse_oper, CeedInt *&dof_map)
Given (fine) CeedOperator, produces everything you need for a coarse level (operator and interpolatio...
int CeedATPMGElemRestriction(int order, int order_reduction, CeedElemRestriction er_in, CeedElemRestriction *er_out, CeedInt *&dof_map)
Take given (high-order) CeedElemRestriction and make a new CeedElemRestriction, which corresponds to ...
OutStream out(std::cout)
Global stream used by the library for standard output. Initially it uses the same std::streambuf as s...
Definition globals.hpp:66