1 // --------------------------------------------------------------------- 2 // 3 // Copyright (C) 2023 by the deal.II authors 4 // 5 // This file is part of the deal.II library. 6 // 7 // The deal.II library is free software; you can use it, redistribute 8 // it, and/or modify it under the terms of the GNU Lesser General 9 // Public License as published by the Free Software Foundation; either 10 // version 2.1 of the License, or (at your option) any later version. 11 // The full text of the license can be found in the file LICENSE.md at 12 // the top level directory of deal.II. 13 // 14 // Authors: Peter Munch, Martin Kronbichler 15 // 16 // --------------------------------------------------------------------- 17 18 #pragma once 19 #ifndef bps_h 20 # define bps_h 21 22 // deal.II includes 23 # include <deal.II/dofs/dof_tools.h> 24 25 # include <deal.II/fe/mapping.h> 26 27 # include <deal.II/lac/la_parallel_vector.h> 28 29 # include <deal.II/matrix_free/fe_evaluation.h> 30 # include <deal.II/matrix_free/matrix_free.h> 31 # include <deal.II/matrix_free/shape_info.h> 32 # include <deal.II/matrix_free/tools.h> 33 34 using namespace dealii; 35 36 37 38 /** 39 * BP types. For more details, see https://ceed.exascaleproject.org/bps/. 40 */ 41 enum class BPType : unsigned int 42 { 43 BP1, 44 BP2, 45 BP3, 46 BP4, 47 BP5, 48 BP6 49 }; 50 51 52 53 /** 54 * Struct storing relevant information regarding each BP. 55 */ 56 struct BPInfo 57 { 58 BPInfo(const BPType type, const int dim, const int fe_degree) 59 : type(type) 60 , dim(dim) 61 , fe_degree(fe_degree) 62 { 63 if (type == BPType::BP1) 64 type_string = "BP1"; 65 else if (type == BPType::BP2) 66 type_string = "BP2"; 67 else if (type == BPType::BP3) 68 type_string = "BP3"; 69 else if (type == BPType::BP4) 70 type_string = "BP4"; 71 else if (type == BPType::BP5) 72 type_string = "BP5"; 73 else if (type == BPType::BP6) 74 type_string = "BP6"; 75 76 this->n_q_points_1d = (type <= BPType::BP4) ? (fe_degree + 2) : (fe_degree + 1); 77 78 this->n_components = 79 (type == BPType::BP1 || type == BPType::BP3 || type == BPType::BP5) ? 1 : dim; 80 } 81 82 83 BPType type; 84 std::string type_string; 85 unsigned int dim; 86 unsigned int fe_degree; 87 unsigned int n_q_points_1d; 88 unsigned int n_components; 89 }; 90 91 92 93 /** 94 * Base class of operators. 95 */ 96 template <typename Number, typename MemorySpace> 97 class OperatorBase 98 { 99 public: 100 /** 101 * deal.II vector type 102 */ 103 using VectorType = LinearAlgebra::distributed::Vector<Number, MemorySpace>; 104 105 /** 106 * Initialize vector. 107 */ 108 virtual void 109 reinit() = 0; 110 111 /** 112 * Perform matrix-vector product 113 */ 114 virtual void 115 vmult(VectorType &dst, const VectorType &src) const = 0; 116 117 /** 118 * Initialize vector. 119 */ 120 virtual void 121 initialize_dof_vector(VectorType &vec) const = 0; 122 123 /** 124 * Compute inverse of diagonal. 125 */ 126 virtual void 127 compute_inverse_diagonal(VectorType &diagonal) const = 0; 128 }; 129 130 #endif 131