Homec4science
Diffusion tamaas (features/declarative) Browse

tamaas/src/surface/isopowerlaw.cpp
3927fa639a53features/declarative

Last Change
Raw File

isopowerlaw.cpp
View Options

/**
* @file
*
* @author Lucas Frérot <lucas.frerot@epfl.ch>
*
* @section LICENSE
*
* Copyright (©) 2016 EPFL (Ecole Polytechnique Fédérale de
* Lausanne) Laboratory (LSMS - Laboratoire de Simulation en Mécanique des
* Solides)
*
* Tamaas 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.
*
* Tamaas 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/>.
*
*/
/* -------------------------------------------------------------------------- */
#include "isopowerlaw.hh"
#include "fftransform.hh"
#include <map>
/* -------------------------------------------------------------------------- */
__BEGIN_TAMAAS__
template <UInt dim>
void Isopowerlaw<dim>::computeFilter(
GridHermitian<Real, dim>& filter_coefficients) const {
auto wavevectors = FFTransform<Real, dim>::template computeFrequencies<true>(
filter_coefficients.sizes());
auto end = wavevectors.end(dim);
#pragma omp parallel for
for (auto it = wavevectors.begin(dim); it < end; ++it) {
VectorProxy<Real> q(&(*it), dim);
UInt i = wavevectors.getNbPoints() - (end - it);
filter_coefficients(i) = (*this)(q);
}
}
/* -------------------------------------------------------------------------- */
template <UInt dim>
Real Isopowerlaw<dim>::rmsHeights() const {
return std::sqrt(M_PI * ((hurst + 1) / hurst * q1 * q1 -
1. / hurst * std::pow(q1, 2 * (hurst + 1)) *
std::pow(q2, -2 * hurst) -
q0 * q0));
}
/* -------------------------------------------------------------------------- */
/**
* Analytical moments, cf. Yastrebov et al. (2015)
* "From infinitesimal to full contact between rough surfaces: Evolution
* of the contact area", appendix A
*/
template <>
std::vector<Real> Isopowerlaw<2>::moments() const {
std::map<UInt, Real> T;
T[0] = 2 * M_PI;
T[2] = M_PI;
T[4] = 3 * M_PI / 4.;
Real xi = q0 / q1;
Real zeta = q2 / q1;
std::vector<Real> moments;
using namespace std;
for (int q : {0, 2, 4}) {
Real m = T[q] * pow(q1, q - 2 * hurst) *
((1 - pow(xi, q + 2)) / (q + 2) +
(pow(zeta, q - 2 * hurst) - 1) / (q - 2 * hurst));
moments.push_back(m);
}
return moments;
}
template <>
std::vector<Real> Isopowerlaw<1>::moments() const {
TAMAAS_EXCEPTION("Moments have not been implemented for 1D surfaces");
}
/* -------------------------------------------------------------------------- */
template <UInt dim>
Real Isopowerlaw<dim>::alpha() const {
std::vector<Real> m = moments();
return m[0] * m[2] / (m[1] * m[1]);
}
/* -------------------------------------------------------------------------- */
template <UInt dim>
Real Isopowerlaw<dim>::rmsSlopes() const {
std::vector<Real> m = moments();
return std::sqrt(2 * m[1]);
}
template class Isopowerlaw<1>;
template class Isopowerlaw<2>;
__END_TAMAAS__
c4science · Help