I think the CVMix user guide will grow out of this document, but for now it just seems like an easy way to share with MOM and MPAS users what I've done in POP to make calls into CVMix.

Background

POP has a module named vertical_mix that is used to call either constant mixing, Richardson mixing, or KPP mixing. These three methods have their own modules: vmix_const, vmix_rich, and vmix_kpp. (It looks like the tidal mixing is in tidal_mixing.F90, which is used by KPP, but I haven't gotten that far yet.) Routines from vertical_mix are called from the baroclinic module, which is the main driver for computing baroclinic velocities and tracer fields. All of the changes I describe below were made in vertical_mix.F90, and I suspect that MPAS and MOM will want to make changes in a similar module.

Constant / Convective Mixing

The first goal for CVMix was to change POP to call CVMix's vmix_background and vmix_convective modules rather than POP's vmix_const module. This was a pretty straight forward process, requiring three steps.

Use Statements

Constant mixing requires three packages from CVMix:

use vmix_kinds_and_types
use vmix_background
use vmix_convection

The first package is needed to set up all the derived data types:

type (vmix_data_type),  allocatable, dimension(:) :: Vmix_vars
type (vmix_global_params_type)                    :: Vmix_params
type (vmix_bkgnd_params_type)                     :: Vmix_bkgnd_params
type (vmix_conv_params_type)                      :: Vmix_conv_params
Initialization

In POP, if you want to use constant and convective mixing, the vertical_mix module reads in various convection-related namelist variables, allocates memory for VDC and VVC (what POP calls the vertical diffusivity and viscosity coefficients, respectively) and then calls init_vmix_const() in the vmix_const module. That routine reads in POP's &vmix_const namelist and set the arrays VDC and VVC to the scalar constants const_vdc and const_vvc, respectively.

allocate (VDC(nx_block,ny_block,km,1,nblocks_clinic),  &
          VVC(nx_block,ny_block,km,  nblocks_clinic))

Note that nblocks_clinic is the number of threads, nx_block and ny_block are the number of longitude / latitude points in each thread, respectively, and km is the number of vertical levels. So there are nx_block*ny_block*nblocks_clinic columns, each with km levels.

Obviously CVMix will not read in a POP namelist, so I moved the input into the vertical_mix module. I also needed to initialize several CVMix variables:

ncol = nx_block*ny_block*nblocks_clinic
allocate(Vmix_vars(ncol))
do bid=1,nblocks_clinic      
   do inx=1,nx_block
      do jny=1,ny_block
         colid = (bid-1)*nx_block*ny_block+(jny-1)*nx_block+inx
         Vmix_vars(colid)%nlev = KMT(inx, jny, bid)
      end do
   end do
end do
Prandtl = 10.0_r8
call vmix_init(Vmix_vars, Vmix_params, ncol, km, Prandtl)

Vmix_vars(:) is the main CVMix data type. It is an array over each column -- as mentioned before, in POP the number of columns is the product of three of the dimensions of VDC and VVC --, colid is a POP-specific renumbering scheme to use natural ordering to get the (nx,ny,nblock) index down to a single number. Besides allocating memory for Vmix_vars, this routine also assigns a value to each Vmix_vars(:)%nlev, which is the number of ocean levels that particular column has. (In POP, this information is stored in the array KMT.) The last step is to call vmix_init, which is currently in the vmix_kinds_and_types module (on Coffee Talk, Linda Richman would probably point out that "vmix_init is neither a kind nor a type, discuss"). This is where memory for each component of Vmix_vars(:) gets allocated and a couple global parameters (maximum number of levels, which is needed for allocating some memory in the Vmix_bkgnd_params_type, and the Prandtl number, which Steve has pointed out should really be tied to a mixing scheme rather than saved as a 'global' parameter).

Once Vmix_vars is set up, we are ready to set up the background and convective mixing parameters / initialize variable values:

call vmix_init_bkgnd(Vmix_bkgnd_params, const_vvc, const_vdc)
if (convection_itype.eq.convect_type_diff) then
   call vmix_init_conv(Vmix_conv_params, convect_diff, convect_visc)
end if

do bid=1,nblocks_clinic      
  do inx=1,nx_block
    do jny=1,ny_block
      VVC(inx,jny,:,  bid) = Vmix_bkgnd_params%static_visc(1,1)
      VDC(inx,jny,:,1,bid) = Vmix_bkgnd_params%static_diff(1,1)
    end do
  end do
end do

vmix_init_bkgnd allocates memory for Vmix_bkgnd_params%static_visc and %static_diff; because POP is sending scalar constants, these will both be 1x1 arrays. If you send vmix_init_bkgnd() 2D arrays (and the number of columns), it will allocate more memory for those two components. Similarly, vmix_init_conv sets the value for the convective mixing in Vmix_conv_params. The second block then populates the POP mixing coefficients with the values set in the CVMix routine.

Updating the mixing coefficients

At each timestep, the call to update the mixing coefficients is straight-forward (note that in POP this is called inside a threaded region and bid is the block number or thread ID):

! Pack into vmix_kind_type
do k=1,km
  kp1 = min(k+1,km)
  call state(k  ,kp1,TMIX(:,:,k  ,1), TMIX(:,:,k  ,2), &
             this_block, RHOOUT=RHOK_MAT(:,:,k,bid) )
  call state(kp1,kp1,TMIX(:,:,kp1,1), TMIX(:,:,kp1,2), &
             this_block, RHOOUT=RHOKP_MAT(:,:,k,bid) )
end do
do inx=1,nx_block
  do jny=1,ny_block
    colid = (bid-1)*nx_block*ny_block+(jny-1)*nx_block+inx
    nlev = Vmix_vars(colid)%nlev
    Vmix_vars(colid)%dens(:)      = RHOK_MAT (inx,jny,1:nlev,bid)
    Vmix_vars(colid)%dens_lwr(:)  = RHOKP_MAT(inx,jny,1:nlev,bid)

    Vmix_vars(colid)%visc_iface(2:nlev+1)   = VVC(inx,jny,1:nlev,  bid)
    Vmix_vars(colid)%diff_iface(2:nlev+1,1) = VDC(inx,jny,1:nlev,1,bid)
  end do
end do

call vmix_coeffs_bkgnd(Vmix_vars, Vmix_bkgnd_params, &
                       (bid-1)*nx_block*ny_block+1,  &
                       bid*nx_block*ny_block)

if (convection_itype.eq.convect_type_diff) then
  call vmix_coeffs_conv(Vmix_vars, Vmix_conv_params, &
                        (bid-1)*nx_block*ny_block+1, &
                        bid*nx_block*ny_block)
end if

! Unpack from vmix_kind_type
do inx=1,nx_block
  do jny=1,ny_block
    colid = (bid-1)*nx_block*ny_block+(jny-1)*nx_block+inx
    nlev = Vmix_vars(colid)%nlev
    VVC(inx,jny,1:nlev,  bid) = Vmix_vars(colid)%visc_iface(2:nlev+1)
    VDC(inx,jny,1:nlev,1,bid) = Vmix_vars(colid)%diff_iface(2:nlev+1,1)
  end do
end do

The two calls to state() are made to get the density at level k and k+1 (note that this takes two calls instead of just using RHOK_MAT(:,:,k,:) and RHOK_MAT(:,:,k+1,:) because the density at level k is calculated for a parcel that is adiabatically displaced down to level k+1). These densities are needed to determine where to use the convective mixing.

The calls to vmix_coeffs_bkgnd() and vmix_coeffs_conv() actually update the mixing coefficients in Vmix_vars(:), and then the data needs to be unpacked back into the POP variables.

  • No labels

3 Comments

  1. Todd Ringler

    Hi All,

    I have a lot of comments here. But I would rather comment on the driver / vmix_module case than a POP / vmix_module case. How soon could I get access to driver that sets up a single column, calls the init vmix procedure (if there is one) and then calls something to compute diffusivities/viscosities. Can we get the driver / vmix module into a clean svn repo?

    The one comment I will offer up here is that I propose there is a single interface to "put" parameters. The code above suggests that there is a separate call (e.g. vmix_init_background and vmix_init_conv) for each type of vertical mixing. This will get very cluttered very quickly. An alternative is have a single put interface, e.g.

    vmix_put ( text string, data to put )

    The idea here is that ALL vmix_params are stored inside the cvmix code. Each parameter has a default value. Each parameter is associated with a string. The ocean model (POP, MOM, MPAS) alters these parameters (at init or at any other time) by calling the vmix_put. The calling code can also retrieve parameter setting via the vmix_get (text string, parameter to get). 

    I would envision the "put" interface would also turn on/off the various choices that exist with vmix. In that, the default setting for all choices within vmix is false and those that are sought are turned on via a call to vmix_put( text string, data to put) where in this case the "data to put" is a logical turn and the text string refers to one of the choices of vertical mixing.

    In this fashion, there is not need whatsoever for the calling model to have to know anything about any of the following:use vmix_kinds_and_typesuse vmix_backgrounduse vmix_convection
    This is all stored inside the cvmix.

    Cheers, Todd

    1. Michael Levy

      Hi Todd,

      I'm hoping to have the stand-alone driver done before Steve leaves this coming Wednesday (Aug 15). I'll look into the separate repository options, but that might take a little more time to set up. I like your "put" / "get" option, I was initially planning to do something similar but with the pointer / mem copy issue still unsettled I figured I'd throw something clunky together for testing purposes.

      It sounds like you want to replace the Vmix_params_type and Vmix_*_params_type data types with global parameters? As in, have the variables vmix_convect_diff and vmix_convect_visc instead of Vmix_conv_params%convect_diff and Vmix_conv_params%convect_visc? I set it up using the latter scheme because Steve had specifically asked us to pass time-independent variables and parameters rather than using "use" statements. I'm open to switching, but I'd like Steve to chime in first -- either to argue his side or to point out that I was misunderstanding him.

      I also want to clarify what you mean by "the put interface would also turn on/off the various choices that exist with vmix." It makes sense from an initialization point of view (rather than call vmix_init_bkgnd and vmix_init_conv, just use the put command to turn on these options), but what do you see for the routines that actually return the mixing coefficients? Steve hints at an order in his CVMix pdf (background, shear, tidal, double diffusion, kpp, convective) -- do we keep this order but skip over any that are not initialized? I was under the impression that we wanted to call each method individually in case the model needs to run some smoothing in between steps.

      ~Mike

      1. Todd Ringler

        Hi Mike,

        This is off the top of my head and would very much like the opportunity for all of us to think it through it, but it might looks something like

        vmix_params % type_of_scheme % attributes

        Unknown macro: {%value, %string}

        so for example, it would be vmix_params % convection % viscosity % character (where the character string might be "convective_viscosity")

        and then vmix_params % convection % viscosity % real (where the real might be set to 1.0 double )

        and there would also be vmix_params % convection % status, where status is a logical with .true. for on and .false. for off

        the convective viscosity is set by "call vmix_put("convection", "convective_viscosity", 1.0)"

        the put command matches the string to the derived data type attributes defined for convection. when a match is found it 1) makes sure that the types match (i.e the attribute is a real and i passed a real) and 2) sets the vmix_params % convection % viscosity % value to 1.0

        with this approach, the calling model never needs to know anything about vmix_params. the model just has to know what parameters it wants to set. the interface to set parameters is completely uniform across the entire vmix system. by construction, this allows "on the fly" reconfiguration of what vmix does.

        i personally would prefer one call into vmix that follows the sequence outlined in steve's notes. i think that this is easily accomplished with a get / put system and logical on/offs. i don't want to handle diffusivities for each process. i want to determine what my mixing scheme is (via puts), call vmix for each column and get out diffusivities/viscosities, then do a solve and move on. if i have to make separate calls into each type of vmix, then my "driver" routine will be littered with if tests. i think in the end we will end up with better, more extensible code if we push that if testing into the top level of vmix.

        while i think i have the concept correct, i am sure that doug and phil can improve on the design and implementation.

        cheers,

        todd.