| 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_POTENTIAL_DISCRETEPOTENTIAL_HPP_ | ||
| 12 | #define BEMBEL_SRC_POTENTIAL_DISCRETEPOTENTIAL_HPP_ | ||
| 13 | |||
| 14 | namespace Bembel { | ||
| 15 | /** | ||
| 16 | * \ingroup Potential | ||
| 17 | * \brief DiscretePotential | ||
| 18 | * \todo Add a documentation | ||
| 19 | */ | ||
| 20 | template <typename Derived, typename LinOp> | ||
| 21 | class DiscretePotential { | ||
| 22 | public: | ||
| 23 | ////////////////////////////////////////////////////////////////////////////// | ||
| 24 | // constructors | ||
| 25 | ////////////////////////////////////////////////////////////////////////////// | ||
| 26 | DiscretePotential() {} | ||
| 27 |
3/5✓ Branch 2 taken 2 times.
✓ Branch 3 taken 4 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
|
6 | explicit DiscretePotential(const AnsatzSpace<LinOp> &ansatz_space) { |
| 28 |
1/2✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
|
6 | init_DiscretePotential(ansatz_space); |
| 29 | 6 | } | |
| 30 | ////////////////////////////////////////////////////////////////////////////// | ||
| 31 | // init_DiscretePotential | ||
| 32 | ////////////////////////////////////////////////////////////////////////////// | ||
| 33 | 6 | void init_DiscretePotential(const AnsatzSpace<LinOp> &ansatz_space) { | |
| 34 |
1/2✓ Branch 2 taken 6 times.
✗ Branch 3 not taken.
|
6 | ansatz_space_ = ansatz_space; |
| 35 |
1/2✓ Branch 2 taken 6 times.
✗ Branch 3 not taken.
|
6 | fun_ev_ = FunctionEvaluator<LinOp>(ansatz_space_); |
| 36 | /** | ||
| 37 | * \todo obtain this from ansatz space | ||
| 38 | */ | ||
| 39 | 6 | deg_ = ansatz_space_.get_polynomial_degree() + 1; | |
| 40 | 6 | return; | |
| 41 | } | ||
| 42 | ////////////////////////////////////////////////////////////////////////////// | ||
| 43 | // compute | ||
| 44 | ////////////////////////////////////////////////////////////////////////////// | ||
| 45 | Eigen::Matrix<typename PotentialReturnScalar< | ||
| 46 | typename LinearOperatorTraits<LinOp>::Scalar, | ||
| 47 | typename PotentialTraits<Derived>::Scalar>::Scalar, | ||
| 48 | Eigen::Dynamic, PotentialTraits<Derived>::OutputSpaceDimension> | ||
| 49 | 6 | evaluate(const Eigen::Matrix<double, Eigen::Dynamic, 3> &points) { | |
| 50 | 6 | auto FunctionSpaceVectorDimension = | |
| 51 | getFunctionSpaceVectorDimension<LinearOperatorTraits<LinOp>::Form>(); | ||
| 52 | 6 | auto OutputDimension = PotentialTraits<Derived>::OutputSpaceDimension; | |
| 53 | |||
| 54 |
1/2✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
|
6 | GaussSquare<Constants::maximum_quadrature_degree> GS; |
| 55 |
2/4✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 6 times.
✗ Branch 5 not taken.
|
6 | auto Q = GS[deg_]; |
| 56 | |||
| 57 | 6 | auto super_space = ansatz_space_.get_superspace(); | |
| 58 | |||
| 59 | 6 | const ElementTree &element_tree = super_space.get_mesh().get_element_tree(); | |
| 60 | 6 | auto number_of_elements = element_tree.get_number_of_elements(); | |
| 61 | |||
| 62 | 6 | auto polynomial_degree = super_space.get_polynomial_degree(); | |
| 63 | 6 | auto polynomial_degree_plus_one_squared = | |
| 64 | 6 | (polynomial_degree + 1) * (polynomial_degree + 1); | |
| 65 | |||
| 66 | Eigen::Matrix<typename PotentialReturnScalar< | ||
| 67 | typename LinearOperatorTraits<LinOp>::Scalar, | ||
| 68 | typename PotentialTraits<Derived>::Scalar>::Scalar, | ||
| 69 | Eigen::Dynamic, | ||
| 70 | PotentialTraits<Derived>::OutputSpaceDimension> | ||
| 71 |
1/2✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
|
6 | potential; |
| 72 |
1/2✓ Branch 2 taken 6 times.
✗ Branch 3 not taken.
|
6 | potential.resize(points.rows(), |
| 73 | PotentialTraits<Derived>::OutputSpaceDimension); | ||
| 74 |
1/2✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
|
6 | potential.setZero(); |
| 75 | |||
| 76 | #pragma omp parallel | ||
| 77 | { | ||
| 78 | Eigen::Matrix<typename PotentialReturnScalar< | ||
| 79 | typename LinearOperatorTraits<LinOp>::Scalar, | ||
| 80 | typename PotentialTraits<Derived>::Scalar>::Scalar, | ||
| 81 | Eigen::Dynamic, | ||
| 82 | PotentialTraits<Derived>::OutputSpaceDimension> | ||
| 83 |
1/2✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
|
6 | my_potential; |
| 84 |
1/2✓ Branch 2 taken 6 times.
✗ Branch 3 not taken.
|
6 | my_potential.resize(points.rows(), |
| 85 | PotentialTraits<Derived>::OutputSpaceDimension); | ||
| 86 |
1/2✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
|
6 | my_potential.setZero(); |
| 87 | 6 | for (auto element = element_tree.cpbegin(); | |
| 88 |
2/2✓ Branch 3 taken 414 times.
✓ Branch 4 taken 6 times.
|
420 | element != element_tree.cpend(); ++element) { |
| 89 | #pragma omp single nowait | ||
| 90 | { | ||
| 91 |
1/2✓ Branch 1 taken 414 times.
✗ Branch 2 not taken.
|
414 | SurfacePoint qp; |
| 92 |
2/2✓ Branch 1 taken 2166 times.
✓ Branch 2 taken 414 times.
|
2580 | for (auto j = 0; j < Q.w_.size(); ++j) { |
| 93 |
2/4✓ Branch 1 taken 2166 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2166 times.
✗ Branch 5 not taken.
|
4332 | super_space.map2surface( |
| 94 | 2166 | *element, Q.xi_.col(j), | |
| 95 |
2/4✓ Branch 5 taken 2166 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 2166 times.
✗ Branch 9 not taken.
|
2166 | element->get_h() * element->get_h() * Q.w_(j), &qp); |
| 96 |
2/2✓ Branch 1 taken 994200 times.
✓ Branch 2 taken 2166 times.
|
996366 | for (auto i = 0; i < points.rows(); ++i) { |
| 97 |
3/6✓ Branch 1 taken 994200 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 994200 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 994200 times.
✗ Branch 8 not taken.
|
1988400 | my_potential.row(i) += pot_.evaluateIntegrand_impl( |
| 98 |
2/4✓ Branch 1 taken 994200 times.
✗ Branch 2 not taken.
✓ Branch 5 taken 994200 times.
✗ Branch 6 not taken.
|
2982600 | fun_ev_, *element, points.row(i), qp); |
| 99 | } | ||
| 100 | } | ||
| 101 | } | ||
| 102 | } | ||
| 103 | #pragma omp critical | ||
| 104 |
1/2✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
|
6 | potential += my_potential; |
| 105 | 6 | } | |
| 106 | 12 | return potential; | |
| 107 | 6 | } | |
| 108 | ////////////////////////////////////////////////////////////////////////////// | ||
| 109 | // setter | ||
| 110 | ////////////////////////////////////////////////////////////////////////////// | ||
| 111 | 6 | void set_cauchy_data( | |
| 112 | const Eigen::Matrix<typename LinearOperatorTraits<LinOp>::Scalar, | ||
| 113 | Eigen::Dynamic, 1> &cauchy_data) { | ||
| 114 |
1/2✓ Branch 2 taken 6 times.
✗ Branch 3 not taken.
|
6 | fun_ev_.set_function(cauchy_data); |
| 115 | 6 | } | |
| 116 | void set_degree(const int °) { deg_ = deg; } | ||
| 117 | ////////////////////////////////////////////////////////////////////////////// | ||
| 118 | // getter | ||
| 119 | ////////////////////////////////////////////////////////////////////////////// | ||
| 120 | 3 | Derived &get_potential() { return pot_; } | |
| 121 | ////////////////////////////////////////////////////////////////////////////// | ||
| 122 | // private member variables | ||
| 123 | ////////////////////////////////////////////////////////////////////////////// | ||
| 124 | private: | ||
| 125 | int deg_; | ||
| 126 | Derived pot_; | ||
| 127 | AnsatzSpace<LinOp> ansatz_space_; | ||
| 128 | FunctionEvaluator<LinOp> fun_ev_; | ||
| 129 | }; // namespace Bembel | ||
| 130 | |||
| 131 | } // namespace Bembel | ||
| 132 | #endif // BEMBEL_SRC_POTENTIAL_DISCRETEPOTENTIAL_HPP_ | ||
| 133 |