MOSD.H

Go to the documentation of this file.

#ifndef MOSD_H
#define MOSD_H
 
#include "AMReX_AmrCore.H"
#include "AMReX.H"
#include "AMReX_REAL.H"
 
using namespace amrex::literals;
namespace amr_wind {
struct MOSD
{
    /*
     * A dynamic wall model that calculates the stress from wave to wind
     * based on geometric information of the wave.
     */
 
    amrex::Real amplitude{0.05_rt};
    amrex::Real wavenumber{4.0_rt};
    amrex::Real omega{0.8_rt};
    amrex::Real time;
 
    AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE amrex::Real get_dyn_tau(
        const amrex::Real u_dx,
        const amrex::Real v_dx,
        const amrex::Real xc,
        const amrex::Real unit_nor) const
    {
 
        // Building the wave surface, gradients, wave velocities and unit normal
        const amrex::Real dx_eta_wave =
            -amplitude * wavenumber * std::sin(wavenumber * xc - omega * time);
        const amrex::Real dy_eta_wave = 0;
        const amrex::Real dt_eta_wave =
            amplitude * omega * std::sin(wavenumber * xc - omega * time);
        const amrex::Real grad_eta_wave =
            std::sqrt(dx_eta_wave * dx_eta_wave + dy_eta_wave * dy_eta_wave);
        const amrex::Real Cx_wave = -dt_eta_wave * dx_eta_wave *
                                    (1.0_rt / (grad_eta_wave * grad_eta_wave));
        const amrex::Real Cy_wave = -dt_eta_wave * dy_eta_wave *
                                    (1.0_rt / (grad_eta_wave * grad_eta_wave));
        const amrex::Real n_x = dx_eta_wave / grad_eta_wave;
        const amrex::Real n_y = dy_eta_wave / grad_eta_wave;
 
        // Calculating the relative velocity, heaviside function
        const amrex::Real u_r = u_dx - Cx_wave;
        const amrex::Real v_r = v_dx - Cy_wave;
        const amrex::Real ur_mag =
            std::sqrt(u_r * u_r * n_x * n_x + v_r * v_r * n_y * n_y);
        const amrex::Real Heavi_arg = (u_r * dx_eta_wave + v_r * dy_eta_wave);
        const amrex::Real Heavi =
            (Heavi_arg + std::abs(Heavi_arg)) / (2.0_rt * Heavi_arg);
 
        // Calculating the magnitude of the stress
        return (1.0_rt / static_cast<amrex::Real>(M_PI)) * ur_mag * ur_mag *
               grad_eta_wave * grad_eta_wave * Heavi *
               (unit_nor == 0 ? n_x : n_y);
    }
};
 
} // namespace amr_wind
#endif /* MOSD_H */