Newer
Older
#include "../includeCU/sfx_surface.cuh"
#include "../includeCU/sfx_esm_compute_subfunc.cuh"
FluxEsmBase<T, memIn, memOut, RunMem>::FluxEsmBase(sfxDataVecTypeC* sfx_in,
meteoDataVecTypeC* meteo_in,
const sfx_esm_param model_param_in,
const sfx_surface_param surface_param_in,
const sfx_esm_numericsTypeC numerics_in,
const sfx_phys_constants phys_constants_in,
const int grid_size_in) : ModelBase<T, memIn, memOut, RunMem>(sfx_in, meteo_in, grid_size_in)
surface_param = surface_param_in;
phys_constants = phys_constants_in;
template<typename T, MemType memIn, MemType memOut >
void FluxEsm<T, memIn, memOut, MemType::CPU>::compute_flux()
T h, U, dT, Tsemi, dQ, z0_m;
T Re, z0_t, B, h0_m, h0_t, u_dyn0, zeta, Rib, zeta_conv_lim, Rib_conv_lim, f_m_conv_lim, f_h_conv_lim, psi_m, psi_h, phi_m, phi_h, Km, Pr_t_inv, Cm, Ct;
int surface_type;
T fval;
for (int step = 0; step < grid_size; step++)
U = meteo.U[step];
Tsemi = meteo.Tsemi[step];
dT = meteo.dT[step];
dQ = meteo.dQ[step];
h = meteo.h[step];
z0_m = meteo.z0_m[step];
surface_type = z0_m < 0.0 ? surface_param.surface_ocean : surface_param.surface_land;
if (surface_type == surface_param.surface_ocean)
{
get_charnock_roughness(z0_m, u_dyn0, h, U, model_param, surface_param, numerics);
h0_m = h / z0_m;
}
if (surface_type == surface_param.surface_land)
{
h0_m = h / z0_m;
u_dyn0 = U * model_param.kappa / log(h0_m);
}
Re = u_dyn0 * z0_m / phys_constants.nu_air;
get_thermal_roughness(z0_t, B, z0_m, Re, surface_param, surface_type);
h0_t = h / z0_t;
Rib = (phys_constants.g / Tsemi) * h * (dT + 0.61e0 * Tsemi * dQ) / (U*U);
get_convection_lim(zeta_conv_lim, Rib_conv_lim, f_m_conv_lim, f_h_conv_lim,
h0_m, h0_t, B,
model_param);
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
if (Rib > 0.0)
{
Rib = std::min(Rib, model_param.Rib_max);
get_psi_stable(psi_m, psi_h, zeta, Rib, h0_m, h0_t, B, model_param);
fval = model_param.beta_m * zeta;
phi_m = 1.0 + fval;
phi_h = 1.0/model_param.Pr_t_0_inv + fval;
}
else if (Rib < Rib_conv_lim)
{
get_psi_convection(psi_m, psi_h, zeta, Rib, h0_m, h0_t, B, zeta_conv_lim, f_m_conv_lim, f_h_conv_lim, model_param, numerics);
fval = pow(zeta_conv_lim / zeta, 1.0/3.0);
phi_m = fval / f_m_conv_lim;
phi_h = fval / (model_param.Pr_t_0_inv * f_h_conv_lim);
}
else if (Rib > -0.001)
{
get_psi_neutral(psi_m, psi_h, zeta, h0_m, h0_t, B, model_param);
phi_m = 1.0;
phi_h = 1.0 / model_param.Pr_t_0_inv;
}
else
{
get_psi_semi_convection(psi_m, psi_h, zeta, Rib, h0_m, h0_t, B, model_param, numerics);
phi_m = pow(1.0 - model_param.alpha_m * zeta, -0.25);
phi_h = 1.0 / (model_param.Pr_t_0_inv * sqrt(1.0 - model_param.alpha_h_fix * zeta));
}
Cm = model_param.kappa / psi_m;
Ct = model_param.kappa / psi_h;
Km = model_param.kappa * Cm * U * h / phi_m;
Pr_t_inv = phi_m / phi_h;
sfx.zeta[step] = zeta;
sfx.Rib[step] = Rib;
sfx.Re[step] = Re;
sfx.B[step] = B;
sfx.z0_m[step] = z0_m;
sfx.z0_t[step] = z0_t;
sfx.Rib_conv_lim[step] = Rib_conv_lim;
sfx.Cm[step] = Cm;
sfx.Ct[step] = Ct;
sfx.Km[step] = Km;
sfx.Pr_t_inv[step] = Pr_t_inv;
}
if(MemType::CPU != memOut)
memproc::memcopy<memOut, MemType::CPU>((void*&)res_sfx->zeta, (void*&)sfx.zeta, new_size);
memproc::memcopy<memOut, MemType::CPU>((void*&)res_sfx->Rib, (void*&)sfx.Rib, new_size);
memproc::memcopy<memOut, MemType::CPU>((void*&)res_sfx->Re, (void*&)sfx.Re, new_size);
memproc::memcopy<memOut, MemType::CPU>((void*&)res_sfx->B, (void*&)sfx.B, new_size);
memproc::memcopy<memOut, MemType::CPU>((void*&)res_sfx->z0_m, (void*&)sfx.z0_m, new_size);
memproc::memcopy<memOut, MemType::CPU>((void*&)res_sfx->z0_t, (void*&)sfx.z0_t, new_size);
memproc::memcopy<memOut, MemType::CPU>((void*&)res_sfx->Rib_conv_lim, (void*&)sfx.Rib_conv_lim, new_size);
memproc::memcopy<memOut, MemType::CPU>((void*&)res_sfx->Cm, (void*&)sfx.Cm, new_size);
memproc::memcopy<memOut, MemType::CPU>((void*&)res_sfx->Ct, (void*&)sfx.Ct, new_size);
memproc::memcopy<memOut, MemType::CPU>((void*&)res_sfx->Km, (void*&)sfx.Km, new_size);
memproc::memcopy<memOut, MemType::CPU>((void*&)res_sfx->Pr_t_inv, (void*&)sfx.Pr_t_inv, new_size);
template class FluxEsm<float, MemType::CPU, MemType::CPU, MemType::CPU>;
template class FluxEsmBase<float, MemType::CPU, MemType::CPU, MemType::CPU>;
template class FluxEsmBase<float, MemType::GPU, MemType::GPU, MemType::GPU>;
template class FluxEsmBase<float, MemType::GPU, MemType::GPU, MemType::CPU>;
template class FluxEsmBase<float, MemType::GPU, MemType::CPU, MemType::GPU>;
template class FluxEsmBase<float, MemType::CPU, MemType::GPU, MemType::GPU>;
template class FluxEsmBase<float, MemType::CPU, MemType::CPU, MemType::GPU>;
template class FluxEsmBase<float, MemType::CPU, MemType::GPU, MemType::CPU>;
template class FluxEsmBase<float, MemType::GPU, MemType::CPU, MemType::CPU>;
template class FluxEsm<float, MemType::GPU, MemType::GPU, MemType::GPU>;
template class FluxEsm<float, MemType::GPU, MemType::GPU, MemType::CPU>;
template class FluxEsm<float, MemType::GPU, MemType::CPU, MemType::GPU>;
template class FluxEsm<float, MemType::CPU, MemType::GPU, MemType::GPU>;
template class FluxEsm<float, MemType::CPU, MemType::CPU, MemType::GPU>;
template class FluxEsm<float, MemType::CPU, MemType::GPU, MemType::CPU>;
template class FluxEsm<float, MemType::GPU, MemType::CPU, MemType::CPU>;