laplace_2d_singular_collocation_integrals.hpp

Go to the documentation of this file.

// This file is a part of NiHu, a C++ BEM template library.
//
// Copyright (C) 2012-2014  Peter Fiala <fiala@hit.bme.hu>
// Copyright (C) 2012-2014  Peter Rucz <rucz@hit.bme.hu>
//
// This program is free software: you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
// GNU General Public License for more details.
//
// You should have received a copy of the GNU General Public License
// along with this program.  If not, see <http://www.gnu.org/licenses/>.
 
 
#ifndef LAPLACE_2D_SINGULAR_COLLOCATION_INTEGRALS_HPP_INCLUDED
#define LAPLACE_2D_SINGULAR_COLLOCATION_INTEGRALS_HPP_INCLUDED
 
#include "../core/match_types.hpp"
#include "../core/singular_integral_shortcut.hpp"
#include "../util/math_functions.hpp"
#include "integral_shortcut_concepts.hpp"
#include "lib_element.hpp"
#include "laplace_kernel.hpp"
#include "normal_derivative_singular_integrals.hpp"
#include "quadrature_store_helper.hpp"
 
#include <boost/math/constants/constants.hpp>
 
#define NIHU_PRINT_DEBUG
 
#ifdef NIHU_PRINT_DEBUG
#include <iostream>
#endif
 
namespace NiHu
{
 
template <class TestField, class TrialField, size_t order>
class laplace_2d_SLP_collocation_curved
{
    typedef TestField test_field_t;
    typedef TrialField trial_field_t;
 
    typedef typename test_field_t::nset_t test_shape_t;
    typedef typename trial_field_t::nset_t trial_shape_t;
 
    static size_t const nTest = test_shape_t::num_nodes;
    static size_t const nTrial = trial_shape_t::num_nodes;
 
    typedef Eigen::Matrix<double, nTest, nTrial> result_t;
    typedef Eigen::Matrix<double, 1, nTrial> result_row_t;
 
    typedef typename trial_field_t::elem_t elem_t;
    typedef typename elem_t::domain_t domain_t;
 
    typedef typename domain_t::xi_t xi_t;
    typedef typename elem_t::x_t x_t;
 
    typedef regular_quad_store<domain_t, order> reg_quadrature_t;
    static size_t const inner_order =
        trial_field_t::nset_t::polynomial_order +   // Ny (exact)
        elem_t::lset_t::jacobian_order * 2;         // Jy (possible overkill for quadratic, exact for linear)
    typedef log_quad_store<inner_order> log_quadrature_t;
 
public:
 
    static result_row_t eval_row(xi_t const &xi0, elem_t const &elem)
    {
        using boost::math::double_constants::two_pi;
 
        result_row_t result_row;
        result_row.setZero();
 
        // the singular location
        x_t x = elem.get_x(xi0);
 
        // traverse integration domains
        // dom = 0: (-1 .. xi0)
        // dom = 1: (xi0 .. +1)
        for (int dom = 0; dom < 2; ++dom) {
            xi_t xi_lim = domain_t::get_corner(dom);
            double jac = std::abs(xi_lim(0) - xi0(0));
 
            // the collocation point is in the corner
            if (jac < 1e-5)
                continue;
 
            // traverse regular quadrature points
            for (auto const &q : reg_quadrature_t::quadrature) {
                // transform quadrature to [0 1];
                double v = (q.get_xi()(0) + 1.) / 2.;
                double wv = q.get_w() / 2;
 
                xi_t eta = xi0 + v * (xi_lim - xi0);
                double r = (elem.get_x(eta) - x).norm();
                double Jeta = elem.get_dx(eta).norm();
 
                // evaluate Green's function
                double Greg = -std::log(r / v) / two_pi;
 
                // integrate difference numerically
                result_row += Greg * Jeta * jac * wv
                    * trial_shape_t::template eval_shape<0>(eta);
            } // loop over regular quadrature
 
            // traverse log quadrature points
            for (auto const &q : log_quadrature_t::quadrature) {
                double v = q.get_xi()(0);
                double wv = q.get_w();
 
                xi_t eta = xi0 + v * (xi_lim - xi0);
                double Jeta = elem.get_dx(eta).norm();
 
                // evaluate Green's function
                double Glog = 1.0 / two_pi;
 
                // integrate difference numerically
                result_row += Glog * Jeta * jac * wv
                    * trial_shape_t::template eval_shape<0>(eta);
            } // loop over log quadrature
        } // loop over subdomains
        return result_row;
    }
 
    static result_t eval(elem_t const &elem)
    {
#ifdef NIHU_PRINT_DEBUG
        static bool printed = false;
        if (!printed) {
            std::cout << "Using laplace_2d_SLP_collocation_curved, inner order: " << inner_order << std::endl;
            printed = true;
        }
#endif
 
        result_t result;
 
        // traverse collocation points (rows)
        for (size_t i = 0; i < nTest; ++i) {
            xi_t const &xi0 = test_shape_t::corner_at(i);
            result.row(i) = eval_row(xi0, elem);
        } // loop over collocation points
 
        return result;
    } // function eval
}; // class laplace_2d_SLP_collocation_curved
 
 
template <class TestField, class TrialField>
class laplace_2d_SLP_collocation_straight
{
    typedef TestField test_field_t;
    typedef TrialField trial_field_t;
 
    typedef typename test_field_t::nset_t test_shape_set_t;
    typedef typename trial_field_t::nset_t trial_shape_set_t;
 
    typedef line_1_elem elem_t;
 
    typedef typename test_shape_set_t::xi_t xi_t;
 
    static size_t const rows = test_shape_set_t::num_nodes;
    static size_t const cols = trial_shape_set_t::num_nodes;
 
    typedef Eigen::Matrix<double, rows, cols> result_t;
    typedef Eigen::Matrix<double, 1, cols> result_row_t;
 
public:
    static result_row_t eval_row(xi_t const &xi0, elem_t const &elem)
    {
        using boost::math::double_constants::two_pi;
 
        // get Jacobian
        double jac = elem.get_normal().norm();
 
        result_row_t result_row;
        result_row.setZero();
 
        for (size_t s = 0; s < 2; s++)  // s == 0 : left side, s == 1 : right side
        {
            double sign = (s == 0 ? 1. : -1.);
            double rho = 1. + sign * xi0(0);    // length of the actual side in intrinsic
            if (rho < 1e-5)
                continue;
 
            double d = jac * rho;               // physical length of the actual side
            double logd = std::log(d);
 
            auto N0 = trial_shape_set_t::template eval_shape<0>(xi0);
            result_row += N0 * d * (1. - logd);
 
            if (trial_shape_set_t::polynomial_order >= 1) {
                auto N1 = (trial_shape_set_t::template eval_shape<1>(xi0) / jac).eval();
                result_row += -sign * N1 * d * d * (1. / 2. - logd) / 2.;
            }
 
            if (trial_shape_set_t::polynomial_order >= 2) {
                auto N2 = (trial_shape_set_t::template eval_shape<2>(xi0) / jac / jac / 2.).eval();
                result_row += N2 * d * d * d * (1. / 3. - logd) / 3.;
            }
        } // loop over two sides of the integration domain
 
        return result_row / two_pi;
    } // function eval
 
 
    static result_t eval(elem_t const &elem)
    {
#ifdef NIHU_PRINT_DEBUG
        static bool printed = false;
        if (!printed) {
            std::cout << "Using laplace_2d_SLP_collocation_straight" << std::endl;
            printed = true;
        }
#endif
 
        result_t result;
        // traverse collocational points
        for (size_t i = 0; i < rows; ++i) {
            xi_t const &xi0 = test_shape_set_t::corner_at(i);
            result.row(i) = eval_row(xi0, elem);
        } // loop over collocation points
        return result;
    } // function eval
}; // class laplace_2d_SLP_collocation_straight
 
 
template <class TestField, class TrialField, size_t order>
class laplace_2d_DLP_collocation_curved
{
    typedef TestField test_field_t;
    typedef TrialField trial_field_t;
 
    typedef typename test_field_t::nset_t test_shape_t;
    typedef typename trial_field_t::nset_t trial_shape_t;
 
    static size_t const nTest = test_shape_t::num_nodes;
    static size_t const nTrial = trial_shape_t::num_nodes;
 
    typedef Eigen::Matrix<double, nTest, nTrial> result_t;
 
    typedef typename trial_field_t::elem_t elem_t;
    typedef typename elem_t::domain_t domain_t;
 
    typedef typename domain_t::xi_t xi_t;
    typedef typename elem_t::x_t x_t;
 
    typedef regular_quad_store<domain_t, order> quadrature_t;
 
    typedef laplace_2d_DLP_kernel kernel_t;
 
public:
    static result_t eval(elem_t const &elem)
    {
#ifdef NIHU_PRINT_DEBUG
        static bool printed = false;
        if (!printed) {
            std::cout << "Using laplace_2d_DLP_collocation_curved" << std::endl;
            printed = true;
        }
#endif
 
        kernel_t kernel;
 
        result_t result = result_t::Zero();
 
        // traverse collocation points
        for (size_t i = 0; i < nTest; ++i) {
            xi_t const &xi0 = test_shape_t::corner_at(i);
 
            // singular point and jacobians
            x_t x = elem.get_x(xi0);
 
            // traverse integration domains
            for (int dom = 0; dom < 2; ++dom) {
                double low = xi0(0), high = domain_t::get_corner(dom)(0);
                xi_t center;
                center << (low + high) / 2.;
 
                // traverse quadrature points
                for (auto const &q : quadrature_t::quadrature) {
                    // transform quadrature to [low high];
                    xi_t xi = q.get_xi() * ((high - low) / 2.) + center;
                    double w = q.get_w() * (std::abs(high - low) / 2.);
 
                    // get trial location, jacobian and normal
                    x_t y = elem.get_x(xi);
                    x_t Jy = elem.get_normal(xi);
 
                    // evaluate Green's function
                    double G = kernel(x, y, Jy);
                    // evaluate shape function
                    auto Ny = trial_shape_t::template eval_shape<0>(xi);
                    // integrate kernel
                    result.row(i) += G * Ny * w;
                } // quadrature points
            } // loop over subdomains
        } // loop over collocation points
 
        return result;
    } // function eval
}; // laplace_2d_DLP_collocation_curved
 
 
template <class TestField, class TrialField, size_t order>
class laplace_2d_DLPt_collocation_curved
{
    typedef TestField test_field_t;
    typedef TrialField trial_field_t;
 
    typedef typename test_field_t::nset_t test_shape_t;
    typedef typename trial_field_t::nset_t trial_shape_t;
 
    static size_t const nTest = test_shape_t::num_nodes;
    static size_t const nTrial = trial_shape_t::num_nodes;
 
    typedef Eigen::Matrix<double, nTest, nTrial> result_t;
 
    typedef typename trial_field_t::elem_t elem_t;
    typedef typename elem_t::domain_t domain_t;
 
    typedef typename domain_t::xi_t xi_t;
    typedef typename elem_t::x_t x_t;
 
    typedef regular_quad_store<domain_t, order> quadrature_t;
 
    typedef laplace_2d_DLPt_kernel kernel_t;
 
public:
    static result_t eval(elem_t const &elem)
    {
#ifdef NIHU_PRINT_DEBUG
        static bool printed = false;
        if (!printed) {
            std::cout << "Using laplace_2d_DLPt_collocation_curved" << std::endl;
            printed = true;
        }
#endif
 
        kernel_t kernel;
 
        result_t result = result_t::Zero();
 
        // traverse collocation points
        for (size_t i = 0; i < nTest; ++i) {
            xi_t const &xi0 = test_shape_t::corner_at(i);
 
            // singular point and jacobians
            x_t x = elem.get_x(xi0);
            x_t nx = elem.get_normal(xi0).normalized();
 
            // traverse integration domains
            for (int dom = 0; dom < 2; ++dom) {
                double low = xi0(0), high = domain_t::get_corner(dom)(0);
                xi_t center;
                center << (low + high) / 2.;
 
                // traverse quadrature points
                for (auto const &q : quadrature_t::quadrature) {
                    // transform quadrature to [low high];
                    xi_t xi = q.get_xi() * ((high - low) / 2.) + center;
                    double w = q.get_w() * (std::abs(high - low) / 2.);
 
                    // get trial location, jacobian and normal
                    x_t y = elem.get_x(xi);
                    double Jy = elem.get_dx(xi).norm();
 
                    // evaluate Green's function
                    double G = kernel(x, y, nx);
                    // evaluate shape function
                    auto Ny = trial_shape_t::template eval_shape<0>(xi);
                    // integrate kernel
                    result.row(i) += G * Jy * Ny * w;
                } // quadrature points
            } // loop over subdomains
        } // loop over collocation points
 
        return result;
    } // function eval
}; // laplace_2d_DLPt_collocation_curved
 
 
template <class TestField, class TrialField, size_t order>
class laplace_2d_HSP_collocation_curved
{
    typedef TestField test_field_t;
    typedef TrialField trial_field_t;
 
    typedef typename test_field_t::nset_t test_shape_t;
    typedef typename trial_field_t::nset_t trial_shape_t;
 
    static size_t const nTest = test_shape_t::num_nodes;
    static size_t const nTrial = trial_shape_t::num_nodes;
 
    typedef Eigen::Matrix<double, nTest, nTrial> result_t;
 
    typedef typename trial_field_t::elem_t elem_t;
    typedef typename elem_t::domain_t domain_t;
 
    typedef typename domain_t::xi_t xi_t;
    typedef typename elem_t::x_t x_t;
 
    typedef regular_quad_store<domain_t, order> quadrature_t;
 
    typedef laplace_2d_HSP_kernel kernel_t;
 
public:
    static result_t eval(elem_t const &elem)
    {
#ifdef NIHU_PRINT_DEBUG
        static bool printed = false;
        if (!printed) {
            std::cout << "Using laplace_2d_HSP_collocation_curved" << std::endl;
            printed = true;
        }
#endif
 
        using namespace boost::math::double_constants;
 
        kernel_t kernel;
 
        result_t result = result_t::Zero();
 
        xi_t const &a = domain_t::get_corner(0);
        xi_t const &b = domain_t::get_corner(1);
 
        // traverse collocation points
        for (size_t i = 0; i < nTest; ++i) {
            xi_t const &xi0 = test_shape_t::corner_at(i);
 
            // trial shape function and its derivative at singular point
            auto N0 = trial_shape_t::template eval_shape<0>(xi0);
            auto N1 = trial_shape_t::template eval_shape<1>(xi0);
 
            // singular point and normal at singular point
            x_t x = elem.get_x(xi0);
            x_t Jxvec = elem.get_normal(xi0);
            x_t nx = Jxvec.normalized();
            double jac0 = Jxvec.norm();
            double twopiJ0 = two_pi * jac0;
 
            for (int dom = 0; dom < 2; ++dom) {
                double high = domain_t::get_corner(dom)(0);
                double low = xi0(0);
                xi_t center;
                center << (low + high) / 2.;
 
                // traverse quadrature points
                for (auto const &q : quadrature_t::quadrature) {
                    // transform quadrature to [xi0 b];
                    xi_t xi = q.get_xi() * ((high - low) / 2.) + center;
                    double w = q.get_w() * (std::abs(high - low) / 2.);
                    double rho = xi(0) - xi0(0);
 
                    // get trial location, jacobian and normal
                    x_t y = elem.get_x(xi);
                    x_t Jyvec = elem.get_normal(xi);
                    x_t ny = Jyvec.normalized();
                    double jac = Jyvec.norm();
 
                    // evaluate Green's function
                    double G = kernel(x, y, nx, ny);
                    // multiply by Shape function and Jacobian
                    auto N = trial_shape_t::template eval_shape<0>(xi);
                    auto F = G * N * jac;
                    // evaluate singular part
                    auto F0 = (N0 / rho + N1) / rho / twopiJ0;
 
                    // integrate difference numerically
                    result.row(i) += (F - F0) * w;
                } // loop over quadrature points
            } // loop over domains
 
            // add analytic integral of singular part
            result.row(i) += ((1. / (a(0) - xi0(0)) - 1. / (b(0) - xi0(0))) * N0 +
                std::log(std::abs((b(0) - xi0(0)) / (a(0) - xi0(0)))) * N1) / twopiJ0;
        } // collocation points
        return result;
    } // function  eval
}; // class laplace_2d_HSP_collocation_curved
 
 
template <class TestField, class TrialField>
class laplace_2d_HSP_collocation_straight
{
    typedef TestField test_field_t;
    typedef TrialField trial_field_t;
 
    typedef typename test_field_t::nset_t test_shape_set_t;
    typedef typename trial_field_t::nset_t trial_shape_set_t;
 
    typedef typename test_field_t::elem_t elem_t;
 
    typedef typename test_shape_set_t::xi_t xi_t;
 
    static size_t const nTest = test_shape_set_t::num_nodes;
    static size_t const nTrial = trial_shape_set_t::num_nodes;
 
    typedef Eigen::Matrix<double, nTest, nTrial> result_t;
 
public:
    static result_t eval(elem_t const &elem)
    {
        using namespace boost::math::double_constants;
 
#ifdef NIHU_PRINT_DEBUG
        static bool printed = false;
        if (!printed) {
            std::cout << "Using laplace_2d_HSP_collocation_straight" << std::endl;
            printed = true;
        }
#endif
 
        double jac = elem.get_normal().norm();
 
        result_t result = result_t::Zero();
 
        // integration limits in intrinsic coordinates
        double const a = -1.;
        double const b = +1.;
 
        for (size_t i = 0; i < nTest; ++i) {
            xi_t xi0 = test_shape_set_t::corner_at(i);
            auto N = trial_shape_set_t::template eval_shape<0>(xi0);
            result.row(i) = N * (-1. / (b - xi0(0)) - -1. / (a - xi0(0)));
 
            if (trial_shape_set_t::polynomial_order >= 1) {
                auto dN = trial_shape_set_t::template eval_shape<1>(xi0);
                result.row(i) += dN * std::log(std::abs((b - xi0(0)) / (a - xi0(0))));
            }
 
            if (trial_shape_set_t::polynomial_order >= 2) {
                auto ddN = trial_shape_set_t::template eval_shape<2>(xi0);
                result.row(i) += ddN / 2. * (b - a);
            }
        } // collocation points
 
        return result / (two_pi * jac);
    } // function eval
}; // class laplace_2d_HSP_collocation_straight
 
template <class TestField, class TrialField>
class singular_integral_shortcut<laplace_2d_SLP_kernel, TestField, TrialField, match::match_1d_type,
    typename std::enable_if<
        is_collocational<TestField, TrialField>::value &&
        is_straight_line<TrialField>::value
    >::type
>
{
public:
    template <class result_t>
    static result_t &eval(
        result_t &result,
        kernel_base<laplace_2d_SLP_kernel> const &,
        field_base<TestField> const &,
        field_base<TrialField> const &trial_field,
        element_match const &)
    {
        result = laplace_2d_SLP_collocation_straight<
            TestField, TrialField
        >::eval(trial_field.get_elem());
        return result;
    }
};
 
 
template <class TestField, class TrialField>
class singular_integral_shortcut<laplace_2d_SLP_kernel, TestField, TrialField, match::match_0d_type,
    typename std::enable_if<
        is_collocational<TestField, TrialField>::value &&
        is_straight_line<TrialField>::value
    >::type
>
{
    typedef typename TrialField::elem_t trial_elem_t;
    typedef typename trial_elem_t::domain_t domain_t;
    typedef typename domain_t::xi_t xi_t;
    static unsigned const order = 8;
    typedef regular_quad_store<domain_t, order> reg_quadrature_t;
    typedef typename TestField::nset_t test_nset_t;
    typedef typename TrialField::nset_t trial_nset_t;
    static size_t const nTest = test_nset_t::num_nodes;
 
public:
    template <class result_t>
    static result_t &eval(
        result_t &result,
        kernel_base<laplace_2d_SLP_kernel> const &kernel,
        field_base<TestField> const &test_field,
        field_base<TrialField> const &trial_field,
        element_match const &match)
    {
        unsigned test_idx = match.get_overlap().get_ind1();
        auto test_singular_xi = TestField::elem_t::lset_t::corner_at(test_idx);
 
        unsigned trial_idx = match.get_overlap().get_ind2();
        auto trial_singular_xi = TrialField::elem_t::lset_t::corner_at(trial_idx);
 
        auto const &test_elem = test_field.get_elem();
        auto const &trial_elem = trial_field.get_elem();
 
        for (size_t i = 0; i < nTest; ++i)
        {
            auto xi0 = TestField::nset_t::corner_at(i);
            if (xi0 == test_singular_xi)
            {
                // singular integral
                result.row(i) = laplace_2d_SLP_collocation_straight<TestField, TrialField>::eval_row(trial_singular_xi, trial_field.get_elem());
            }
            else
            {
                result.row(i).setZero();
                for (auto const &q : reg_quadrature_t::quadrature) {
                    xi_t xi = q.get_xi();
                    double w = q.get_w();
                    double jac = trial_elem.get_normal(xi).norm();
                    double K = kernel(test_elem.get_x(xi0), trial_elem.get_x(xi));
                    result.row(i) += (w * jac * K) * trial_nset_t::template eval_shape<0>(xi);
                }
            }
        }
 
        return result;
    }
};
 
 
 
template <class TestField, class TrialField>
class singular_integral_shortcut<laplace_2d_DLP_kernel, TestField, TrialField, match::match_0d_type,
    typename std::enable_if<
        is_collocational<TestField, TrialField>::value &&
        is_straight_line<TrialField>::value
    >::type
>
{
    typedef typename TrialField::elem_t trial_elem_t;
    typedef typename trial_elem_t::x_t x_t;
    typedef typename trial_elem_t::domain_t domain_t;
    typedef typename domain_t::xi_t xi_t;
    static unsigned const order = 8;
    typedef regular_quad_store<domain_t, order> reg_quadrature_t;
    typedef typename TestField::nset_t test_nset_t;
    typedef typename TrialField::nset_t trial_nset_t;
    static size_t const nTest = test_nset_t::num_nodes;
 
public:
    template <class result_t>
    static result_t &eval(
        result_t &result,
        kernel_base<laplace_2d_DLP_kernel> const &kernel,
        field_base<TestField> const &test_field,
        field_base<TrialField> const &trial_field,
        element_match const &match)
    {
        unsigned test_idx = match.get_overlap().get_ind1();
        auto test_singular_xi = TestField::elem_t::lset_t::corner_at(test_idx);
 
        // unsigned trial_idx = match.get_overlap().get_ind2();
        // auto trial_singular_xi = TrialField::elem_t::lset_t::corner_at(trial_idx);
 
        auto const &test_elem = test_field.get_elem();
        auto const &trial_elem = trial_field.get_elem();
 
        for (size_t i = 0; i < nTest; ++i)
        {
            auto xi0 = TestField::nset_t::corner_at(i);
            if (xi0 == test_singular_xi)
            {
                // singular integral
                result.row(i).setZero();
            }
            else
            {
                result.row(i).setZero();
                for (auto const &q : reg_quadrature_t::quadrature) {
                    xi_t xi = q.get_xi();
                    double w = q.get_w();
                    x_t normal = trial_elem.get_normal(xi);
                    x_t unit_normal = normal.normalized();
                    double jac = normal.norm();
                    double K = kernel.derived()(test_elem.get_x(xi0), trial_elem.get_x(xi), unit_normal);
                    result.row(i) += (w * jac * K) * trial_nset_t::template eval_shape<0>(xi);
                }
            }
        }
 
        return result;
    }
};
 
 
 
 
template <class TestField, class TrialField>
// requires Collocational<TestField, TrialField> && NotStraightLine<TrialField>
class singular_integral_shortcut<laplace_2d_SLP_kernel, TestField, TrialField, match::match_1d_type,
    typename std::enable_if<
        is_collocational<TestField, TrialField>::value &&
        !is_straight_line<TrialField>::value
    >::type
>
{
public:
    template <class result_t>
    static result_t &eval(
        result_t &result,
        kernel_base<laplace_2d_SLP_kernel> const &,
        field_base<TestField> const &,
        field_base<TrialField> const &trial_field,
        element_match const &)
    {
        result = laplace_2d_SLP_collocation_curved<
            TestField, TrialField, 8
        >::eval(trial_field.get_elem());
        return result;
    }
};
 
 
template <class TestField, class TrialField>
// requires Collocational<TestField, TrialField> && NotStraightLine<TrialField>
class singular_integral_shortcut<laplace_2d_DLP_kernel, TestField, TrialField, match::match_1d_type,
    typename std::enable_if<
        is_collocational<TestField, TrialField>::value &&
        !is_straight_line<TrialField>::value
    >::type
>
{
public:
    template <class result_t>
    static result_t &eval(
        result_t &result,
        kernel_base<laplace_2d_DLP_kernel> const &,
        field_base<TestField> const &,
        field_base<TrialField> const &trial_field,
        element_match const &)
    {
        result = laplace_2d_DLP_collocation_curved<TestField, TrialField, 10>::eval(
            trial_field.get_elem());
        return result;
    }
};
 
 
template <class TestField, class TrialField>
// requires Collocational<TestField, TrialField> && NotStraightLine<TrialField>
class singular_integral_shortcut<laplace_2d_DLPt_kernel, TestField, TrialField, match::match_1d_type,
    typename std::enable_if<
        is_collocational<TestField, TrialField>::value &&
        !is_straight_line<TrialField>::value
    >::type
>
{
public:
    template <class result_t>
    static result_t &eval(
        result_t &result,
        kernel_base<laplace_2d_DLPt_kernel> const &,
        field_base<TestField> const &,
        field_base<TrialField> const &trial_field,
        element_match const &)
    {
        result = laplace_2d_DLPt_collocation_curved<TestField, TrialField, 10>::eval(
            trial_field.get_elem());
        return result;
    }
};
 
 
template <class TestField, class TrialField>
// requires Collocational<TestField, TrialField> && NotStraightLine<TrialField>
class singular_integral_shortcut<laplace_2d_HSP_kernel, TestField, TrialField, match::match_1d_type,
    typename std::enable_if<
        is_collocational<TestField, TrialField>::value &&
        !is_straight_line<TrialField>::value
    >::type
>
{
public:
    template <class result_t>
    static result_t &eval(
        result_t &result,
        kernel_base<laplace_2d_HSP_kernel> const &,
        field_base<TestField> const &,
        field_base<TrialField> const &trial_field,
        element_match const &)
    {
        result = laplace_2d_HSP_collocation_curved<
            TestField, TrialField, 10
        >::eval(trial_field.get_elem());
        return result;
    }
};
 
 
template <class TestField, class TrialField>
// requires Collocational<TestField, TrialField> && StraightLine<TrialField>
class singular_integral_shortcut<laplace_2d_HSP_kernel, TestField, TrialField, match::match_1d_type,
    typename std::enable_if<
        is_collocational<TestField, TrialField>::value &&
        is_straight_line<TrialField>::value
    >::type
>
{
public:
    template <class result_t>
    static result_t &eval(
        result_t &result,
        kernel_base<laplace_2d_HSP_kernel> const &,
        field_base<TestField> const &,
        field_base<TrialField> const &trial_field,
        element_match const &)
    {
        result = laplace_2d_HSP_collocation_straight<TestField, TrialField>::eval(
            trial_field.get_elem());
        return result;
    }
};
 
} // namespace NiHu
 
#endif