Newer
Older
#include "../includeC/sfx_data.h"
#include "../includeCU/sfx_math.cuh"
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
template<typename T>
FUCNTION_DECLARATION_SPECIFIER void get_thermal_roughness(T &z0_t, T &B,
const T z0_m, const T Re,
const struct sfx_surface_param& param,
const int surface_type)
{
// --- define B = log(z0_m / z0_t)
B = (Re <= param.Re_rough_min) ? param.B1_rough * log(param.B3_rough * Re) + param.B2_rough :
param.B4_rough * (pow(Re, param.B2_rough));
// --- apply max restriction based on surface type
if (surface_type == param.surface_ocean) B = sfx_math::min(B, param.B_max_ocean);
if (surface_type == param.surface_lake) B = sfx_math::min(B, param.B_max_land);
if (surface_type == param.surface_land) B = sfx_math::min(B, param.B_max_lake);
// --- define roughness [thermal]
z0_t = z0_m / exp(B);
}
template<typename T>
FUCNTION_DECLARATION_SPECIFIER void get_charnock_roughness(T &z0_m, T &u_dyn0,
const T h, const T U,
const sfx_esm_param& model_param,
const sfx_surface_param& surface_param,
const sfx_esm_numericsTypeC& numerics)
{
T Uc, a, b, c, c_min, f;
Uc = U;
a = 0.0;
b = 25.0;
c_min = log(surface_param.h_charnock) / model_param.kappa;
for (int i = 0; i < numerics.maxiters_charnock; i++)
{
f = surface_param.c1_charnock - 2.0 * log(Uc);
for (int j = 0; j < numerics.maxiters_charnock; j++)
{
c = (f + 2.0 * log(b)) / model_param.kappa;
if (U <= 8.0e0)
a = log(1.0 + surface_param.c2_charnock * ( pow(b / Uc, 3) ) ) / model_param.kappa;
c = sfx_math::max(c - a, c_min);
b = c;
}
z0_m = surface_param.h_charnock * exp(-c * model_param.kappa);
z0_m = sfx_math::max(z0_m, T(0.000015e0));
Uc = U * log(surface_param.h_charnock / z0_m) / log(h / z0_m);
}
u_dyn0 = Uc / c;
}