Homec4science
Diffusion expolit (master) Browse

expolit/integration.hh
046a43293eecmaster

Last Change
Raw File

integration.hh
View Options

/**
* @file
*
* @author Lucas Frérot <lucas.frerot@epfl.ch>
*
* @section LICENSE
*
* Copyright (©) 2018 EPFL (Ecole Polytechnique Fédérale de
* Lausanne) Laboratory (LSMS - Laboratoire de Simulation en Mécanique des
* Solides)
*
* Expolit is free software: you can redistribute it and/or modify it under the
* terms of the GNU Lesser General Public License as published by the Free
* Software Foundation, either version 3 of the License, or (at your option) any
* later version.
*
* Expolit 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 Lesser General Public License for
* more details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with Tamaas. If not, see <http://www.gnu.org/licenses/>.
*
*/
#ifndef INTEGRATION_HH
#define INTEGRATION_HH
#include "differentiation.hh"
#include "exponential.hh"
#include "polynomial.hh"
#include "types.hh"
namespace expolit {
template <typename T, UInt order, typename CT>
constexpr auto integrate(const Polynomial<T, order, CT>& p) {
Polynomial<T, order + 1, CT> ip;
for (UInt i = 1; i < ip.terms; ++i)
ip.coeffs[i] = p.coeffs[i - 1] / static_cast<Real>(i);
return ip;
}
template <typename T, typename CT, typename Derived>
constexpr auto integrate(const Product<Constant<T, CT>, Derived>& p) {
return p.operands.first * integrate(p.operands.second);
}
template <typename T, typename CT, typename Derived>
constexpr auto integrate(const Product<Derived, Constant<T, CT>>& p) {
return p.operands.second * integrate(p.operands.first);
}
template <typename T, typename CT>
constexpr auto integrate(const Exponential<Polynomial<T, 1, CT>>& e) {
return Constant<T, CT>({1 / e.expression.coeffs.back()}) * e;
}
template <typename Derived, typename T2, typename CT2, UInt order>
constexpr auto
integrate(const Product<Derived, Polynomial<T2, order, CT2>>& p) {
auto u = integrate(p.operands.first);
return u * p.operands.second -
integrate(u * differentiate(p.operands.second));
}
template <typename Derived, typename T2, typename CT2, UInt order>
constexpr auto
integrate(const Product<Polynomial<T2, order, CT2>, Derived>& p) {
return integrate(p.commute());
}
template <typename T, typename Derived>
constexpr auto
definite_integral(const std::pair<T, T>& bounds, const Expression<Derived>& exp) {
auto integ = integrate(static_cast<const Derived&>(exp));
return integ(bounds.second) - integ(bounds.first);
}
} // namespace expolit
#endif
c4science · Help