cvmix / cvmix-src Goto Github PK
View Code? Open in Web Editor NEWCVMix source code (including protex documentation) as well as build system and examples / tests
License: Other
CVMix source code (including protex documentation) as well as build system and examples / tests
License: Other
Need a general cvmix_finalize to deallocate memory (or clean up pointers?) as well as method-specific cvmix_MODULE_finalize routines for methods that allocate memory (background mixing is the only one for now).
The Makefiles I've written get the job done, though occasionally lack of maintenance has led to issues with some of the compilers / targets. It would probably be cleaner and easier to manage if we switched to something like CMake
When computing the bulk Richardson number, we want
Rib = (d-dr)*(Br-B(d))/( |Vr-V(d)|^2 + Vt^2(d) )
Where d is the current cell-center depth and dr is the center of the surface layer. As @toddringler pointed out, dr = 0.5_surf_layer_ext_d, so instead of adding an additional argument to cvmix_kpp_compute_bulk_Richardson() we'll compute
Rib = (1-0.5_surf_layer_ext)_d*(Br-B(d))/( |Vr-V(d)|^2 + Vt^2(d) )
In cvmix_kpp_coeffs_low, the indexing in variables like Gat1, DGat1, visc_at_OBL isn't very clear... I think I might put in a module-wide parameter index (CVMIX_KPP_MDIFF_IND, etc) to clarify.
The non-local transport returned from cvmix_coeffs_kpp_low() contains a NaN in the first element (surface). Bug appeared in commit 4d2d987. Bug is still manifest in 4f3943a.
A possibly related change is that the diffusivity at the surface is no longer zero but a smallish (not tiny 1e-5) value.
Configuration uses MatchTechnique='SimpleShapes'.
The last working commit was 7ec9961.
See https://github.com/CVMix/CVMix-src/tree/feature/kpp_optimized_for_openmp_no_allocates for tips on removing allocate calls that slow down the code (GFDL reported 30% speed up using this branch rather than corresponding trunk)
static_diff(:,2) and static_visc(:) don't mesh with the renamed Tdiff, Sdiff, and Mdiff.
Per a discussion on September 23, 2013 -- should we add new interface for compute_turbulent_shear() to allow for array of OBL_depth (used for computing Ri_bulk)? In more general terms, should more routines have interfaces to handle multiple columns at once?
Shear mixing needs two background mixing parameters for Pacanowski-Philander. There are three options for how to handle this:
Currently option (3) is in use, but I think it's the worst option of the group. I'm not sure which of the first two would be the best to use, though. I'm leaning towards (2) because that already has all the needed support for 1D or 2D variations in the background field.
We need to use N / N^2 at the bottom interface of a cell when computing unresolved shear because the parameter Cv (and the term Cv*N) is defined in this manner. Computing N at the cell center would require recomputing Cv.
We cleaned up KPP before the beta release, but other modules still have allocatable arrays that should have a declared size...
I believe cvmix_kpp requires Sdiff at this time, but the actual KPP routine doesn't require double diffusion and therefore flops / memory can be saved for users who only have momentum and tracer diffusion coefficients.
While updating cvmix_shear to separate inputs and outputs / have a low-level and a wrapped version of cvmix_coeffs_shear, I don't know that I handled the "no_diff" operation correctly. I'll go back and check once the rest of the modules have separate inputs / outputs and low-level coeff routines.
The original PP scheme includes a background term that is added on.
The implementation in CVMix is the same and, as a result, breaks the "a al carte" approach that we strive for in CVMix. In order to run the CVMix shear PP module we have to initialize CVMix background module, i.e. shear depends on background.
I propose that we drop the background part of PP. Most everyone who uses CVMix shear will be using some form of background mixing that is set outside of PP. For example, right now MPAS-O has to test if shear mixing is using PP and, if so, subtract out the background from the result (because we have added it elsewhere.)
It will lead to clearer code all around to simply remove the addition of background from CVMix shear PP. (Note that it is not added into CVMix shear KPP).
v0.53:
cvmix_kinds_and_types.F90 says:
! For KPP we need buoyancy (as opposed to buoyancy frequency)
then defines:
real(cvmix_r8), dimension(:), pointer :: buoyancy_cntr => NULL() ! units: m s^-2
but I can not find where/how buoyancy_cntr is used in the code. Is buoyancy_cntr used in CVMix/KPP?
a "grep buoyancy_cntr *.F90" turns up nothing except puts/gets.
In LMD94, the nonlocal term is defined (by manipulating Eqs (10) and (20) by
K*gamma = -Cs * G(sigma) * Q
(the minus sign is the CVMix convention and CVMix returns K*gamma/Q to the base model).
@StephenGriffies pointed out that LMD assumes G(0) = 0 => the nonlocal term is 0 at the surface... but it may be desirable to use a shape function that allows K*gamma = -Q at the surface.
I think there are two ways to accomplish this. Both start with removing cvmix_kpp_compute_nonlocal() [which returns G(sigma)*Cs] and simply computing G(sigma). From there, we can either:
As Steve and I discussed over the phone, we can add Cs to the cvmix_kpp_params_type. When you call cvmix_init_kpp(), if you omit the Cs argument then it computes Cs = C* * kappa * (__)^1/3 or you can set it to 1
Something else that I just thought of would be to create Cs as a module-wide parameter. If you choose MatchTechnique = 'ParabolicNonlocal' then Cs = 1, otherwise Cs = the value above.
In either case, you multiply Tnonlocal and Snonlocal by Cs after computing all the G(sigma) values in the column.
Relatedly, Steve has suggested two other shape functions for the nonlocal term, a linear function and a cubic function. I think the best way to implement those would be to have another MatchTechnique parameter (so MatchTechnique controls the local shape function and MatchTechnique2 or MatchTechniqueNonlocal controls the nonlocal shape function). I'll tackle that after fixing this issue with Cs, regardless of whether we go with option (1) or option (2) it should be a straight-forward update.
A few functions within the CVMix KPP interface take the argument max_nlev to allow KPP to be computed over a subset of the entire column (e.g. cvmix_coeffs_kpp).
This option should be available to all routines (e.g. compute_turbulent_scales, compute_bulk_richardson).
Even when the max_nlev option is included in a routine, there are portions of that routine that still compute over nlevs For example, within coeffs_kpp, compute_turbulent_scales is called it is over nlevs (determined by the size of zw_ctr). Thus, if max_nlev != nlevs and all of zw_ctr is sent, an issue can arise.
Testing should be split into two directories:
Also, should provide some light-weight examples for both stand-alone single column cases and (if possible) an outline of steps needed to import into a pre-existing ocean model.
Random thoughts:
Instead of computing the enhanced diffusivity when using "MatchBoth" (and no other time), introduce a new parameter [cvmix_kpp_params%enhance_diff = .true. / .false.]. For now, we can default to enhance_diff = .true. but that may be changed later (or depending on choice of MatchTechnique).
Right now, POP+CVMix can only run with a single thread per MPI task (by default, the NCAR machine runs POP with 2 threads per MPI task). At this point, I don't know if the issue is on the POP side or the CVMix side.
The cvmix_put_kpp_real routine has a different structure from every other cvmix_put_[module] routine... it uses an optional cvmix_kpp_params_type argument and otherwise allows the user to edit variables in the private cvmix_kpp_params_saved variable. We want these routines to be consistent from module to module, and I like the optional argument structure.
Obtained with intel, g77 and pgi compilers using -g. I suspect that OBL_Mdiff, OBL_Tdiff and OBL_Sdiff should be allocated with one extra slot.
Image PC Routine Line Source
MOM6 00000000097E558A cvmix_kpp_mp_cvmi 797 cvmix_kpp.F90
MOM6 000000000689F111 mom_kpp_mp_kpp_ca 617 MOM_KPP.F90
MOM6 0000000008D6016C mom_diabatic_driv 473 MOM_diabatic_driver.F90
MOM6 000000000B99E375 mom_mp_step_mom_ 803 MOM.F90
MOM6 0000000009B52933 ocean_model_mod_m 368 ocean_model_MOM.F90
MOM6 000000000154C038 MAIN__ 695 coupler_main.F90
forrtl: severe (408): fort: (2): Subscript #1 of the array OBL_MDIFF has value 6 which is greater than the upper bound of 5
'''
We need to pass either N_iface or Nsqr_iface to the routine but we want to do computations at cell center. As discussed on the weekly telecon, we want to do all averaging to N_iface the latter is not appropriate for grids with variable vertical resolutions... this means changing
N_cntr(kt)=sqrt((max(Nsqr_iface(kt),cvmix_zero) + &
max(Nsqr_iface(kt+1),cvmix_zero)) * &
0.5_cvmix_r8)
to
N_cntr(kt)=(sqrt(max(Nsqr_iface(kt),cvmix_zero)) + &
sqrt(max(Nsqr_iface(kt+1),cvmix_zero))) * &
0.5_cvmix_r8
Also, we want the default value for lavg_N_or_Nsqr to be .true. rather than .false.
Each module should have a private saved cvmix_[module]_params_type, and each public routine that takes that type as an argument should do so optionally - if an argument is present, it provides whatever parameters are needed otherwise the module-wide structure is used.
This option currently exists in
but is missing in
The variable naming convention should be the following:
Exceptions to the above - diffusivities will be slightly abbreviated:
Same standards should be applied to intent(in)s and intent(out)s.
The Pacanowski-Philander (PP) version of shear mixing is based on the gradient Richardson number (RICH). The scheme proposed in PP (1981) clearly assumes RICH >= 0. But within CVMix, RICH is not tested to be 0 or positive when using PP. Note that when using shear KPP, RICH<0 is tested, so CVMix is assuming RICH can be either positive or negative.
I only ran into this issue when using KPP, but the problem ends up being a combination of lines:
451 and 600-602 in https://github.com/CVMix/CVMix-src/blob/master/src/shared/cvmix_kpp.F90
Basically, the arrays new_Mdiff (and such) are allocated using nlev, which is the number of active levels in a column, not necessarily the number of total levels in the column. Later a direct copy happens moving old_Mdiff into Mdiff_out, but they have different sizes (one is active, one is total) which causes an array out of bounds with some compilers. Ideally, the arrays would be allocated using the total number of levels whether or not they are used.
Many cvmix_kpp routines were made public for testing, but really should only be callable from the cvmix_coeffs_kpp() routine [shape functions, kOBL_depth, enhanced_diff, compute_nonlocal, compute_nu_at_OBL_depth]
As I have been looking through CVMIX, I have an issue with the following code in cvmix_kpp_compute_turbulent_scales_1d. In it, zeta (sigma * h / L) is computed as follows
if (surf_fric_vel.ne.cvmix_zero) then
do kw=1,n_sigma
! compute scales at sigma if sigma < surf_layer_ext, otherwise compute
! at surf_layer_ext
zeta(kw) = min(surf_layer_ext, sigma_coord(kw)) * OBL_depth * &
surf_buoy_force*vonkar/(surf_fric_vel**3)
end do
I am wondering why the convective case limitation on sigma is included for stable buoyancy forcing as well? Large et al. (1994), Troen and Mahrt (1986), and the CVMIX description document assume that the similarity functions for the surface layer extend from 0 < sigma < 1
I have run a test case with positive buoyancy forcing and a constant wind stress in a SCM framework. Below is a plot of temperature profiles through time, with the kpp diagnosed boundary layer depth in black.
The corresponding result from a similarly forced LES is below (the line is a diagnosed mixed layer depth from a temperature threshold criterion)
I have then run a test that changed the CVMIX code above to
if (surf_fric_vel.ne.cvmix_zero) then
do kw=1,n_sigma
! compute scales at sigma if sigma < surf_layer_ext, otherwise compute
! at surf_layer_ext
if(surf_buoy_force .ge. cvmix_zero) then
zeta(kw) = sigma_coord(kw) * OBL_depth * &
surf_buoy_force*vonkar/(surf_fric_vel**3)
else
zeta(kw) = min(surf_layer_ext, sigma_coord(kw)) * OBL_depth * &
surf_buoy_force*vonkar/(surf_fric_vel**3)
endif
end do
with this change, the temperature profile plot becomes
which is closer to the LES result.
It seems that the in stable conditions, sigma should not be confined less than the surface layer extent.
Need physical parameters for testing KPP. Ideally, would have multiple sets of surface buoyancy force, surface friction velocity, OBL_depths, and interior mixing values to exercise as many branches of the KPP code as possible. Specifically, for computing turbulent velocity scales, would like to test the case where surf_fric_vel = 0 & surf_buoy_force < 0 as well other general cases.
@StephenGriffies pointed out two parameters that might be useful:
A minimum OBL_depth
A prescribed OBL_depth -- this could be helpful in diagnosing differences in the WSwPSBF tests
Need clear instructions for members of community that want to add mixing techniques to the code. For example:
We want each module to have a clear owner so that users can ask the author directly instead of going through us.
The Simmons scheme uses a static tidal energy field, and the vertical deposition is also static in time, so we only need to compute energy_flux*vert_dep/rho_fw once
(rather than every time step, as is currently the case). I'd like to refactor the tidal mixing module to include a compute_tidal_diffusivity()
routine that can be during the initialization phase of the GCM to compute this value and store it in the cvmix_data_type
structure, then cvmix_coeffs_tidal_low()
can reference the stored data rather than computing it every time.
If you want to write a variable to netCDF, you should be able to specify one of many names (for example, "SqrBuoyancy", "SqrBuoyancyFreq", "SqrBuoyancyFreq_iface", "buoyancy_iface", etc) but those should all result in a variable with the same name ("SqrBuoyancyFreq"?) being made in the netCDF file. Similarly, you should be able to edit attributes using a wide range of names.
I think the best way to do this will be to create a simplify_name routine so that the multiple names only need to be handled in a single routine.
Nothing in the reg_test/ directory exercises code in cvmix_convection module.
I think the code is more readable if all routines are written in the same style, and the put routines exist for a reason.
Need a clean way to have CVMix optionally output diagnostic variables.
Need to make sure parameters are set with reasonable values (in cvmix_[module]_init routines?). This is already done in some situations: kpp makes sure c_s and c_m are positive and that surf_layer_ext is in [0,1].
Add parameter for min_Vtsqr (default 1e-10 m^2/s^2), have cvmix_kpp_compute_unresolved_shear return
min(computed Vt^2, min_Vtsqr)
Currently cvmix_kpp is the only module that allows users to pass in arrays rather than requiring users to pack data in to the CVmix_data_type structure. The current cvmix_[module]coeffs should be renamed cvmix[module]coeffs_wrap and it should call cvmix[module]_coeffs_low, where computation is actually done.
Per a discussion on September 23, 2013 -- the cvmix_put and cvmix_get routines improve the readability of the code but result in a performance hit. Is that a worthwhile trade-off?
We should make sure the shape of any arrays being stored via cvmix_put are the same shape as the data type attribute they are being stored in.
Example:
allocate(depth(5))
cvmix_put(CVMix_vars, 'nlev', 5)
cvmix_put(CVMix_vars, 'zt', depth)
allocate(depth(8))
cvmix_put(CVMix_vars, 'nlev', 5)
cvmix_put(CVMix_vars, 'zt', depth)
allocate(depth(5))
cvmix_put(CVMix_vars, 'nlev', 8)
cvmix_put(CVMix_vars, 'zt', depth)
Obviously some system is needed to ensure code on the master branch compiles and runs as expected... ideally with little or no human intervention.
Necessary Features:
Some discussion has taken place on the developer's mailing list:
Besides needing to create this testing infrastructure, I think we need to completely overhaul reg_tests/ so that it will play nice with the test system and be more amenable to introducing new regression tests.
There have been a lot of changes to cvmix_coeffs_kpp_low between June 18 (4d2d987) and June 27 (61292f1), and it has changed answers in the GFDL single-column test case... Before making the beta tag, we want to make sure we understand what is happening in the routine with the GFDL default options and that everything looks okay.
Shear mixing will error out if you do not set all the necessary parameters while every other module sets the parameters to default values if they are not provided by the user.
The way CVMix handles the nonlocal term doesn't account for enhanced diffusivity... need to multiply by Kenh/Kold at bottom of boundary layer
After first public release, we should look into vectorizing the code. (This is a pretty vague ticket, details will be added via comments.)
Low-level routines (cvmix_[module]_coeffs) should have separate inputs and outputs for the diffusivities so that the user can determine what to do with new diffusivities (overwrite old, add to old, etc). The interfaces using the CVMix wrapper should have some parameter flag to determine what to return (see cvmix_background and the 'handle_old_vals' parameter).
POP's convective mixing scheme is only applied below the boundary layer. I have a branch of the repository (on my personal fork) where I've added a new option to the convective mixing initialization named lnoOBL
, with a default value of .true.
implying we do not apply in the boundary layer. I don't know if
(a) this is a feature we want in CVMix (alternatively, GCMs can zero out the convective mixing contribution in the OBL), or
(b) if this is the best implementation.
Some shortfallings -- currently, lnoOBL
is only checked if we are using Brunt-Vaisala convective mixing with BVsqr_convect
>= 0. Once I have a better idea on the desired application / defaults / variable names I'll clean up the code and apply it regardless of convective mixing scheme.
We've talked via email about adding an error code as an optional argument (to allow modules to handle errors on their own rather than getting a stop call which is sure to cause issues with MPI). There would also be a cvmix_error module in src/ to handle the errors for the stand-alone drivers... likely just printing a statement followed by a stop call.
On the same note, we should put some thought into handling warnings - I think we should print warnings in a variety of cases (over-writing existing parameters, using cvmix_put to write a scalar into an array, etc)... but there should be an easy way to suppress the warnings. I'm leaning towards a namelist option that sets the global variable lwarn, but I'm open to other ideas - a namelist variable stored somewhere else? Copious ifdefs along with a compile-time setting? Going with the namelist option, I'm picturing a cvmix_warn routine that takes a string as an argument and prints "WARNING: " // warn_string... with the check against lwarn or whatever in the cvmix_warn routine. Maybe something like
module cvmix_warn
public :: cvmix_warn
public :: cvmix_set_lwarn
logical, private, save :: lwarn = .true.
...
end module cvmix_warn
So by default we print all warnings, but the user can call cvmix_set_lwarn(.false.) to suppress them (and then later call cvmix_set_lwarn(.true.) to turn them back on).
A declarative, efficient, and flexible JavaScript library for building user interfaces.
๐ Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
An Open Source Machine Learning Framework for Everyone
The Web framework for perfectionists with deadlines.
A PHP framework for web artisans
Bring data to life with SVG, Canvas and HTML. ๐๐๐
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
Some thing interesting about web. New door for the world.
A server is a program made to process requests and deliver data to clients.
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
Some thing interesting about visualization, use data art
Some thing interesting about game, make everyone happy.
We are working to build community through open source technology. NB: members must have two-factor auth.
Open source projects and samples from Microsoft.
Google โค๏ธ Open Source for everyone.
Alibaba Open Source for everyone
Data-Driven Documents codes.
China tencent open source team.