ROMS/TOMS Developers

1. I was asked to help Ken Coyle with a problem in the GOANPZ code. It turned out to be a missing chunk in inp_par.F, converting the boundary nudging timescale from units of days to sec-1. This code is duplicated for each biology model and was missing for the GOA version – could this be cleaned up so that such duplication isn’t needed?

2. I was getting odd differences in the ice depending on my domain tiling. I tracked it down to code in the ICE part of bulk_flux.F, setting the surface albedo if there is ice, leaving it random otherwise. I also took out the fice variable since it was just used the once (and not being set right either).

3. Of more general interest, I tried out the new code for interpolating from coarse forcing files. I was surprised to get stress fields such as:

The operation order is that the coarse grid data is read in and interpolated to the ROMS grid during get_data.F. It is then time-interpolated, exchanged, and finally rotated in set_data.F. It seems to me that it would be more efficient to rotate it in get_data.F since that is only happening every time you read the (six hourly) winds, not every timestep. I’m off to try to fix this.

Just a nit on the CPP option around the rotation code: CURVGRID is for when the dmde and dndx terms are nonzero and require the extra advection terms to be computed. It is not a flag for whether your grid is aligned north-south and east-west. In fact, north-south grids generally have nonzero dmde on a sphere.

Edit: I think it’s working! Just a note that if you look at my get_data.F, it is just a bit different because Don Stark asked me to change it (and change who calls it) for the CCSM coupling effort at NCAR. Most important, he asked me to put the initialization stuff elsewhere (now in get_data_init.F).

Corrected a couple of bugs in the W4DPSAS driver:

• Remove the writing of the nonlinear model initial conditions at the end of the inner loop in w4dpsas_ocean.h. The correction due to data assimilation is added in main3d.F to the first guess state. The nonlinear model is always started from the first guess. Otherwise, the solution dynamical range with accumulates in each outer loop. The final (optimal) nonlinear initial conditions are written after done with data assimilation.
• Modified def_rst.F, wrt_rst.F, tl_def_his.F, and tl_wrt_his.F to write the basic state vertical mixing coefficients when any of the vertical parameterization are activated. This will allow to control instabilities due to data assimilation. Recall that the finite amplitude tangent linear model (RPM) is used as the basic state for the next outer loop iteration in the W4DVAR algorithm. Since we are not using the tangent linear of these vertical parameterizations, the nonlinear model values are a good approximation in the Picard iterations.

In order to facilitate coupling to other models the directory tree of ROMS was changed to:

    src/
       /Compilers
       /Lib
       /Master
       /ROMS/
            /Adjoint
            /Bin
            /Drivers
            /External
            /Include
            /Modules
            /Nonlinear
            /Obsolete
            /Programs
            /Representer
            /SeaIce
            /Tangent
            /Utility

The makefile is now located in the top directory src. This allows us to used the same makefile in other components of the coupling sytem. Notice the a new directory ROMS was created with contains all previous directories except Lib and Compilers. They are now located in the top directory src so it can be used by other components of the coupling system.

A new directory Master was created and contains the main drivers of the coupled system. It currently contains the files air_ocean.h, master.F, ocean_control.F and propagator.F. These files were located previously in the Drivers directory.

All the Module.mk in each directory were modified to add the ROMS subdirectory. The file Compilers/MakeDepend was changed to conform with the new directory structure. Recall that this file is created automatically by executing a make depend from the top directory src.

The model is now compiled from the top directory src.

The user can add more directories in src for other coupling systems like that for atmospheric, wave models, etc.

Following files were modified to consolidate sediment branch code to Hernan’s latest version.

• master.F – added waves_ocean coupling and call to waves_ocean.h.
• nl_ocean.h – added waves_ocean coupling and calls to initialize coupling for atmospheric and waves models.
• Drivers/Module.mk – change CYGWIN special attribute to accomodate df and ifort.
• mod_arrays.F – added rho point source option Q_PSOURCE.
• mod_averages.F – updated sediment averages.
• mod_forces.F – added wave input parameters.
• globaldefs.h – create short notation definitions for combinations of wave data for HWAVE, LWAVE, DWAVE, PWAVE_TOP, PWAVE_BOT, and for BEDLOAD.
• mod_grid.F – added bed_thick arrays to track changes of bed elevation for morphology.
• mod_ncparam.F – added NetCDF cases for waves and radiation stress I/O.
• mod_ocean.F – added stokes velocities.
• mod_parallel.F – added WAVES_OCEAN parallel nodes variables.
• mod_param.F – added setup parameters for new test cases: INLET_TEST, VISSER, SHOREFACE, TEST_CHAN, TRENCH, SED_TOY, BEVANO. Increase unmber of diagnostic terms to include radiation stresses.
• mod_scalars.F – define diagnostics for radiation stresses. Define variables for model coupling.
• mod_sediment.F – define new variables for bed_layer_thickness, morph_fac, bed biodiffusivity.
• mod_sources.F – add Q_PSOURCE.
• mod_stepping.F – add Q_PSOURCE.
• checkdefs.F – add new applications and CPP definitions.
• def_avg.F – change 365.25 to Julian for time description.
• def_diags.F – change 365.25 to Julian for time description.
• def_floats.F – change 365.25 to Julian for time description.
• def_his.F – change 365.25 to Julian for time description, added nearshore radiation stress terms, allow h to vary for morphology, write out rmask for wet_dry, write out waves data.
• def_info.F – change 365.25 to Julian for time description, prevent writing of static h for morphologic simulations.
• def_ini.F – change 365.25 to Julian for time description.
• def_rst.F – change 365.25 to Julian for time description, read in bathy for morphology.
• get_state.F – including read for time dependent h when using sediment and morphology, include read for TKE and GLS when
• inp_par.F – updated read of model coupling input file, added write out of raditaion stress terms, GLS wave breaking criteria, wet_dry, and h for morphology, update reports for stations file.
• coupling.in – update methodology to input coupling parameters for node allocations, time steps between coupling, and reads for input file names.
• waves_coupler.F – new file added that initializes, runs, and finalizes MCT coupling.
• waves_ocean.h – updated logic to determine colors to spawn new file added.
• metrics.F – Use MAX(h) to calculate Courant number for wet_dry.
• mp_routines.F – Added USE ifport statements for CYGWIN-ifort.mk.
• set_scoord.F – Added option for log profile distribution of vertical levels.
• wrt_his.F – Add write of time dependent bathy for morphology, write out radiation stress terms, stokes velocities, wave fields, and write out wet_dry mask. When using Nearshore_Mellor would like to write out velocities in Eulerian form by subtracting Stokes velocities but this may conflict with Adjoint (needs revisiting !!)
• wrt_info.F – suppress write bathy for time dependent h with morphology.
• wrt_rst.F – Add write of time dependent bathy for morphology, Stokes velocities, and write out wet_dry mask.
• analytical.F – Added initial and forcings for various test cases.
• step2d.F – Added contributions from nearshore_mellor radiation stresses, included wet_dry.
• u2dbc_im.F – Update h to be h+zeta for shallow water (SQRT(g h+zeta)).
• v2dbc_im.F – Update h to be h+zeta for shallow water (SQRT(g h+zeta)).
• zetabc.F – Update h to be h+zeta for shallow water (SQRT(g h+zeta)).
• step_floats.F – Modify float types. Float types are 1 = neutral, 2 = constant distance below moving sea surface.
• step3d_uv.F – Include methodology for wet_dry.
• step3d_t.F – Include poitn Q sources and add vertical velowity omega to bottom face for morphology.
• ssw_bbl.h – update routine to account for nearshore velocities, set formdrag as a cpp definition.
• sg_bbl.h – update routine to account for nearshore velocities.
• mb_bbl.h – update routine to account for nearshore velocities.
• rhs3d.F – update to include radiation stresses and include omega BC for morphology.
• pre_step3d.F – update to include radiation stresses diagnostics and include omega BC for morphology.
• omega.F – include bottom boundary condition for morphology, include Q_PSOURCE.
• main2d.F – include coupler to wave model.
• main3d.F – include coupler to wave model.
• sediment.F – Update bedload routines, reform susp load.

Made few changes to allow model-observation comparison statistics after the 4DVAR data assimilation. Many thanks to Andy and Brian for their suggestions.

• Modified is4dvar_ocean.h so the nonlinear model is run at the end after data assimilation, when using IS4DVAR. The nonlinear model is initialized with the estimated initial conditions from data assimilation. The nonlinear model is then interpolated at the observation locations for comparison statistics.
• Added a new routine Utility/stats_4dvar.F to compute model-observation comparisons. The following measurements are computed for each state variable:

  1) Observations mean (obs_mean) and standard deviation (obs_std).
  2) Model mean (model_mean) and standard deviation (model_std).
  3) Model bias (model_bias).
  4) Model-Observations standard deviation error (SDE).
  5) Model-Observations cross-correlation (CC).
  6) Model-Observations mean squared error (MSE).

  All these statistical quantities are written into 4DVAR output NetCDF file (MODname) and standard output file.

• Added stats_4dvar calls to all 4DVAR data assimilation drivers: is4dvar_ocean.h, s4dvar_ocean.h, w4dpsas_ocean.h, and w4dvar_ocean.h.
• Fixed an out-range index bound in arrayd A2d and A3d in the exact background error normalization coefficients routine normalization.F. This is a benign bug.
• Made changes to inp_par.F to insure that there is only one record in the restart file after the nonlinear model is run at the end of an assimilation cycle. The parameter nRST is set to ntimes. The restart file can be used as the first guess for the next assimilation file in sequential data assimilation.