Custom Query (986 matches)

I have been working on a major overhaul of Matlab scripts and tools in the ROMS repository for several weeks. Since all the data in ROMS is managed with NetCDF files, some of these scripts use a NetCDF interface to Matlab to process the data. There are several interfaces for Matlab available from third parties. The most widely used interfaces are MEXNC and SNCTOOLS.

However, starting with Matlab version 2012a, released on Feb 9, 2012, the native interface to NetCDF is the preferred method for processing NetCDF data in the scripts distributed in the ROMS repository matlab and described here. The native interface was introduced in Matlab version 2008b for NetCDF-3 type files. The NetCDF-4 support started in version 2010b. The support for HDF5 files was completed in version 2011a. The OpenDAP support began in version 2012a. If your Matlab version is older than 2008b, we highly recommend that you update to the newest version. However, in the basic generic scripts we have switches for older versions to activate either the MEXNC interface for standard NetCDF files and the SNCTOOLS interface to process NetCDF files on an OpenDAP server.

Nesting Scripts:

This update also includes several new scripts to process ROMS nesting contact regions and contact points. The contact points are processed outside of ROMS and all the connectivity within nested grids is read from an input NetCDF file. This facilitates the configuration of various grid Sub-Classes. It tremendously simplifies the processing of such points in parallel computations. It also gives the user full editing control of the contact points. The contact points need to be set-up once for a particular application.

The information about these new scripts is very technical and difficult to explain with plain words. We need definitions, equations, and a lot of diagrams to explain and illustrate the complexity of ROMS nesting capabilities. Therefore, we decided to put all that information in ​WikiROMS instead of a very long message in this forum. In particular, we added two new menus to the documentation portal:

​Nested Grids

​Matlab Scripts

This update includes a full revision of ROMS shared-memory pragma directives using OpenMP standard. This is a very important and delicate update that requires expertise. Luckly, I doubth that will affect your customized code.

All the parallel loops of ROMS are modified to simpler directives. For example, the old strategy:

!$OMP PARALLEL DO PRIVATE(thread,subs,tile) SHARED(numthreads)
            DO thread=0,numthreads-1
              subs=NtileX(ng)*NtileE(ng)/numthreads
              DO tile=subs*thread,subs*(thread+1)-1,+1
                ...
              END DO
            END DO
!$OMP END PARALLEL DO

is replaced with:

            DO tile=first_tile(ng),last_tile(ng),+1
              ...
            END DO
!$OMP BARRIER

In shared-memory, the parallel threads are spawn in the higher calling routines (drivers). For example, we now have:

!$OMP PARALLEL
      CALL main3d (RunInterval)
!$OMP END PARALLEL

This directive is less restrictive and allows MASTER, BARRIER, and other useful OpenMP pragmas inside the parallel regions. If you are interested, please see the following discussion in the ​Forum.

This change cleans the code and facilitates parallelization of tricky algorithms for nesting, MPDATA, random number generation, point-sources, etc using the shared-memory paradigm.

WARNINGS:

• The values of NtileX(ng) and NtileE(ng) are no longer equal to one in distributed-memory (MPI). They have the same values as the ones specified in standard input: NtileI(ng) and NtileJ(ng). Notice that in the critical regions for global reduction operations, we now use the following code instead:

  #ifdef DISTRIBUTE
        NSUB=1                             ! distributed-memory
  #else
        IF (DOMAIN(ng)%SouthWest_Corner(tile).and.                        &
       &    DOMAIN(ng)%NorthEast_Corner(tile)) THEN
          NSUB=1                           ! non-tiled application
        ELSE
          NSUB=NtileX(ng)*NtileE(ng)       ! tiled application
        END IF
  #endif
  

  That is, we do a special exception for distribute-memory cases. This change is necessary in your customized versions of ana_grid.h and ana_psource.h.

• Notice that we no longer use the TILE (uppercase) as argument to the kernel routines. We use tile (lowercase) instead. This was important in previous versions of the distributed-memory code where TILE was replaced with MyRank during C-preprocessing. Be careful with this one...

• Notice that few important variables in mod_scalars.F and mod_stepping.F use the THREADPRIVATE directive in shared-memory applications so all the parallel threads have a private copy of such variables to avoid parallel collisions.

• Two new variables (first_tile(ng) and last_tile(ng)) are introduced to specify the tile range in each parallel region:

        integer, allocatable :: first_tile(:)
        integer, allocatable :: last_tile(:)
  
  !$OMP THREADPRIVATE (first_tile, last_tile)
  

  These variables are set during the initialization of ROMS kernel using:

  !$OMP PARALLEL
  #if defined _OPENMP
        MyThread=my_threadnum()
  #elif defined DISTRIBUTE
        MyThread=MyRank
  #else
        MyThread=0
  #endif
        DO ng=1,Ngrids
          chunk_size=(NtileX(ng)*NtileE(ng)+numthreads-1)/numthreads
          first_tile(ng)=MyThread*chunk_size
          last_tile (ng)=first_tile(ng)+chunk_size-1
        END DO
  !$OMP END PARALLEL
  

Many thanks to Sasha shchepetkin for ​suggesting this strategy. Also many thanks to Mark Hadfield for his persistence and testing.

A new option FLOAT_OYSTER has been added to model the swimming, sinking, and settlement of planktonic oyster larvae. This development includes the infrastructure to add different behaviors to the Lagrangian drifters (FLOATS). The design is similar to that of the ecosystem models in ROMS and includes several header files, as described below. This option is under the umbrella of internal C-preprocessing flag FLOAT_BIOLOGY, which is defined in globaldefs.h.

The biological behavior for the floats is achieved by interpolating additional variables using interp_floats in step_floats.F and then calling biology_floats in the predictor and corrector steps:

# ifdef FLOAT_BIOLOGY
!
!-----------------------------------------------------------------------
!  Compute biological float behavior, corrector step.
!-----------------------------------------------------------------------
!
      CALL biology_floats (ng, Lstr, Lend, .FALSE., my_thread)
# endif

The planktonic behavior is based on the model of Dekshenieks et al. (1993; 1996; 1997). It calculates the size (length) and development of oyster larvae. This model was implement in ROMS with the help of Diego Narvaez as part of his Ph.D. dissertation. Many thanks to John Klinck for his suggestions about the numerical implementation using the predictor/corrector scheme of the floats. Results about this oyster model can be found in Narvaez et al. 2012 a,b (J. Mar. Res., in review).

This model includes the following files:

• ROMS/Modules/mod_behavior.F: module declaring all the biological behavior model internal parameters.

        MODULE mod_behavior
  
  #if defined FLOATS && defined FLOAT_BIOLOGY
  # if defined FLOAT_OYSTER
  #  include <oyster_floats_mod.h>
  # endif
  #endif
  
        END MODULE mod_behavior
  

• ROMS/Nonlinear/Biology/biology_floats.F: Lagrangian drifters behavior model driver.

        MODULE biology_floats_mod
  #if defined NONLINEAR && defined FLOATS && defined FLOAT_BIOLOGY
  
        implicit none
  
        PRIVATE
        PUBLIC  :: biology_floats
  
        CONTAINS
  
  # ifdef FLOAT_OYSTER
  #  include <oyster_floats.h>
  # endif
  
  #endif
  
        END MODULE biology_floats_mod
  

• ROMS/Nonlinear/Biology/oyster_floats.h: planktonic oyster larvae model using Lagrangian particles.

• ROMS/Nonlinear/Biology/oyster_floats_mod.h: oyster model declaration and allocation of internal parameters.

• ROMS/Nonlinear/Biology/oyster_floats_inp.h: code to read input model parameters which is included read_fltbiopar.F.

• ROMS/Nonlinear/Biology/oyster_floats_def.h: code to define input oyster model parameters in all output NetCDF files which is included in ROMS/Utility/def_info.F.

• ROMS/Nonlinear/Biology/oyster_floats_wrt.h: code to write out input model parameters in all output NetCDF files which is included in ROMS/Utility/wrt_info.F.

• ROMS/External/behavior_oyster.in: oyster model standard input parameters.

• ROMS/Utility/read_fltbiopar.F: routine to read floats biological behavior parameters.

  #include "cppdefs.h"
  #if defined FLOATS && defined FLOAT_BIOLOGY
  
  # if defined FLOAT_OYSTER
  #  include <oyster_floats_inp.h>
  # endif
  #else
        SUBROUTINE read_FltBioPar
        END SUBROUTINE read_FltBioPar
  #endif
  

NOTES:

• Mostly all the floats parameter are moved from mod_scalars.F to mod_floats.F. I think that it is better to have everything in a single module.

• All the header files (*.h) above are included within <...> to allow the user to customize any of them in the project directory while keeping the distributed code intact (check the build script for details).

• All the floats input scripts (floats_*.in) have two new keyword entries: A switch Fprint to control the printing of floats initial positions to standard out file. It can be used to turn off the printing of information when thousands of floats are released. This information is still in output floats NetCDF file.

        Fprint == T
  

  And the float's behavior model input script file name:

       FBIONAM =  behavior_oyster.in
  

• Currently, the oyster model assumes constant food (mg Carbon/l) and turbidity (g/l). This can be changed in the future with spatial and temporal variability for either food and turbidity. If this is the case, we need to have something like:

            IF (Predictor) THEN
              Lfood=track(ifood,nf,l)
              Lturb=track(iturb,nf,l)
              ...
            ELSE
              Lfood=track(ifood,nfp1,l)
              Lturb=track(iturb,nfp1,l)
              ...
            END IF
  

  in oyster_floats.h. The food variability may be from data, ecosystem model, or analytical functions. Similarly, the turbidity may be from data, sediment model, or analytical functions.

• A new script ROMS/External/floats_upwelling.in was added to test the oyster model in the UPWELLING example.

• Notice that behavior_oyster.in includes three look tables to: (i) compute larvae swimming speed as a function of larvae size and temperature, (ii) larvae growth rate factor as function of salinity and temperature, and (iii) larvae growth rate as function of food supply and larvae size. I recommend extreme caution when editing this part of the input script. Otherwise, ROMS may read these value incorrectly!!!

This capability illustrates how the biological behavior kernel can be implemented in ROMS floats. It is up to the user to validate and test their own customized behavior models.

Good luck...