| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | // This file is part of Bembel, the higher order C++ boundary element library. | ||
| 2 | // | ||
| 3 | // Copyright (C) 2022 see <http://www.bembel.eu> | ||
| 4 | // | ||
| 5 | // It was written as part of a cooperation of J. Doelz, H. Harbrecht, S. Kurz, | ||
| 6 | // M. Multerer, S. Schoeps, and F. Wolf at Technische Universitaet Darmstadt, | ||
| 7 | // Universitaet Basel, and Universita della Svizzera italiana, Lugano. This | ||
| 8 | // source code is subject to the GNU General Public License version 3 and | ||
| 9 | // provided WITHOUT ANY WARRANTY, see <http://www.bembel.eu> for further | ||
| 10 | // information. | ||
| 11 | #ifndef BEMBEL_SRC_ANSATZSPACE_SUPERSPACE_HPP_ | ||
| 12 | #define BEMBEL_SRC_ANSATZSPACE_SUPERSPACE_HPP_ | ||
| 13 | namespace Bembel { | ||
| 14 | /** | ||
| 15 | * \ingroup AnsatzSpace | ||
| 16 | * \brief The superspace manages local polynomial bases on each element of the | ||
| 17 | * mesh and provides an interface to evaluate them. | ||
| 18 | * | ||
| 19 | * | ||
| 20 | */ | ||
| 21 | template <typename Derived> | ||
| 22 | struct SuperSpace { | ||
| 23 | typedef typename LinearOperatorTraits<Derived>::Scalar Scalar; | ||
| 24 | ////////////////////////////////////////////////////////////////////////////// | ||
| 25 | // constructors | ||
| 26 | ////////////////////////////////////////////////////////////////////////////// | ||
| 27 | /** | ||
| 28 | * \brief Default constructor for the SuperSpace class. | ||
| 29 | * | ||
| 30 | * This constructor creates a SuperSpace object with default parameters. | ||
| 31 | */ | ||
| 32 | 868 | SuperSpace() {} | |
| 33 | /** | ||
| 34 | * \brief Parameterized constructor for the SuperSpace class. | ||
| 35 | * | ||
| 36 | * This constructor initializes a SuperSpace object with the provided | ||
| 37 | * parameters. | ||
| 38 | * | ||
| 39 | * \param geom The geometry object defining the space. | ||
| 40 | * \param M The refinement level of the space. | ||
| 41 | * \param P The degree of polynomials used in the space. | ||
| 42 | */ | ||
| 43 |
1/2✓ Branch 2 taken 24 times.
✗ Branch 3 not taken.
|
48 | SuperSpace(Geometry& geom, int M, int P) { init_SuperSpace(geom, M, P); } |
| 44 | /** | ||
| 45 | * \brief Copy constructor for the SuperSpace class. | ||
| 46 | * | ||
| 47 | * This constructor initializes a SuperSpace object by copying another | ||
| 48 | * SuperSpace object. | ||
| 49 | * | ||
| 50 | * \param other The SuperSpace object to copy from. | ||
| 51 | */ | ||
| 52 | 617 | SuperSpace(const SuperSpace& other) { | |
| 53 | 617 | mesh_ = other.mesh_; | |
| 54 | 617 | phi = other.phi; | |
| 55 | 617 | phiDx = other.phiDx; | |
| 56 | 617 | phiPhi = other.phiPhi; | |
| 57 | 617 | phiPhiDx = other.phiPhiDx; | |
| 58 | 617 | phiPhiDy = other.phiPhiDy; | |
| 59 | 617 | phiTimesPhi = other.phiTimesPhi; | |
| 60 | // vPhiScalVPhi = other.vPhiScalVPhi; | ||
| 61 | 617 | divPhiTimesDivPhi = other.divPhiTimesDivPhi; | |
| 62 | 617 | polynomial_degree = other.polynomial_degree; | |
| 63 | 617 | polynomial_degree_plus_one_squared = | |
| 64 | 617 | other.polynomial_degree_plus_one_squared; | |
| 65 | 617 | } | |
| 66 | /** | ||
| 67 | * \brief Move constructor for the SuperSpace class. | ||
| 68 | * | ||
| 69 | * This constructor initializes a SuperSpace object by moving from another | ||
| 70 | * SuperSpace object. | ||
| 71 | * | ||
| 72 | * \param other The SuperSpace object to move from. | ||
| 73 | */ | ||
| 74 | SuperSpace(SuperSpace&& other) { | ||
| 75 | mesh_ = other.mesh_; | ||
| 76 | phi = other.phi; | ||
| 77 | phiDx = other.phiDx; | ||
| 78 | phiPhi = other.phiPhi; | ||
| 79 | phiPhiDx = other.phiPhiDx; | ||
| 80 | phiPhiDy = other.phiPhiDy; | ||
| 81 | phiTimesPhi = other.phiTimesPhi; | ||
| 82 | // vPhiScalVPhi = other.vPhiScalVPhi; | ||
| 83 | divPhiTimesDivPhi = other.divPhiTimesDivPhi; | ||
| 84 | polynomial_degree = other.polynomial_degree; | ||
| 85 | polynomial_degree_plus_one_squared = | ||
| 86 | other.polynomial_degree_plus_one_squared; | ||
| 87 | } | ||
| 88 | /** | ||
| 89 | * \brief Assignment operator for the SuperSpace class. | ||
| 90 | * | ||
| 91 | * This operator assigns the contents of another SuperSpace object to this | ||
| 92 | * one. | ||
| 93 | * | ||
| 94 | * \param other The SuperSpace object to copy from. | ||
| 95 | * \return A reference to the updated SuperSpace object. | ||
| 96 | */ | ||
| 97 | 596 | SuperSpace& operator=(SuperSpace other) { | |
| 98 | 596 | mesh_ = other.mesh_; | |
| 99 | 596 | phi = other.phi; | |
| 100 | 596 | phiDx = other.phiDx; | |
| 101 | 596 | phiPhi = other.phiPhi; | |
| 102 | 596 | phiPhiDx = other.phiPhiDx; | |
| 103 | 596 | phiPhiDy = other.phiPhiDy; | |
| 104 | 596 | phiTimesPhi = other.phiTimesPhi; | |
| 105 | // vPhiScalVPhi = other.vPhiScalVPhi; | ||
| 106 | 596 | divPhiTimesDivPhi = other.divPhiTimesDivPhi; | |
| 107 | 596 | polynomial_degree = other.polynomial_degree; | |
| 108 | 596 | polynomial_degree_plus_one_squared = | |
| 109 | 596 | other.polynomial_degree_plus_one_squared; | |
| 110 | 596 | return *this; | |
| 111 | } | ||
| 112 | ////////////////////////////////////////////////////////////////////////////// | ||
| 113 | // getters | ||
| 114 | ////////////////////////////////////////////////////////////////////////////// | ||
| 115 | 10836533 | int get_polynomial_degree() const { return polynomial_degree; } | |
| 116 | int get_polynomial_degree_plus_one_squared() const { | ||
| 117 | return polynomial_degree_plus_one_squared; | ||
| 118 | } | ||
| 119 | 2572209 | int get_refinement_level() const { return mesh_->get_max_level(); } | |
| 120 | int get_number_of_elements() const { return mesh_->get_number_of_elements(); } | ||
| 121 | 1944 | int get_number_of_patches() const { return mesh_->get_geometry().size(); } | |
| 122 | 7680 | const PatchVector& get_geometry() const { return mesh_->get_geometry(); } | |
| 123 | 55640 | const ClusterTree& get_mesh() const { return *mesh_; } | |
| 124 | ////////////////////////////////////////////////////////////////////////////// | ||
| 125 | // init_SuperSpace | ||
| 126 | ////////////////////////////////////////////////////////////////////////////// | ||
| 127 | 296 | void init_SuperSpace(const Geometry& geom, int M, int P) { | |
| 128 | 296 | polynomial_degree = P; | |
| 129 | 296 | polynomial_degree_plus_one_squared = | |
| 130 | 296 | (polynomial_degree + 1) * (polynomial_degree + 1); | |
| 131 | 296 | phi = (Basis::BasisHandler<Scalar>::funPtrPhi(P)); | |
| 132 | 296 | phiDx = (Basis::BasisHandler<Scalar>::funPtrPhiDx(P)); | |
| 133 | 296 | phiPhi = (Basis::BasisHandler<Scalar>::funPtrPhiPhi(P)); | |
| 134 | 296 | phiPhiDx = (Basis::BasisHandler<Scalar>::funPtrPhiPhiDx(P)); | |
| 135 | 296 | phiPhiDy = (Basis::BasisHandler<Scalar>::funPtrPhiPhiDy(P)); | |
| 136 | 296 | phiTimesPhi = (Basis::BasisHandler<Scalar>::funPtrPhiTimesPhi(P)); | |
| 137 | // vPhiScalVPhi = (Basis::BasisHandler<typename | ||
| 138 | // LinearOperatorTraits<Derived>::Scalar>::funPtrVPhiScalVPhi(P)); | ||
| 139 | 296 | divPhiTimesDivPhi = | |
| 140 | 296 | (Basis::BasisHandler<Scalar>::funPtrDivPhiTimesDivPhi(P)); | |
| 141 | 296 | mesh_ = std::make_shared<ClusterTree>(); | |
| 142 | 296 | mesh_->init_ClusterTree(geom, M); | |
| 143 | 296 | mesh_->checkOrientation(); | |
| 144 | 296 | return; | |
| 145 | } | ||
| 146 | ////////////////////////////////////////////////////////////////////////////// | ||
| 147 | // map2surface | ||
| 148 | ////////////////////////////////////////////////////////////////////////////// | ||
| 149 | /** | ||
| 150 | * \brief Evaluation of a point in the element and its Jacobian matrix. | ||
| 151 | * | ||
| 152 | * This function performs the affine transformation of an element to the | ||
| 153 | * reference domain of the patch and returns the output in a surface point. | ||
| 154 | * | ||
| 155 | * \param e : Element to be evaluated, | ||
| 156 | * \param xi : Point in [0, 1]^2 of the element | ||
| 157 | * \param w : Quadrature weight | ||
| 158 | * \param surf_pt : Evaluated point and its jacobian | ||
| 159 | */ | ||
| 160 | 99820505 | void map2surface(const ElementTreeNode& e, const Eigen::Vector2d& xi, | |
| 161 | double w, SurfacePoint* surf_pt) const { | ||
| 162 |
3/6✓ Branch 2 taken 99820505 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 99820505 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 99820505 times.
✗ Branch 9 not taken.
|
99820505 | Eigen::Vector2d st = e.llc_ + e.get_h() * xi; |
| 163 |
1/2✓ Branch 4 taken 99820505 times.
✗ Branch 5 not taken.
|
99820505 | mesh_->get_geometry()[e.patch_].updateSurfacePoint(surf_pt, st, w, xi); |
| 164 | 199641010 | return; | |
| 165 | } | ||
| 166 | ////////////////////////////////////////////////////////////////////////////// | ||
| 167 | // Methods | ||
| 168 | ////////////////////////////////////////////////////////////////////////////// | ||
| 169 | /** | ||
| 170 | * \brief Compute all products of local shape functions on the unit square at | ||
| 171 | * coordinates s,t, scale by w and add to intval. | ||
| 172 | */ | ||
| 173 | 8189792 | void addScaledBasisInteraction( | |
| 174 | Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic>* intval, | ||
| 175 | typename LinearOperatorTraits<Derived>::Scalar w, | ||
| 176 | const Eigen::Vector2d& s, const Eigen::Vector2d& t) const { | ||
| 177 |
2/4✓ Branch 2 taken 8189792 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 8189792 times.
✗ Branch 6 not taken.
|
8189792 | phiTimesPhi(intval, w, s, t); |
| 178 | 8189792 | } | |
| 179 | /** | ||
| 180 | * \brief Compute all products of local shape functions on the unit square at | ||
| 181 | * coordinates s,t. | ||
| 182 | */ | ||
| 183 | 2563584 | Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> basisInteraction( | |
| 184 | const Eigen::Vector2d& s, const Eigen::Vector2d& t) const { | ||
| 185 | 2563584 | Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> intval( | |
| 186 | 2563584 | polynomial_degree_plus_one_squared, polynomial_degree_plus_one_squared); | |
| 187 |
1/2✓ Branch 1 taken 2563584 times.
✗ Branch 2 not taken.
|
2563584 | intval.setZero(); |
| 188 |
3/6✓ Branch 1 taken 2563584 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2563584 times.
✗ Branch 5 not taken.
✓ Branch 8 taken 2563584 times.
✗ Branch 9 not taken.
|
2563584 | phiTimesPhi(&intval, 1., s, t); |
| 189 | 2563584 | return intval; | |
| 190 | ✗ | } | |
| 191 | |||
| 192 | /** | ||
| 193 | * \brief Compute all products of surface curls of local shape functions | ||
| 194 | * on the unit square at coordinates s,t. | ||
| 195 | */ | ||
| 196 | 54 | void addScaledSurfaceCurlInteraction( | |
| 197 | Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic>* intval, Scalar w, | ||
| 198 | const SurfacePoint& p1, const SurfacePoint& p2) const { | ||
| 199 | // surface measures | ||
| 200 |
4/8✓ Branch 1 taken 54 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 54 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 54 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 54 times.
✗ Branch 11 not taken.
|
54 | double kappa1 = p1.segment<3>(6).cross(p1.segment<3>(9)).norm(); |
| 201 |
4/8✓ Branch 1 taken 54 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 54 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 54 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 54 times.
✗ Branch 11 not taken.
|
54 | double kappa2 = p2.segment<3>(6).cross(p2.segment<3>(9)).norm(); |
| 202 | // compute basis functions's surface curl. Each column of s_curl is a basis | ||
| 203 | // function's surface curl at point s. | ||
| 204 |
4/8✓ Branch 1 taken 54 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 54 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 54 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 54 times.
✗ Branch 11 not taken.
|
108 | Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> s_curl = |
| 205 |
1/2✓ Branch 1 taken 54 times.
✗ Branch 2 not taken.
|
54 | (1.0 / kappa1) * |
| 206 |
5/10✓ Branch 1 taken 54 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 54 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 54 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 54 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 54 times.
✗ Branch 14 not taken.
|
108 | (-p1.segment<3>(6) * basisDy(p1.segment<2>(0)).transpose() + |
| 207 |
6/12✓ Branch 1 taken 54 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 54 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 54 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 54 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 54 times.
✗ Branch 14 not taken.
✓ Branch 16 taken 54 times.
✗ Branch 17 not taken.
|
108 | p1.segment<3>(9) * basisDx(p1.segment<2>(0)).transpose()); |
| 208 |
4/8✓ Branch 1 taken 54 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 54 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 54 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 54 times.
✗ Branch 11 not taken.
|
108 | Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> t_curl = |
| 209 |
1/2✓ Branch 1 taken 54 times.
✗ Branch 2 not taken.
|
54 | (1.0 / kappa2) * |
| 210 |
5/10✓ Branch 1 taken 54 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 54 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 54 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 54 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 54 times.
✗ Branch 14 not taken.
|
108 | (-p2.segment<3>(6) * basisDy(p2.segment<2>(0)).transpose() + |
| 211 |
6/12✓ Branch 1 taken 54 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 54 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 54 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 54 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 54 times.
✗ Branch 14 not taken.
✓ Branch 16 taken 54 times.
✗ Branch 17 not taken.
|
108 | p2.segment<3>(9) * basisDx(p2.segment<2>(0)).transpose()); |
| 212 | // inner product of surface curls of any two basis functions | ||
| 213 |
2/2✓ Branch 0 taken 216 times.
✓ Branch 1 taken 54 times.
|
270 | for (int j = 0; j < polynomial_degree_plus_one_squared; ++j) |
| 214 |
2/2✓ Branch 0 taken 864 times.
✓ Branch 1 taken 216 times.
|
1080 | for (int i = 0; i < polynomial_degree_plus_one_squared; ++i) |
| 215 | 1728 | (*intval)(j * polynomial_degree_plus_one_squared + i) += | |
| 216 |
4/8✓ Branch 1 taken 864 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 864 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 864 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 864 times.
✗ Branch 11 not taken.
|
864 | w * s_curl.col(i).dot(t_curl.col(j)); |
| 217 | 54 | } | |
| 218 | |||
| 219 | /** | ||
| 220 | * \brief Compute all products of surface gradients of local shape functions | ||
| 221 | * on the unit square at coordinates s,t. | ||
| 222 | */ | ||
| 223 | 54 | void addScaledSurfaceGradientInteraction( | |
| 224 | Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic>* intval, Scalar w, | ||
| 225 | const SurfacePoint& p1, const SurfacePoint& p2) const { | ||
| 226 | // inner product of surface gradients of any two basis functions equals to | ||
| 227 | // inner product of surface curls of any two basis functions | ||
| 228 | 54 | addScaledSurfaceCurlInteraction(intval, w, p1, p2); | |
| 229 | 54 | } | |
| 230 | |||
| 231 | /** | ||
| 232 | * \brief Compute all scalar products of vector valued local shape functions | ||
| 233 | * on the surface points with reference coordinates s,t, scale by w and add to | ||
| 234 | * intval. | ||
| 235 | */ | ||
| 236 | 2563584 | void addScaledVectorBasisInteraction( | |
| 237 | Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic>* intval, Scalar w, | ||
| 238 | const Eigen::Vector2d& s, const Eigen::Vector2d& t, | ||
| 239 | const Eigen::Vector3d x_f_dx, const Eigen::Vector3d x_f_dy, | ||
| 240 | const Eigen::Vector3d y_f_dx, const Eigen::Vector3d y_f_dy) const { | ||
| 241 |
1/2✓ Branch 1 taken 2563584 times.
✗ Branch 2 not taken.
|
2563584 | auto basis_interaction = basisInteraction(s, t); |
| 242 | ✗ | intval->block(0, 0, polynomial_degree_plus_one_squared, | |
| 243 |
3/6✓ Branch 1 taken 2563584 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2563584 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 2563584 times.
✗ Branch 8 not taken.
|
2563584 | polynomial_degree_plus_one_squared) += |
| 244 |
1/2✓ Branch 1 taken 2563584 times.
✗ Branch 2 not taken.
|
2563584 | w * x_f_dx.dot(y_f_dx) * basis_interaction; |
| 245 | ✗ | intval->block(0, polynomial_degree_plus_one_squared, | |
| 246 | 2563584 | polynomial_degree_plus_one_squared, | |
| 247 |
3/6✓ Branch 1 taken 2563584 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2563584 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 2563584 times.
✗ Branch 8 not taken.
|
2563584 | polynomial_degree_plus_one_squared) += |
| 248 |
1/2✓ Branch 1 taken 2563584 times.
✗ Branch 2 not taken.
|
2563584 | w * x_f_dx.dot(y_f_dy) * basis_interaction; |
| 249 | ✗ | intval->block(polynomial_degree_plus_one_squared, 0, | |
| 250 | 2563584 | polynomial_degree_plus_one_squared, | |
| 251 |
3/6✓ Branch 1 taken 2563584 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2563584 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 2563584 times.
✗ Branch 8 not taken.
|
2563584 | polynomial_degree_plus_one_squared) += |
| 252 |
1/2✓ Branch 1 taken 2563584 times.
✗ Branch 2 not taken.
|
2563584 | w * x_f_dy.dot(y_f_dx) * basis_interaction; |
| 253 | ✗ | intval->block(polynomial_degree_plus_one_squared, | |
| 254 | 2563584 | polynomial_degree_plus_one_squared, | |
| 255 | 2563584 | polynomial_degree_plus_one_squared, | |
| 256 |
3/6✓ Branch 1 taken 2563584 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2563584 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 2563584 times.
✗ Branch 8 not taken.
|
2563584 | polynomial_degree_plus_one_squared) += |
| 257 |
1/2✓ Branch 1 taken 2563584 times.
✗ Branch 2 not taken.
|
2563584 | w * x_f_dy.dot(y_f_dy) * basis_interaction; |
| 258 | 2563584 | } | |
| 259 | /** | ||
| 260 | * \brief Compute all scalar products of vector valued local shape functions | ||
| 261 | * on the surface points with reference coordinates s,t. | ||
| 262 | */ | ||
| 263 | Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> vectorBasisInteraction( | ||
| 264 | const Eigen::Vector2d& s, const Eigen::Vector2d& t, | ||
| 265 | const Eigen::Vector3d x_f_dx, const Eigen::Vector3d x_f_dy, | ||
| 266 | const Eigen::Vector3d y_f_dx, const Eigen::Vector3d y_f_dy) const { | ||
| 267 | Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> intval( | ||
| 268 | 2 * polynomial_degree_plus_one_squared, | ||
| 269 | 2 * polynomial_degree_plus_one_squared); | ||
| 270 | intval.setZero(); | ||
| 271 | addScaledVectorBasisInteraction(&intval, 1., s, t, x_f_dx, x_f_dy, y_f_dx, | ||
| 272 | y_f_dy); | ||
| 273 | return intval; | ||
| 274 | } | ||
| 275 | /** | ||
| 276 | * \brief Compute all products of divergences of local shape functions on the | ||
| 277 | * unit square at coordinates s,t, scale by w and add to intval. | ||
| 278 | */ | ||
| 279 | 2563584 | void addScaledVectorBasisDivergenceInteraction( | |
| 280 | Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic>* intval, Scalar w, | ||
| 281 | const Eigen::Vector2d& s, const Eigen::Vector2d& t) const { | ||
| 282 |
2/4✓ Branch 2 taken 2563584 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 2563584 times.
✗ Branch 6 not taken.
|
2563584 | divPhiTimesDivPhi(intval, w, s, t); |
| 283 | 2563584 | } | |
| 284 | /** | ||
| 285 | * \brief Compute all products of divergences of local shape functions on the | ||
| 286 | * unit square at coordinates s,t. | ||
| 287 | */ | ||
| 288 | Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> | ||
| 289 | vectorBasisDivergenceInteraction(const Eigen::Vector2d& s, | ||
| 290 | const Eigen::Vector2d& t) const { | ||
| 291 | Eigen::Matrix<Scalar, Eigen::Dynamic, Eigen::Dynamic> intval( | ||
| 292 | 2 * polynomial_degree_plus_one_squared, | ||
| 293 | 2 * polynomial_degree_plus_one_squared); | ||
| 294 | intval.setZero(); | ||
| 295 | divPhiTimesDivPhi(&intval, 1., s, t); | ||
| 296 | return intval; | ||
| 297 | } | ||
| 298 | /** | ||
| 299 | * \brief Evaluate local shape functions on the unit square at coordinate s, | ||
| 300 | * scale by w and add to intval. | ||
| 301 | */ | ||
| 302 | 1410 | void addScaledBasis(Eigen::Matrix<Scalar, Eigen::Dynamic, 1>* intval, | |
| 303 | Scalar w, const Eigen::Vector2d& s) const { | ||
| 304 |
1/2✓ Branch 2 taken 1410 times.
✗ Branch 3 not taken.
|
1410 | phiPhi(intval, w, s); |
| 305 | 1410 | } | |
| 306 | /** | ||
| 307 | * \brief Evaluate local shape functions on the unit square at coordinate s. | ||
| 308 | */ | ||
| 309 | 1003145 | Eigen::Matrix<Scalar, Eigen::Dynamic, 1> basis( | |
| 310 | const Eigen::Vector2d& s) const { | ||
| 311 | 1003145 | Eigen::Matrix<Scalar, Eigen::Dynamic, 1> intval( | |
| 312 | 1003145 | polynomial_degree_plus_one_squared); | |
| 313 |
1/2✓ Branch 1 taken 999305 times.
✗ Branch 2 not taken.
|
1003145 | intval.setZero(); |
| 314 |
3/5✓ Branch 1 taken 999305 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 868186 times.
✓ Branch 5 taken 131119 times.
✗ Branch 6 not taken.
|
1003145 | phiPhi(&intval, 1., s); |
| 315 | 1003145 | return intval; | |
| 316 | ✗ | } | |
| 317 | /** | ||
| 318 | * \brief Evaluate derivatives in x direction of local shape functions on the | ||
| 319 | * unit square at coordinate s, scale by w and add to intval. | ||
| 320 | */ | ||
| 321 | void addScaledBasisDx(Eigen::Matrix<Scalar, Eigen::Dynamic, 1>* intval, | ||
| 322 | typename LinearOperatorTraits<Derived>::Scalar w, | ||
| 323 | const Eigen::Vector2d& s) const { | ||
| 324 | phiPhiDx(intval, w, s); | ||
| 325 | } | ||
| 326 | /** | ||
| 327 | * \brief Evaluate derivatives in x direction of local shape functions on the | ||
| 328 | * unit square at coordinate s. | ||
| 329 | */ | ||
| 330 | 91908 | Eigen::Matrix<Scalar, Eigen::Dynamic, 1> basisDx( | |
| 331 | const Eigen::Vector2d& s) const { | ||
| 332 | 91908 | Eigen::Matrix<Scalar, Eigen::Dynamic, 1> intval( | |
| 333 | 91908 | polynomial_degree_plus_one_squared); | |
| 334 |
1/2✓ Branch 1 taken 91908 times.
✗ Branch 2 not taken.
|
91908 | intval.setZero(); |
| 335 |
3/5✓ Branch 1 taken 91908 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 108 times.
✓ Branch 5 taken 91800 times.
✗ Branch 6 not taken.
|
91908 | phiPhiDx(&intval, 1., s); |
| 336 | 91908 | return intval; | |
| 337 | ✗ | } | |
| 338 | /** | ||
| 339 | * \brief Evaluate derivatives in y direction of local shape functions on the | ||
| 340 | * unit square at coordinate s, scale by w and add to intval. | ||
| 341 | */ | ||
| 342 | void addScaledBasisDy(Eigen::Matrix<Scalar, Eigen::Dynamic, 1>* intval, | ||
| 343 | typename LinearOperatorTraits<Derived>::Scalar w, | ||
| 344 | const Eigen::Vector2d& s) const { | ||
| 345 | phiPhiDy(intval, w, s); | ||
| 346 | } | ||
| 347 | /** | ||
| 348 | * \brief Evaluate derivatives in y direction of local shape functions on the | ||
| 349 | * unit square at coordinate s. | ||
| 350 | */ | ||
| 351 | 91908 | Eigen::Matrix<Scalar, Eigen::Dynamic, 1> basisDy( | |
| 352 | const Eigen::Vector2d& s) const { | ||
| 353 | 91908 | Eigen::Matrix<Scalar, Eigen::Dynamic, 1> intval( | |
| 354 | 91908 | polynomial_degree_plus_one_squared); | |
| 355 |
1/2✓ Branch 1 taken 91908 times.
✗ Branch 2 not taken.
|
91908 | intval.setZero(); |
| 356 |
3/5✓ Branch 1 taken 91908 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 108 times.
✓ Branch 5 taken 91800 times.
✗ Branch 6 not taken.
|
91908 | phiPhiDy(&intval, 1., s); |
| 357 | 91908 | return intval; | |
| 358 | ✗ | } | |
| 359 | /** | ||
| 360 | * \brief Evaluate local shape functions on the unit interval at coordinate s, | ||
| 361 | * scale by w and add to intval. | ||
| 362 | */ | ||
| 363 | 13320 | void addScaledBasis1D(Eigen::Matrix<Scalar, Eigen::Dynamic, 1>* intval, | |
| 364 | Scalar w, double s) const { | ||
| 365 | 13320 | phi(intval, w, s); | |
| 366 | 13320 | } | |
| 367 | /** | ||
| 368 | * \brief Evaluate local shape functions on the unit interval at coordinate s. | ||
| 369 | */ | ||
| 370 | Eigen::Matrix<Scalar, Eigen::Dynamic, 1> basis1D(double s) const { | ||
| 371 | Eigen::Matrix<Scalar, Eigen::Dynamic, 1> intval(polynomial_degree + 1); | ||
| 372 | intval.setZero(); | ||
| 373 | phi(&intval, 1., s); | ||
| 374 | return intval; | ||
| 375 | } | ||
| 376 | /** | ||
| 377 | * \brief Evaluate derivatives of local shape functions on the unit interval | ||
| 378 | * at coordinate s, scale by w and add to intval. | ||
| 379 | */ | ||
| 380 | 180 | void addScaledBasis1DDx(Eigen::Matrix<Scalar, Eigen::Dynamic, 1>* intval, | |
| 381 | Scalar w, double s) const { | ||
| 382 | 180 | phiDx(intval, w, s); | |
| 383 | 180 | } | |
| 384 | /** | ||
| 385 | * \brief Evaluate derivatives of local shape functions on the unit interval | ||
| 386 | * at coordinate s. | ||
| 387 | */ | ||
| 388 | Eigen::Matrix<Scalar, Eigen::Dynamic, 1> basis1DDx(double s) const { | ||
| 389 | Eigen::Matrix<Scalar, Eigen::Dynamic, 1> intval(polynomial_degree + 1); | ||
| 390 | intval.setZero(); | ||
| 391 | phiDx(&intval, 1., s); | ||
| 392 | return intval; | ||
| 393 | } | ||
| 394 | ////////////////////////////////////////////////////////////////////////////// | ||
| 395 | // member variables | ||
| 396 | ////////////////////////////////////////////////////////////////////////////// | ||
| 397 | private: | ||
| 398 | std::shared_ptr<ClusterTree> mesh_; | ||
| 399 | Basis::funptr_phi<Scalar> phi; | ||
| 400 | Basis::funptr_phidx<Scalar> phiDx; | ||
| 401 | Basis::funptr_phiphi<Scalar> phiPhi; | ||
| 402 | Basis::funptr_phiphidx<Scalar> phiPhiDx; | ||
| 403 | Basis::funptr_phiphidy<Scalar> phiPhiDy; | ||
| 404 | Basis::funptr_phitimesphi<Scalar> phiTimesPhi; | ||
| 405 | // Basis::funptr_vphiscalvphi<typename LinearOperatorTraits<Derived>::Scalar> | ||
| 406 | // vPhiScalVPhi; | ||
| 407 | Basis::funptr_divphitimesdivphi<Scalar> divPhiTimesDivPhi; | ||
| 408 | int polynomial_degree; | ||
| 409 | int polynomial_degree_plus_one_squared; | ||
| 410 | }; // namespace Bembel | ||
| 411 | } // namespace Bembel | ||
| 412 | #endif // BEMBEL_SRC_ANSATZSPACE_SUPERSPACE_HPP_ | ||
| 413 |