Neko 1.99.3
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
most_cpu Module Reference

Implements the CPU kernel for the most_t type.

Data Types

interface  corr_h_interface
 
interface  corr_m_interface
 
interface  dfdl_interface
 
interface  f_interface
 
interface  slaw_h_interface
 
interface  slaw_m_interface
 

Functions/Subroutines

subroutine select_bc_operators (bc_type, bc_value, q, ts, ti, kappa, utau, z0h, hi, pr)
 Selects different expressions for the similarity functions in MOST based on the type of bottom boundary condition for temperature.
 
subroutine compute_ri_b (bc_type, g_dot_n, hi, ti, ts, magu, kappa, q, pr, ri_b)
 Computes the Richardson number.
 
subroutine set_stability_regime (ri_b, ri_threshold)
 Sets the stability regime based on the Richardson number value (quite arbitrary).
 
subroutine, public most_compute_cpu (u, v, w, temp, ind_r, ind_s, ind_t, ind_e, n_x, n_y, n_z, h, tau_x, tau_y, tau_z, n_nodes, lx, nelv, kappa, mu_w, rho_w, g_vec, pr, z0, z0h_in, bc_type, bc_value, tstep, ri_b_diagn, l_ob_diagn, utau_diagn, magu_diagn, ti_diagn, ts_diagn, q_diagn, h_x_idx, h_y_idx, h_z_idx)
 Main routine to compute the surface stresses based on MOST.
 
real(kind=rp) function slaw_m_stable (z, l_ob, z0)
 Similarity laws and corrections for the STABLE regime: REFERENCE: Holtslag, A. A. M., & De Bruin, H. A. R. (1988). Applied Modeling of the Nighttime Surface Energy Balance over Land. Journal of Applied Meteorology, 27(6), 689–704. NOTE: This formulation is chosen for its superior behavior in very stable conditions (large z/L), avoiding the numerical decoupling found in older linear functions (e.g., Dyer).
 
real(kind=rp) function slaw_h_stable (z, l_ob, z0h)
 
real(kind=rp) function corr_m_stable (z, l_ob)
 
real(kind=rp) function corr_h_stable (z, l_ob)
 
real(kind=rp) function slaw_m_convective (z, l_ob, z0)
 Similarity laws and corrections for the UNSTABLE (convective) regime: REFERENCE: Dyer, A. J. (1974), A review of flux-profile relationships, Bound.-Layer Meteorol., 7, 363-372. INTEGRATION: Paulson, C. A. (1970), The mathematical representation of wind speed and temperature profiles in the unstable atmospheric surface layer, J. Appl. Meteorol., 9, 857-861.
 
real(kind=rp) function slaw_h_convective (z, l_ob, z0h)
 
real(kind=rp) function corr_m_convective (z, l_ob)
 
real(kind=rp) function corr_h_convective (z, l_ob)
 
real(kind=rp) function slaw_m_neutral (z, l_ob, z0)
 Similarity laws and corrections for the NEUTRAL regime:
 
real(kind=rp) function slaw_h_neutral (z, l_ob, z0h)
 
real(kind=rp) function f_neumann (ri_b, z, z0, z0h, pr, l_ob, slaw_m, slaw_h)
 Simialrity laws (different for neumann and dirichlet bc's)
 
real(kind=rp) function dfdl_neumann (l_upper, l_lower, z, z0, z0h, pr, l_ob, slaw_m, slaw_h, fd_h)
 
real(kind=rp) function f_dirichlet (ri_b, z, z0, z0h, pr, l_ob, slaw_m, slaw_h)
 
real(kind=rp) function dfdl_dirichlet (l_upper, l_lower, z, z0, z0h, pr, l_ob, slaw_m, slaw_h, fd_h)
 

Variables

procedure(slaw_m_interface), pointer slaw_m_ptr => null()
 
procedure(slaw_h_interface), pointer slaw_h_ptr => null()
 
procedure(corr_m_interface), pointer corr_m_ptr => null()
 
procedure(corr_h_interface), pointer corr_h_ptr => null()
 
procedure(f_interface), pointer f_ptr => null()
 
procedure(dfdl_interface), pointer dfdl_ptr => null()
 

Function/Subroutine Documentation

◆ compute_ri_b()

subroutine most_cpu::compute_ri_b ( character(len=*), intent(in bc_type,
real(kind=rp), intent(in g_dot_n,
real(kind=rp), intent(in hi,
real(kind=rp), intent(in ti,
real(kind=rp), intent(in ts,
real(kind=rp), intent(in magu,
real(kind=rp), intent(in kappa,
real(kind=rp), intent(inout q,
real(kind=rp), intent(in pr,
real(kind=rp), intent(inout ri_b 
)
private

Definition at line 123 of file most_cpu.f90.

Here is the caller graph for this function:

◆ corr_h_convective()

real(kind=rp) function most_cpu::corr_h_convective ( real(kind=rp), intent(in z,
real(kind=rp), intent(in l_ob 
)
private

Definition at line 396 of file most_cpu.f90.

Here is the caller graph for this function:

◆ corr_h_stable()

real(kind=rp) function most_cpu::corr_h_stable ( real(kind=rp), intent(in z,
real(kind=rp), intent(in l_ob 
)
private

Definition at line 350 of file most_cpu.f90.

Here is the caller graph for this function:

◆ corr_m_convective()

real(kind=rp) function most_cpu::corr_m_convective ( real(kind=rp), intent(in z,
real(kind=rp), intent(in l_ob 
)
private

Definition at line 384 of file most_cpu.f90.

Here is the caller graph for this function:

◆ corr_m_stable()

real(kind=rp) function most_cpu::corr_m_stable ( real(kind=rp), intent(in z,
real(kind=rp), intent(in l_ob 
)
private

Definition at line 336 of file most_cpu.f90.

Here is the caller graph for this function:

◆ dfdl_dirichlet()

real(kind=rp) function most_cpu::dfdl_dirichlet ( real(kind=rp), intent(in l_upper,
real(kind=rp), intent(in l_lower,
real(kind=rp), intent(in z,
real(kind=rp), intent(in z0,
real(kind=rp), intent(in z0h,
real(kind=rp), intent(in pr,
real(kind=rp), intent(in l_ob,
procedure(slaw_m_interface slaw_m,
procedure(slaw_h_interface slaw_h,
real(kind=rp), intent(in fd_h 
)
private

Definition at line 454 of file most_cpu.f90.

Here is the caller graph for this function:

◆ dfdl_neumann()

real(kind=rp) function most_cpu::dfdl_neumann ( real(kind=rp), intent(in l_upper,
real(kind=rp), intent(in l_lower,
real(kind=rp), intent(in z,
real(kind=rp), intent(in z0,
real(kind=rp), intent(in z0h,
real(kind=rp), intent(in pr,
real(kind=rp), intent(in l_ob,
procedure(slaw_m_interface slaw_m,
procedure(slaw_h_interface slaw_h,
real(kind=rp), intent(in fd_h 
)
private

Definition at line 433 of file most_cpu.f90.

Here is the caller graph for this function:

◆ f_dirichlet()

real(kind=rp) function most_cpu::f_dirichlet ( real(kind=rp), intent(in ri_b,
real(kind=rp), intent(in z,
real(kind=rp), intent(in z0,
real(kind=rp), intent(in z0h,
real(kind=rp), intent(in pr,
real(kind=rp), intent(in l_ob,
procedure(slaw_m_interface slaw_m,
procedure(slaw_h_interface slaw_h 
)
private

Definition at line 445 of file most_cpu.f90.

Here is the caller graph for this function:

◆ f_neumann()

real(kind=rp) function most_cpu::f_neumann ( real(kind=rp), intent(in ri_b,
real(kind=rp), intent(in z,
real(kind=rp), intent(in z0,
real(kind=rp), intent(in z0h,
real(kind=rp), intent(in pr,
real(kind=rp), intent(in l_ob,
procedure(slaw_m_interface slaw_m,
procedure(slaw_h_interface slaw_h 
)
private

Definition at line 424 of file most_cpu.f90.

Here is the caller graph for this function:

◆ most_compute_cpu()

subroutine, public most_cpu::most_compute_cpu ( real(kind=rp), dimension(lx, lx, lx, nelv), intent(in u,
real(kind=rp), dimension(lx, lx, lx, nelv), intent(in v,
real(kind=rp), dimension(lx, lx, lx, nelv), intent(in w,
real(kind=rp), dimension(lx, lx, lx, nelv), intent(in temp,
integer, dimension(n_nodes), intent(in ind_r,
integer, dimension(n_nodes), intent(in ind_s,
integer, dimension(n_nodes), intent(in ind_t,
integer, dimension(n_nodes), intent(in ind_e,
real(kind=rp), dimension(n_nodes), intent(in n_x,
real(kind=rp), dimension(n_nodes), intent(in n_y,
real(kind=rp), dimension(n_nodes), intent(in n_z,
real(kind=rp), dimension(n_nodes), intent(in h,
real(kind=rp), dimension(n_nodes), intent(inout tau_x,
real(kind=rp), dimension(n_nodes), intent(inout tau_y,
real(kind=rp), dimension(n_nodes), intent(inout tau_z,
integer, intent(in n_nodes,
integer, intent(in lx,
integer, intent(in nelv,
real(kind=rp), intent(in kappa,
real(kind=rp), dimension(n_nodes), intent(in mu_w,
real(kind=rp), dimension(n_nodes), intent(in rho_w,
real(kind=rp), dimension(3), intent(in g_vec,
real(kind=rp), intent(in pr,
real(kind=rp), intent(in z0,
real(kind=rp), intent(in z0h_in,
character(len=*), intent(in bc_type,
real(kind=rp), intent(in bc_value,
integer, intent(in tstep,
real(kind=rp), dimension(n_nodes), intent(inout ri_b_diagn,
real(kind=rp), dimension(n_nodes), intent(inout l_ob_diagn,
real(kind=rp), dimension(n_nodes), intent(inout utau_diagn,
real(kind=rp), dimension(n_nodes), intent(inout magu_diagn,
real(kind=rp), dimension(n_nodes), intent(inout ti_diagn,
real(kind=rp), dimension(n_nodes), intent(inout ts_diagn,
real(kind=rp), dimension(n_nodes), intent(inout q_diagn,
integer, dimension(n_nodes), intent(in h_x_idx,
integer, dimension(n_nodes), intent(in h_y_idx,
integer, dimension(n_nodes), intent(in h_z_idx 
)
Parameters
tstepThe current time-step

Definition at line 161 of file most_cpu.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ select_bc_operators()

subroutine most_cpu::select_bc_operators ( character(len=*), intent(in bc_type,
real(kind=rp), intent(in bc_value,
real(kind=rp), intent(inout q,
real(kind=rp), intent(inout ts,
real(kind=rp), intent(in ti,
real(kind=rp), intent(in kappa,
real(kind=rp), intent(in utau,
real(kind=rp), intent(in z0h,
real(kind=rp), intent(in hi,
real(kind=rp), intent(in pr 
)
private

Definition at line 101 of file most_cpu.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ set_stability_regime()

subroutine most_cpu::set_stability_regime ( real(kind=rp), intent(in ri_b,
real(kind=rp), intent(in ri_threshold 
)
private

Definition at line 140 of file most_cpu.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ slaw_h_convective()

real(kind=rp) function most_cpu::slaw_h_convective ( real(kind=rp), intent(in z,
real(kind=rp), intent(in l_ob,
real(kind=rp), intent(in z0h 
)
private

Definition at line 377 of file most_cpu.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ slaw_h_neutral()

real(kind=rp) function most_cpu::slaw_h_neutral ( real(kind=rp), intent(in z,
real(kind=rp), intent(in l_ob,
real(kind=rp), intent(in z0h 
)
private

Definition at line 416 of file most_cpu.f90.

Here is the caller graph for this function:

◆ slaw_h_stable()

real(kind=rp) function most_cpu::slaw_h_stable ( real(kind=rp), intent(in z,
real(kind=rp), intent(in l_ob,
real(kind=rp), intent(in z0h 
)
private

Definition at line 329 of file most_cpu.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ slaw_m_convective()

real(kind=rp) function most_cpu::slaw_m_convective ( real(kind=rp), intent(in z,
real(kind=rp), intent(in l_ob,
real(kind=rp), intent(in z0 
)
private

Definition at line 370 of file most_cpu.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

◆ slaw_m_neutral()

real(kind=rp) function most_cpu::slaw_m_neutral ( real(kind=rp), intent(in z,
real(kind=rp), intent(in l_ob,
real(kind=rp), intent(in z0 
)
private

Definition at line 409 of file most_cpu.f90.

Here is the caller graph for this function:

◆ slaw_m_stable()

real(kind=rp) function most_cpu::slaw_m_stable ( real(kind=rp), intent(in z,
real(kind=rp), intent(in l_ob,
real(kind=rp), intent(in z0 
)
private

Definition at line 322 of file most_cpu.f90.

Here is the call graph for this function:
Here is the caller graph for this function:

Variable Documentation

◆ corr_h_ptr

procedure(corr_h_interface), pointer most_cpu::corr_h_ptr => null()
private

Definition at line 93 of file most_cpu.f90.

◆ corr_m_ptr

procedure(corr_m_interface), pointer most_cpu::corr_m_ptr => null()
private

Definition at line 92 of file most_cpu.f90.

◆ dfdl_ptr

procedure(dfdl_interface), pointer most_cpu::dfdl_ptr => null()
private

Definition at line 95 of file most_cpu.f90.

◆ f_ptr

procedure(f_interface), pointer most_cpu::f_ptr => null()
private

Definition at line 94 of file most_cpu.f90.

◆ slaw_h_ptr

procedure(slaw_h_interface), pointer most_cpu::slaw_h_ptr => null()
private

Definition at line 91 of file most_cpu.f90.

◆ slaw_m_ptr

procedure(slaw_m_interface), pointer most_cpu::slaw_m_ptr => null()
private

Definition at line 90 of file most_cpu.f90.