cosine_wave.cc

#include <vtestbed/config/real_type.hh>
#include <vtestbed/core/cosine_wave.hh>

template <class T,int N>
T
vtb::core::Standing_cosine_wave<T,N>::velocity_potential(
    const vec& x,
    T t,
    T z
) {
    using blitz::dot;
    using std::sin;
    using std::sinh;
    const T k_length = length(this->wave_number());
    const T& angfreq = this->angular_frequency();
    return this->amplitude() * angfreq
           * cos(dot(this->wave_number(), x)) * sin(angfreq*t)
           * sinh(k_length * z) / k_length;
}

template <class T,int N>
auto
vtb::core::Standing_cosine_wave<T,N>::max_velocity() const -> vec {
    using std::abs;
    const T& a = this->amplitude();
    const auto& k = this->wave_number();
    const T k_length = length(k);
    const T phi = abs(a * this->angular_frequency()
        * sinh(k_length * a) / k_length);
    return {k(0)*phi, k(1)*phi};
}

template <class T,int N>
vtb::core::Grid<T,3>
vtb::core::Standing_cosine_wave<T,N>::acf_grid() const {
    static_assert(N == 3, "bad N");
    using blitz::ceil;
    using blitz::max;
    using blitz::where;
    using std::ceil;
    auto length = this->length();
    const auto& period = this->period();
    const auto& mask = isfinite(length);
    T length_max = max(where(mask, length, T{0}));
    length = T{0.75}*where(mask, length, length_max);
    blitz::TinyVector<T,2> npoints = ceil(length) + 1;
    Grid<T,N> grid{
        {0,0,0},
        {period,length(0),length(1)},
        {2,int(npoints(0)),int(npoints(1))}
    };
    T dt = T{0.25} / courant_number_sum(*this, grid);
    grid.end(0) = int(ceil(period / dt)) + 1;
    std::clog << "dt=" << dt << std::endl;
    std::clog << "dx=" << grid.delta(1) << std::endl;
    std::clog << "dy=" << grid.delta(2) << std::endl;
    std::clog << "grid=" << grid << std::endl;
    T courant = courant_number(*this, grid);
    std::clog << "courant=" << courant << std::endl;
    return grid;
}

template class vtb::core::Standing_cosine_wave<VTB_REAL_TYPE,3>;
#if defined(VTB_REAL_TYPE_FLOAT)
template class vtb::core::Standing_cosine_wave<double,3>;
#endif