Home

Technical Manual for a Coupled Sea-Ice/Ocean Circulation

1. 2 3 CPP CPPFLAGS MY_CPP_FLAGS lt gt 0 CLEAN endef compile rules define compile rules foreach f local_src call one compile rule call source to object f call 90 source f f endef 1 We define a function to convert the path from the source directory to the Build directory called source dir to binary dir Note that the Build directory is called SCRATCH_ DIR here All it does is strip off the leading directory with the the built in function notdir then paste on the Build directory 2 Next comes source to object which calls the function above to return the object filename when given the source filename It assumes that all sources have a F extension 141 3 A similar function is f90 source which returns the name of the intermediate source which is created by cpp from our F file 4 The Module mk fragment in each library source directory invokes make library which takes the library name and the list of sources as its arguments The function adds its library to the global list of libraries and provides rules for building itself The double dollar signs are to delay the variable substitution Note that we call source dir to binary dir instead of source to ob ject this is a work around for a make bug 5 The next one compile rule takes three arguments the o filename the f90 filename and the F filename From these it produces the make rules for running cpp a
2. t Gm There is a related set of secondary weights bm used as a Sa bu In order to maintain constancy preservation this relation must hold OT Ory mn jAt Cn zy il zy o y 27 2 gt 37 Shchepetkin and McWilliams 73 introduce a range of possible weights but the ones used here have a shape function qye q T AA as where p q are parameters and Ap To and r are chosen to satisfy normalization consistency and second order accuracy conditions i T A T dr 1 n 0 1 2 39 0 using Newton iterations 7 is the upper limit of 7 with A 7 gt 0 In practice we initially set p 2 p q 2 Ap 1 r 0 and r p D p 4 1 compute A T using eq 38 normalize using M M i Y am l Y any 1 40 m 1 m 1 and adjust r iteratively to satisfy the n 2 condition of 39 We are using values of p 2 q 4 and r 0 284 This form allows some negative weights for small m allowing M to be less than 1 5M ROMS also supports an older cosine weighting option which isn t recommended since it is only first order accurate 4 6 Density in the mode coupling Equation contains the term R computed as the difference between the 3 D right hand side and the 2 D right hand side The pressure gradient therefore has the form gD OC gD o mo n F 41 where the term in square brackets is the mode coupling term and is held fixed over all the barotr
3. The C preprocessor will also delete C language comments starting with and ending with as in Hendif MASKING When mixed with Fortran code it is safer to use a Fortran comment F 5 A note on style In ROMS the style adopted is as shown above and here define NEP5 Hif defined NEP5 defined BERING ifdef MASKING endif Hendif The gnu community has instead adopted this style define NEP5 1 if NEPS BERING if MASKING endif Hendif Both styles work and you should be aware that there s more than one way to do it F 6 Potential problems The use of the C preprocessor is not entirely free of problems but many can be worked around or avoided by using the gnu version of cpp with the traditional flag 1 Apostrophes in Fortran comments cpp does not know that it is in a comment and some versions will complain about unmatched apostrophes in the following c Some useful comment about Green s functions 131 2 C comments Some of the newer versions of cpp will remove C comments which go from to the end of the line Some perfectly reasonable Fortran lines contain two consecutive slashes such as common vari var2 44 format and Fortran 90 string concatenation mystring Hello World 3 Macro replacement One feature of cpp is that you can define macros and have it perform replacements The code define REAL double precision REAL really_long_vari
4. defined BIO_GOANPZ defined PASSIVE_TRACERS undef BIOFLUX sum Nitrogen fluxes between boxes define ANA_BIOLOGY analytical biology initial conditions define ANA_BPFLUX analytical bottom passive tracers fluxes define ANA_SPFLUX analytical surface passive tracers fluxes define DIAPAUSE Enable Neocalanus seasonal vertical migration undef FLOAT_VWALK define IRON_LIMIT Add iron as passive 11th tracer undef TCLM_NUDGING Nudging of tracer climatology for iron endif if defined NEMURO undef ANA_BIOLOGY analytical biology initial conditions define ANA_BPFLUX analytical bottom passive tracers fluxes define ANA_SPFLUX analytical surface passive tracers fluxes define IRON_LIMIT Add iron as passive 11th tracer 109 define IRON_RELAX undef IRON_RSIN define BIO_SEDIMENT define HOLLING_GRAZING undef IVLEV_EXPLICIT undef ANA_BIOSWRAD undef DIAGNOSTICS_BIO undef BIO_SEDIMENT endif HE HH HH H OH OF Here we have declared that we want two open sides and two closed walls For each open side we re using a Chapman Flather combination on the 2D flow and a radiation nudging combination on the 3D flow We ve got masking salinity and the non linear equation of state We want Lapla cian viscosity on o surfaces diffusion along constant z surfaces and the full non linear curvilinear momentum equations Sea ice is turned on as is the NEMURO ecosystem
5. oooooocooococcommoomo 1 923 0 3453 Model 2D kernel fae ww eai NARA 365 560 65 6278 2D 3D coupling vertical metrics 1 950 0 3500 Omega vertical velocity oooooocmoommoomo eee 3 135 0 5628 Equation of state for seawater 04 1 645 0 2954 3D equations right side terms 14 268 2 5615 3D equations predictor step 00 00 31 518 5 6583 Pressure gradient oi dicc ienris iiaii eee 10 157 1 8235 Harmonic mixing of tracers S surfaces 4 130 0 7414 Harmonic stress tensor S surfaces 8 211 1 4741 Corrector time step for 3D momentum 27 976 5 0225 Corrector time step for tracers 19 558 3 5112 Total 548 216 98 4193 All percentages are with respect to total time 557 021 ROMS TOMS Output NetCDF summary for Grid 01 number of time records written in HISTORY file 00000021 98 number of time records written in RESTART file 00000002 number of time records written in AVERAGE file 00000020 Analytical header files used ROMS Functionals ana_btflux h ROMS Functionals ana_grid h ROMS Functionals ana_initial h ROMS Functionals ana_smflux h ROMS Functionals ana_stflux h ROMS Functionals ana_vmix h ROMS TOMS DONE Friday October 30 2009 9 20 38 AM Note that it ends by printing out a profile of where the time was used the ratios will vary depending
6. 1 Cou T mS 228 The ocean is checked at each depth k and at each timestep for supercooling If the water is below freezing the temperature and salinity are adjusted as in equations and and the ice above is thickened by the amount Ah mazos 229 2 5 2 3 Differences from Mellor and Kantha We have tried to modify the hakkis model to more closely follow MK89 However there are also ways in which we have deviated from it e Add advection of snow e Add lateral melting of snow when ice is melting laterally We iterate on the solution of T3 We added various limiters Ice concentration Amin lt A lt 1 0 Amin 0 0 Ice thickness h gt 0 0 Snow thickness hs gt 0 0 l Brine fraction r lt fmax Tmax 0 2 Surface water 0 0 lt Ws lt Ws max Ws max 0 1 45 6 Details of the Code 6 1 Directory structure The directory structure is as shown in Fig 9 with the ability to run the ocean alone or coupled to atmospheric and or wave models If running just the ocean the model can be run forward in time the nonlinear model or as an adjoint tangent linear or representer model for data assimilation purposes This document describes the uncoupled forward model only specifically the version used for the Northeast Pacific domain containing sea ice and other changes from the main trunk code Details are subject to change without notice check your own source code for specific de
7. 47 48 49 50 51 Y Kanarska A F Shchepetkin and J C McWilliams Algorithm for non hydrostatic dy namics in roms Ocean Modeling 18 143 174 2007 R F Keeling B B Stephens R G Najjar S C Doney D Archer and M Heimann Seasonal variations in the atmospherice 02 n2 ratio in relation to kinetics of air sea gas exchange Global Biogeochem Cycles 12 141 163 1998 B W Kernighan and D M Ritchie The C Programming Language Prentice Hall Englewood Cliffs New Jersey 07632 second edition 1988 M J Kishi M Kashiwai D M Ware B A Megrey D L Eslinger F E Werner M Noguchi Aita T Azumaya M Fujii S Hashimoto D Huang H lizumi Y Ishida S Kang G A Kantakov H C Kim K Komatsu V V Navrotsky S Lan Smith K Tadokoro A Tsuda O Yamamura Y Yamanaka K Yokouchi N Yoshie J Zhang Y I Zuenko and V I Zvalin sky Nemuro a lower trophic level model for the north pacific marine ecosystem Ecological Modelling 2002 12 25 2007 H Lamb Aydrodynamcis Cambridge University Press 6 edition 1932 1945 Dover Publica tions reproduction W G Large Modeling and parameterization ocean planetary boundary layers In E P Chas signet and J Verron editors Ocean Modeling and Parameterization pages 81 120 Kluwer Academic Publishers 1998 W G Large J C McWilliams and S C Doney Oceanic vertical mixing a review and a model with a nonlocal boundary laye
8. FFLAGS TARGET_ARCH c Printing the rules database will show variables that make is picking up from the environment from the Makefile and from its built in rules and which of these sources is providing each value G 2 2 Assignments In the old days I only used one kind of assignment equals sign Gmake has several kinds of assignment other makes might as well but I no longer have to know or care An example of the power of gnu make is shown by an example from my Cray X1 Makefile There is a routine which runs much more quickly if a short function in another file is inlined The way to accomplish this is through the O inlinefrom file directive to the compiler I can t add this option to FFLAGS since the inlined routine won t compile with this directive there is only the one file that needs it I had FFLAGS 0 3 aggress e I e m FFLAGS2 0 3 aggress 0 inlinefrom lmd_wscale f90 e I e m lmd_skpp o lt TAB gt FC c FFLAGS2 90 The first change I can make to this using other assignments is FFLAGS 0 3 aggress e I e m FFLAGS2 FFLAGS 0 inlinefrom lmd_wscale f90 The assignment means to evaluate the right hand side immediately In this case there is no reason not to as long as the second assigment follows the first one since it s using the value of FFLAGS For the plain equals make doesn t evaluate the right hand side until its second pass through the Makefile However gnu ma
9. include COMPILERS 0S strip FORT mk We pick one include file at compile time here picking Linux ftn mk containing the Cray cross compiler information The value Linux comes from the uname command the ftn comes from the user and the two are concatenated The sed command will turn the slash in UNICOS mp into a dash the native Cray include file is UNICOS mp ftn mk Strip is a built in function to strip away any extra white space Another tricky system is CYGWIN which puts a version number in the uname output such as CYGWIN_NT 5 1 Gnu make has quite a few built in functions one of which allows us to do string substitution OS patsubst CYGWIN_ CYGWIN 0S In make the symbol is a sort of wild card much like in the shell Here we match CYGWIN followed by an underscore and anything else replacing the whole with simply CYGWIN Another example of a built in function is the substitution we do in objects subst F 0 sources where we build our list of objects from the list of sources There are quite a few other functions plus the user can define their own From the book 52 137 GNU make supports both built in and user defined functions A function invocation looks much like a variable reference but includes one or more parameters separated by commas Most built in functions expand to some value that is then assigned to a variable or passed to a subshell A user defined function is stored in a variable
10. Ah 2074 184 where the nonlinear viscosities are given by P c 2 2 1 2 ie 2 1 1 e 4e7 e 2e11 22 1 1 e n 3 186 We would also like to have an explicit model that can be solved efficiently on parallel computers The EVP rheology has a tunable coefficient E the Young s modulus which can be chosen to make the elastic term small compared to the other terms We rearrange the VP rheology E Ps ij i 187 m Anc _ 7 O kk055 E 4c bij ij then add the elastic term 1 00 n C PS EB a To ou inc OrkOig x big sj 188 37 Much like the ocean model the ice model has a split timestep The internal ice stress term is updated on a shorter timestep so as to allow the elastic wave velocity to be resolved Once the new ice velocities are computed the ice tracers can be advected using the MPDATA scheme 77 The tracers in this case are the ice thickness ice concentration snow thickness internal ice temperature and surface melt ponds The continuity equations describing the evolution of these parameters equations 189 191 also include thermodynamic terms Sp Ss and 4 which will be described in DAh O uAh O vAh a gy FSA RD 189 OAhs O UAhs vAhs a jy Sete 190 OA O uA O vA lt AS lt l 191 The first two equations represent the conservation of ice and snow Equation is discussed in some detail in MK89 but represents the advect
11. ImaxS JminS JmaxS RETURN END SUBROUTINE set_data Here there are two sets of array lower and upper bounds those in the LBi family and those in the IminS family Both depend on the Istr family shown in Fig The IminS family is for work arrays that are local to an MPI process or to an OpenMP thread also local to a _ tile routine They are initialized IminS BOUNDS ng Istr tile 3 ImaxS BOUNDS ng Iend tile 3 JminS BOUNDS ng Jstr tile 3 JmaxS BOUNDS ng Jend tile 3 and used real r8 dimension IminS ImaxS JminS JmaxS worki real r8 dimension IminS ImaxS JminS JmaxS work2 The Istr and LBi families are dimensioned by the number of tiles once it is known by inp_ par DO ng 1 Ngrids Ntiles Ntilel ng NtileJ ng 1 allocate BOUNDS ng LBi 1 Ntiles allocate BOUNDS ng Jend 1 Ntiles END DO 70 E A E E a oo P malo Figure 16 A tiled grid with out of date halo regions shown in grey and the interior points color coded by tile a before an exchange and b after an exchange 71 They are then initialized in calls to the routines in get__bounds F If a tile is on the western edge Itile 0 then UBi is set to LOWER_BOUND_ I Imin LOWER_BOUND_1 IF Itile eq 1 or Itile eq 0 THEN LBi Imin ELSE LBi Istr Nghost END IF Tracking the origins of LOWER_BOUND_I we find it in globaldefs h ifdef EW_PERIODIC define LOWER_BOUND_I GHOST_POINTS else
12. Include as a template for your application you must rename it If you don t change the name ROMS will use the one in ROMS Tnclude and your file will be ignored during the build procedure 7 1 3 Functionals Some of the cpp Options have names beginning with ANA __ For each one of these you will be expected to provide an analytic expression for the field in question in the corresponding include file These files are listed in 6 5 and their location is determined by MY_ANALYTICAL_ DIR You may chose to copy those from User Functionals to some new directory and place your version of the assignments within ifdef WIKI_TEST Set weird and wonderful winds Hendif This makes it easy to search for later if nothing else 7 1 4 checkdefs F For each new cpp variable other than your application name it is recommended that you also add the appropriate code to checkdefs F such as ifdef SLEET IF Master WRITE stdout 20 SLEET amp amp Sleet falling on the ice option is lenstr Coptions 1 Coptions is is 7 SLEET endif ICE Note that the number 7 on the Coptions line must be set according to the length of the string you are adding In this case 7 is for SLEET including the comma and the space Again you do not need to do this for your application name WIKI_ TEST here since checkdefs will print whatever is in My AppCPP 7 1 5 Model domain One of the first things the user must decide is
13. Though this doesn t seem to work on the Mac 3 The very first makefile I showed had a set list of files to remove on make clean We now build a list called clean_ list clean_list core ipo SCRATCH_DIR ifeq strip SCRATCH_DIR clean_list core o 00 mod f90 lib a bak clean_list CURDIR ipo endif In other words we want to clean up the Build directory unless it happens to be the top level directory in which case we only want to remove specific files there 4 all is the first target that gets seen by make making it the default target In this case we know there is only the one binary whose name we know the book 52 shows what to do with more than one binary Both all and clean are phony targets in that no files of those names get generated make has the PHONY designation for such targets Also the clean target doesn t require any compiler information so the compiler include doesn t happen if the target is clean ifneq MAKECMDGOALS clean include COMPILERS 0S strip FORT mk endif MAKECMDGOALS is a special variable containing the current make target 5 We ll be creating different executable names depending on which options we ve picked BIN BINDIR oceanS ifdef USE_DEBUG BIN BINDIR oceanG else 144 ifdef USE_MPI BIN BINDIR oceanM endif ifdef USE_OpenMP BIN BINDIR oceanO endif endif The NetCDF l
14. define LOWER_BOUND_I O endif In the case of set__data we are simply passing array indices for the tiled arrays To access the tiled arrays from within set_data_tile we need to use the relevant modules and then refer to the array with its full name USE mod_forces CALL set_2dfld_tile ng tile iNLM idCfra LBi UBi LBj UBj FORCES ng cloudG FORCES ng cloud update SV PYPP S amp F 2 amp In other cases the parent routine would have the use then would pass the relevant array to the _ tile routine USE mod_grid CALL prsgrd_tile ng tile amp amp GRID ng Hz amp Sinan prsgrd_tile ng tile amp amp Hz z_r ZW amp real r8 intent in Hz LBi LBj This allows the _ tile routine to use Hz with the same syntax as the pre parallel pre module code once had 6 9 4 Input output In ROMS the distributed memory I O is all happening on the master process 0 unless you specifically ask it to use MPI I O which requires both the NETCDF4 and PARALLEL_IO cpp flags to be defined If you do this you will be reading and writing HDF5 files and will need 72 to update your pre and post processing tools accordingly I have tentatively tried the parallel I O and found it to be exceedingly slow I ve been told since that this is the fault of the NetCDF 4 layer sitting on top of HDF5 HDF5 alone should be fast In the case of having all the I O pass through the master process we can still read
15. epineutral sur faces MIX_5S_ TS Define for diffusion along constant s surfaces CLIMA_ TS_ MIX Define for diffusion of tracer perturbation T Telm Vertical mixing BVF_MIXING Define to activate Brunt Vaisala frequency mixing GLS_ MIXING Define for Generic Length Scale mixing 60 CANUTO_A Define for Canuto A stability function formulation CANUTO_B Define for Canuto B stability function formulation CHARNOK Define for Charnok surface roughness from wind stress CRAIG_BANNER Define for Craig and Banner wave breaking surface flux KANTHA_CLAYSON Define for Kantha and Clayson stability function K C2ADVECTION Define for 2th order centered advection K C4ADVECTION Define for 4th order centered advection N2S2_ HORAVG Define for horizontal smoothing of buoyancy shear ZOS_HSIG Define for surface roughness from wave amplitude TKE_WAVEDISS Define for wave breaking surface flux from wave amplitude LMD_MIXING Define to activate Large McWilliams Doney interior closure LMD_BKPP Define to add a bottom boundary layer from a local K Profile Parameterization KPP MD_CONVEC Define to add convective mixing due to shear instabilities MD_DDMIX Define to add double diffusive mixing MD_NONLOCAL Define to add convective nonlocal transport MD_RIMIX Define to add diffusivity due to shear instabilities MD_SHAPIRO Define to Shapiro filtering boundary layer depths MD_SKPP Define to add a surface boundary layer from a local K Profile Parameterizatio
16. for some runs at least Notice that comments are allowed in this file as long as they are in the C C syntax Also conditionals are fine some parts are only wanted when SOLVESD is on There is evidence of options having been tried and rejected such as the quadratic bottom drag in the presence of tides The Eastern side is actually open to the North up in the Arctic but it is set to be a closed boundary Instead there is an imposed flow through Bering Strait by using UV_ PSOURCE and ANA_ PSOURCE 7 3 2 NEP5 code chunks As mentioned above we are specifying point sources in ana_ psource h The NEP5 section of the copy in Apps NEP is Hif defined NEP5 Nsrc 149 DO is 1 Nsrc Dsrc is 0 0_r8 Isrc is 225 Jsrc is 472 i8 Lsrc is itemp FALSE Lsrc is isalt FALSE END DO elif defined BERING_10K if defined UV_PSOURCE defined Q_PSOURCE ifdef SOLVE3D I If appropriate set up nondimensional shape function to distribute mass point sources sinks vertically It must add to unity I if defined NEP5S DO k 1 N ng DO is 1 Nsrc Qshape is k 1 0_r8 REAL N ng r8 END DO 110 END DO endif endif if defined NEP5 cff 800000 _r8 Nsrc DO is 1 Nsrc Qbar is cff END DO else ana_psource h No values provided for Qbar endif The total transport is set to 0 8 Sv and is divided among the Nsrc sources Doing a search for NEP5 shows up a few other instances One is in output F forcing more digits
17. on the application NetCDF history and restart files are also created containing the model fields at the requested times We have asked it to save restart records every day In this case the restart file has been told to cycle or to write over the second last record The restart file at the end of the run contains the fields at day 4 and day 5 The history file contains records for 0 1 4 1 2 3 4 and 1 day etc while the averages and diagnostics files are at 1 8 3 8 5 8 and 7 8 day etc Plots can be made from any one of these files using the plotting software described in Selected frames from such plots are shown in Fig 17 to 20 7 3 Northeast Pacific example The upwelling downwelling examples is one in which all the start up fields are defined analytically The other extreme is one in which everything is read from files as in our Northeast Pacific simu lations Figure 21 shows the bathymetry and the extent of this domain which is rectangular in a conic map projection 7 3 1 nep5 h The C preprocessor variable NEP5 has been chosen to label our Northeast Pacific domain using the fifth generation grid file The header include file then becomes nep5 h Options for Northeast Pacific NEP5 simulation undef NETCDF4 undef PARALLEL_IO undef OFFLINE_FLOATS general define CURVGRID define MASKING define NONLIN_EOS define SOLVE3D define SALINITY ifdef SOLVE3D 99 ROMS 3 2 Upw
18. t stress tensor Ta air stress Tw water stress u v the x y components of ice velocity U Vio Fw 10 meter air and surface water velocities Pa Pw 1 3kgm 1025kgm air and water densities C x y t nonlinear bulk viscosity Culo y t height of the ocean surface Table 4 Variables used in the ice momentum equations 39 The ice contains brine pockets for a total ice salinity of 5 The surface ocean temperature and salinity is half a dz below the surface The water right below the surface is assumed to be at the freezing temperature a logarithmic boundary layer is computed having the temperature and salinity matched at freezing Here the W variables are the freeze or melt rates as shown in Fig 7 and Table 5 The frazil ice growth Wy will be discussed further in 5 2 2 note that it contributes to changes in A as well as to changes in hi The other term that contributes to A is Wao This term includes a factor which Mellor and Kantha set to different values depending on whether ice is melting or freezing 4 0 Wao 0 195 b 0 5 Wao lt 0 196 197 Wir Wi Figure 7 Diagram of the different locations where ice melting and freezing can occur Figure 8 shows the locations of the ice and snow temperatures and the heat fluxes The tem perature profile is assumed to be linear between adjacent temperature points The interior of the ice contains brine pockets leading to a prognostic
19. this will provide the storage for the surface run ning means of the fields you are averaging mod_ bbl F If BBL_ MODEL is defined this will provide the storage for the bottom boundary fields mod _ _biology F If BIOLOGY is defined this will provide the storage for the biology interaction parameters mod_ boundary F This contains the storage for the open boundary conditions If they aren t provided analytically this will also provide the storage for fields read from a file that need to be time interpolated mod_ clima F If one of CLIMATOLOGY or several other options is defined this will provide the storage for the climatology fields mod__coupler F If either MODEL _ COUPLING or ESMF_ LIB is defined this will set up the requisite fields and data structures for the coupling mod__coupling F If SOLVE3D is defined this will provide the storage for the fields used in coupling the 2 D and 3 D components of the simulation mod_ diags F If DIAGNOSTICS is defined this will provide the storage for the various ten dency terms mod_ eclight F If both BIOLOGY and ECOSIM are defined this will set up the spectral irradiance variables mod_ eoscoef F If NONLIN_ EOS is defined this will provide the polynomial expansion coef ficients for the nonlinear equation of state for sea water mod_ filter F If FILTERED is defined this will provide the storage for the weighted means used in detiding the averages mod_ floats F If FLOATS is define
20. 10 000 MAX DEPTH 7402 3 Ao 2 BR 20 30 40 50 bD T ODEH2 2 006402 4 0DEH2 b ODE R B ODEH 1 006403 1 506403 2 00 E 03 2 50 03 5 0E 03 4 0E 03 Figure 21 Bathymetry of the Northeast Pacific domain NEP5 104 undef SPLINES endif define FLOATS define STATIONS undef WET_DRY undef T_PASSIVE ifdef T_PASSIVE define ANA_PASSIVE endif ice ifdef SOLVE3D define ICE_MODEL ifdef ICE_MODEL define ICE_THERMO define ICE_MK undef ICE_ALB_EC92 undef ICE_SMOOTH define ICE_MOMENTUM define ICE_MOM_BULK define ICE_EVP define ICE_ADVECT define ICE_SMOLAR define ICE_UPWIND define ICE_BULK_FLUXES define ANA_AIOBC define ANA_HIOBC define ANA_HSNOBC endif endif output stuff define NO_WRITE_GRID undef OUT_DOUBLE define RST_SINGLE define AVERAGES undef AVERAGES2 ifdef SOLVE3D undef AVERAGES_DETIDE define AVERAGES_AKT define AVERAGES_AKS define AVERAGES_AKV define AVERAGES_FLUXES undef AVERAGES_QUADRATIC undef DIAGNOSTICS_TS endif undef DIAGNOSTICS_UV HH HF H HF OF 105 advection dissipation pressure grad etc ifdef SOLVE3D define DJ_GRADPS Hendif define UV_ADV define UV_COR define UV_SADVECTION tifdef SOLVE3D define TS_U3HADVECTION define TS_C4VADVECTION undef TS_MPDATA endif define UV_VIS2 define UV_SMAGORINSKY define VISC_3DCOEF define MIX_S_UV define VISC_GRID defin
21. 12 maintains third order accuracy The most accurate choices for 6 and e are B 17 120 and e 11 20 which yields a fourth order scheme with low numerical dispersion and diffusion and Qmax 1 851640 75 A 7 Generalized FB with an AB3 AM4 Step One drawback of the previous two schemes is the need to evaluate the right hand side r h s terms twice per time step The original forward backward scheme manages to achieve Qmax 2 while only evaluating each r h s term once It is possible to contruct a robust scheme using a combination of a three time AB3 like step for with a four time AM4 like step for u Sia p 3 F u G 28 Flu a u At on unti u G ay 2e got G 2y 3e GCC VGC GC At This is second order accurate for any set of values for P y and e It is third order accurate if B 5 12 However it has a wider stability range for 0 281105 Shchepetkin and McWilliams 75 choose to use this scheme for the barotropic mode in their model with 6 0 281105 y 0 0880 and e 0 013 to obtain Qmax 1 7802 close to the value of 2 for pure FB but with better stability properties for the combination of waves advection and Coriolis terms 122 B The vertical v coordinate Following Song and Haidvogel 79 but modified by Shchepetkin and McWilliams 73 the vertical coordinate has been chosen to be 2 1 0 h 0 h h Clo 1 lt 0 lt 0 243 where he is either the minimum depth or
22. 262 An estimate of the cloud fraction c will be provided by Jennifer Francis 17 E 2 Longwave radiation The clear sky formula for incoming longwave radiation is given by F oT 1 0 261 exp 7 77 x 107 273 T 263 while the cloud correction is given by LW 1 0 2750 F 264 E 3 Sensible heat The sensible heat is given by the standard aerodynamic formula H paCpCuVug Ta Toso 265 E 4 Latent heat The latent heat depends on the vapor pressure and the saturation vapor pressure given by es 611 x 10 7se7273 16 Tsfe b 267 The vapor pressures are used to compute specific humidities according to Ee 268 s 7 269 ds p 1 eJes 127 Variable Value Description a b 9 5 7 66 vapor pressure constants over ice a b 7 5 35 86 vapor pressure constants over water c cloud cover fraction CE 1 75 x 1078 transfer coefficient for latent heat CH 1 75 x 1078 transfer coefficient for sensible heat Cp 1004 J kg K7 specific heat of dry air o declination e vapor pressure in pascals s saturation vapor pressure 0 622 ratio of molecular weight of water to dry air HA hour angle L 2 5 x 10 J kg latent heat of vaporization L 2 834 x 10 J kg latent heat of sublimation latitude Qo incoming radiation for cloudless skies ds surface specific humidity diom 10 meter specific humidity Pa air density S 1353 W m solar constant o 5 67 x 1078
23. 7 3 2 NEP5 code chunks 0 0 00 000 eee eee 7 3 3 Model dama 7 3 4 Initial and boundary conditions 02 02 2004 7 3 0 Forcing 7 3 6 ocean in 7 3 7 Output 8 Plotting Programs for Postprocessing A Model Time stepping Schemes A l Euler A 2 Leapfrog A 3 Third order Adams Bashforth ABS oo o EROS A A AA A A eee oo A 5 Forward Backward Feedback RK2 FB o o bs io ADRs Se eae a ee A 7 Generalized FB with an AB3 AM4 Step 2 0 202 20 e e B The vertical c coordinate C_ Horizontal curvilinear coordinates D Viscosity and Diffusion D 1 Horizontal viscosityl ee D 2 Horizontal Diffusion ee D 3 Vertical Viscosity and Diffusion 0 0 0 0 000002 eee eee E Radiant heat fluxes E l Shortwave radiation sa s so so ma a oaoa sra a a iao uoga a i a a a e aa E 2 Longwave radiation ooa aaa a E 3 Sensible heat E 4 Latent heat F The C preprocessor F 1 File inclusion E2 Macro substituto F 3 Conditional inclusion 0 a a a a F 4 C comments F 5 A note on style F 6 Potential problems a F 7 Modern Fortran 117 120 120 120 120 121 121 122 122 123 125 126 126 126 126 127 127 127 127 127 ue He ey ee eee ee eae ors Se a reese ye oe a Dia hye eee bans de oe Sees a oe ee ee E eee a Sind pid Ge Hema see bdo wes A B
24. 88 define SALINITY define SOLVE3D define SPLINES define AVERAGES define DIAGNOSTICS_TS define DIAGNOSTICS_UV define EW_PERIODIC define ANA_GRID define ANA_INITIAL define ANA_SMFLUX define ANA_STFLUX define ANA_SSFLUX define ANA_BTFLUX define ANA_BSFLUX if defined GLS_MIXING defined MY25_MIXING define KANTHA_CLAYSON define N252_HORAVG else define ANA_VMIX endif if defined BIO_FENNEL defined ECOSIM defined NPZD_POWELL defined NEMURO A defined BIO_UMAINE define ANA_BIOLOGY define ANA_SPFLUX define ANA_BPFLUX define ANA_SRFLUX endif if defined NEMURO define HOLLING_GRAZING undef IVLEV_EXPLICIT endif ifdef BIO_FENNEL define CARBON define DENITRIFICATION define BIO_SEDIMENT define DIAGNOSTICS_BIO endif ifdef BIO_UMAINE define OXYGEN undef CARBON endif ifdef PERFECT_RESTART 89 undef AVERAGES undef DIAGNOSTICS_BIO undef DIAGNOSTICS_TS undef DIAGNOSTICS_UV define OUT_DOUBLE endif Here we have declared that we want a periodic channel EW__PERIODIC but no masking There is salinity but we re using a linear equation of state The momentum equations have advection Coriolis force and pressure gradients There is both horizontal viscosity and diffusion but they are along constant o surfaces and if you check the input file you find that the horizontal diffusion is set to zero There are ifdefs for various bio
25. Almost as simple is setting the boundary value to a known exterior value e o 170 4 13 4 Flather boundary condition For the normal component of the barotropic velocity one option is to radiate out deviations from exterior values at the speed of the external gravity waves Flather 14 T a gD 171 The exterior values are often used to provide tidal boundary contitions to the barotropic mode However there are times when only the tidal elevation is known A reduced physics option is available for estimating U in that case 4 13 5 Chapman boundary condition The corresponding condition for surface elevation was investigated by Chapman 8 assuming all outgoing signals leave at the shallow water wave speed of gD This can be useful when using the Flather condition on the 2 D momentum equations 0 4 12 The time derivative here can be handled either explicitly or implicitly The model uses an implicit timestep with the term E being evaluated at the new timestep 4 13 6 Radiation boundary condition In realistic domains open boundary conditions can be extremely difficult to get right There are situations in which incoming flow and outgoing flow happen along the same boundary or even at different depths at the same horizontal location Orlanski proposed a radiation scheme in which a local normal phase velocity is computed and used to radiate things out if it is indeed going out 35 This works
26. Assuming transverse isotropy as in Sadourny and Maynard 1997 and Griffies and Hallberg 2000 23 the deviatoric stress tensor can be split into vertical and horizontal sub tensors The horizontal or transverse sub tensor is symmetric it has a null trace and it possesses axial sym metry in the local vertical direction Then transverse stress tensor can be derived from eq and yielding E E w A Bonn ME where F L Ay no z ae 93 pun Ay E poate l 94 Fe L Ay ao x no 95 pon 1 Ang tm a 96 Notice the flux form of eq and and the symmetry between the F S and F terms which are defined at density points on a C grid Similarly the FY and FY terms are symmetric and defined at vorticity points These staggering positions are optimal for the discretization of the tensor it has no computational modes and satisfies first moment conservation The biharmonic friction operator can be computed by applying twice the tensor operator eq and 92 but with the squared root of the biharmonic viscosity coefficient Griffies and Hall berg 2000 23 For simplicity and momentum balance the thickness H appears only when computing the second harmonic operator as in Griffies and Hallberg 2000 4 10 3 Rotated Transverse Stress Tensor In some applications with tall and steep topography it will be advantageous to reduce substantially the contribution of the stress tensor eq and to the vertica
27. General circulation experiments with the primitive equations i the basic experiment Mon Wea Rev 91 99 164 1963 P K Smolarkiewicz and W W Grabowski The multidimensional positive definite advection transport algorithm non oscillatory option J Comp Phys 86 355 375 1990 Y Song A general pressure gradient formulation for ocean models part i Scheme design and diagnostic analysis Mon Wea Rev 126 3213 3230 1998 Y Song and D B Haidvogel A semi implicit ocean circulation model using a generalized topography following coordinate system J Comp Phys 115 1 228 244 1994 M Steele G L Mellor and M G McPhee Role of the molecular sublayer in the melting or freezing of sea ice J Phys Oceanogr 19 139 147 1989 Travis Swicegood Pragmatic Version Control Using Git Pragmatic Bookshelf Raleigh North Carolina and Dallas Texas first edition 2008 P K Taylor and M A Yelland The dependence of sea surface roughness on the height and steepness of the waves J Phys Oceanogr 31 572 590 2001 J Thuburn Multidimensional flux limited advection schemes J Comp Phys 123 74 83 1996 I B Troen and L Mahrt A simple model of the atmospheric boundary layer sensitivity to surface evaporation Boundary Layer Meteor 37 129 148 1986 L Umlauf and H Burchard A generic length scale equation for geophysical turbulence models J Marine Res 61 235 265 2001 157 86 R C Wajsowicz A co
28. In this case the coefficients become AM 0 159 A 2 N sn 160 BQ Hat oe 161 BIN Ant at 163 C 1 Nm at 164 C N 0 165 Di H Atmnky 43 7 n 166 D 2 Nm H 0 AtmnRg 167 a He ke 88 ya OR 168 DIN HRON Atmn Roy TNO din Kote 169 This is a standard tridiagonal system for which the solution procedure can be found in any standard reference such as Press et al 65 34 4 13 Boundary Conditions ROMS comes with a variety of boundary conditions including open closed and periodic See Marchesiello et al for a more thorough exploration of the options Some options require a value for the boundary points from either an included analytic expression 86 5 or from an external NetCDF file Here represents the exterior value of a quantity 4 13 1 Gradient boundary condition This boundary condition is extremely simple and consists of setting the gradient of a field to zero at the edge The outside value is set equal to the closest interior value 4 13 2 Wall boundary condition ROMS now assumes a wall condition if no other boundary condition is chosen This is a zero gradient condition for tracers and the surface elevation and zero flow for the normal velocity For tangential velocities the wall is treated as either no slip or free slip depending on the value of gamma2 chosen by the user 4 13 3 Clamped boundary condition
29. M do i 2 Lm Uflux i j Uflux i j pmask i j enddo enddo endif MASKING will not be in the Fortran source code if MASKING has not been defined Likewise ifndef tests for a macro being undefined ifndef RMDOCINC c rmask Mask at RHO points O Land 1 Sea c pmask Slipperiness mask at PSI points O Land 1 Sea 1 gamma2 boundary c umask Mask at U points O Land 1 Sea c vmask Mask at V points O Land 1 Sea c CF22 endif There are also else and elif else if statements although elif is newer and is not sup ported by all versions of cpp An example using else and elif is shown ifdef SOLITON real r8 g 1 0_r8 non dimensional elif defined WBC_1 defined WBC_2 defined WBC_3 real r8 g 9 8_r8 m s2 elif defined CIRCLE real r8 g 3 92e 2_r8 m s2 else real r8 g 9 81_r8 m s2 endif Actually ifdef is a restricted version of the more general test if expression where expression is a constant integer value Zero evaluates to false and everything else is considered true Compound expressions may be built using the C logical operators amp amp logical and logical or logical not 130 These symbols would be used as in the following example Hif defined CANYON_A defined CANYON_B do j 0 M do i 0 L yc c32000 c16000 sin pix xr i j x1 24 h i j c20 p5 hmax c20 c1 tanh yr i j yc c10000 enddo enddo endif F 4 C comments
30. Sets up analytic fields as if they came from NCEP ana_nudgcoef Sets up spatially dependent nudging coefficients for nudging to a climatology ana_ pair Computes analytic sea level air pressure ana__passive Computes analytic initial conditions for passive tracers ana_perturb Computes analytic perturbations to the initial conditions ana__psource Computes analytic point source fluxes ana_rain Computes analytic rainfall ana_ scope Sets adjoint sensitivity spatial scope masking arrays ana_sediment Computes analytic initial conditions for the sediment tracers ana_smflux Computes analytic kinematic surface momentum flux wind stress ana_specir Sets surface solar downwelling spectral irradiance at just beneath the sea surface ana_ spinning Sets time variable rotation force as the sum of Coriolis and Centripetal accelera tions This is used in polar coordinate applications annulus grid ana_srflux Computes analytic kinematic surface shortwave radiation ana_ssh Computes analytic sea surface height ana_sss Computes analytic sea surface salinity ana_sst Computes analytic sea surface temperature and dQdSST which are used in the surface heat flux correction ana_stflux Computes analytic kinematic surface flux of tracer type variables 55 ana_tair Computes analytic air temperature ana_tclima Computes analytic tracer climatology fields ana_tobc Computes analytic open boundary conditions for all tracers active passive biology
31. TCLINE 10 0d0 surface stretching parameter bottom stretching parameter critical depth m 112 7 3 4 Initial and boundary conditions The best fields we know of at the moment are known as SODA a reanalysis from 1958 through 2005 by Carton et al 6 There are Matlab scripts to create initial and boundary conditions from the SODA files 7 3 5 Forcing We are using the CORE forcing files 42 and computing the momentum heat and salt fluxes from the atmospheric conditions and the model s surface temperature There are two options one provided by bulk_flux F and the other provided by ccsm_flux F We are opting to use the second With minimal fussing the CORE files can be used as is on their native grid then interpolated by ROMS internally to the domain at hand 7 3 6 ocean in We use an internal time step of 200 s and an external time step of 10 s The horizontal viscosity and diffusion is TNU2 5 0d0 5 0d0 m2 s VISC2 25 0d0 m2 s plus the sponge layers as shown above along the SE and SW boundaries 7 3 7 Output The model writes voluminous information to standard out as shown for the UP WELLING case It also writes out NetCDF files for restart history daily averages stations and floats Plots can be made from any of these files some examples are shown in Fig z 113 ROMS 3 0 NEP5 200 00 Day Min 4 1923E 00 Max 9 1601E 00 Free Surface m Figure 22 Surface elevation after 200 days s
32. The overall grid is shown in Fig The thick outer line shows the position of the model boundary The points inside this boundary are those which are advanced in time using the model physics The points on the boundary and those on the outside must be supplied by the boundary conditions The three dimensional model fields are carried in four dimensional arrays where the fourth array index refers to one of two or three time levels The tracers have a fifth array index telling which tracer is being referred to For instance itemp 1 refers to potential temperature while isalt 2 refers to salinity The integers i j and k are used throughout the model to index the three spatial dimensions a Index variable for the direction j Index variable for the 7 direction k Index variable for the o direction k 1 refers to the bottom while k N refers to the surface 7 1 8 Initial conditions The initial values for the model fields are provided by either ana_ initial or get_ state get_ state is also used to read a restart file if the model is being restarted from a previous run Also in initial rho_ eos is called to initialize the density field rho__eos also computes rhoA the vertically averaged density and rhoS the density perturbation Both rhoA and rhoS are used in the barotropic pressure gradient TT 7 1 9 Equation of state The equation of state is defined in the subroutine rho__eos Two versions are provided in ROMS a nonlinear
33. activate one or more cpp options to run different variants of the same application without modifying its header file If this is the case specify each option here using the D syntax Notice that you need to use the shell s quoting syntax either single or double quotes to enclose the definition if you are using one of the build scripts below NestedGrids Integer number of grids in the setup usually 1 Compiler specific Options These flags are used by the files inside the Compilers directory USE_DEBUG Set this to on to turn off optimization and turn on the g flag for debugging USE_MPI Set this if running an MPI parallel job USE_OpenMP Set this if running an OpenMP parallel job USE_MPIF90 I m frankly not sure about this one I suppose if you have both mpich and some other MPI for a given compiler system pair this could be used to switch between them USE_ LARGE Some systems support both 32 bit and 64 bit options Select this to get 64 bit addressing usually used for programs need more than 2 GB of memory NETCDF_INCDIR The location of the netcdf mod and typesizes mod files NETCDF_LIBDIR The location of the NetCDF library USE_NETCDF4 Set this if linking against the NetCDF4 library which needs the HDF5 library and therefore HDF5_LIBDIR The location of the HDF5 library FORT A shorthand name for the compiler to be used when selecting which system compiler file is to be included from the Compilers directory See section G 2 3 a
34. and sediment ana_vmix Computes analytic vertical mixing coefficients ana_ winds Computes analytic winds ana__wwave Computes analytic wind induced wave amplitude direction and period 6 6 Other subroutines and functions NetCDF I O def_ Creates the ROMS NetCDF file of the appropriate type including dimensions attributes and variables def_info Adds some standard scalar variables to any NetCDF file get_srflux Reads shortwave radiation flux wrt_ Writes to the ROMS NetCDF file of the appropriate type 6 7 C preprocessor variables Before it can be compiled the model must be run through the C preprocessor cpp as described in Appendix F The C preprocessor has its own variables which may be defined either with an explicit define command or with a command line option to cpp We have chosen to define these variables in an application specific include file except for some machine dependent ones which are defined in the makefile These variables allow you to conditionally compile sections of the code For instance if MASKING is not defined the masking code will not be seen by the compiler and the masking variables will not be declared These cpp variables can be grouped into several categories Momentum terms The default horizontal advection is 3rd order upstream bias for 3D momen tum and 4th order centered for 2D momentum The default vertical advection is 4th order centered for 3D momentum If this is the case no flags for
35. coefficient in the linear equation of state Slipperiness parameter gamma2 Slipperiness variable either 1 0 free slip or 1 0 no slip Adjoint sensitivity time parameters DstrS Starting day DendS Ending day Adjoint sensitivity vertical level parameters KstrS Starting level KendS Ending level Adjoint sensitivity logical parameters Lstate isFsur Free surface Lstate isUbar 2D u momentum Lstate isVbar 2D v momentum Lstate isUvel 3D u momentum Lstate isVvel 3D v momentum Lstate isTvar Tracers NT values expected Stochastic optimals time scale SO_decay Stochastic optimals time decorrelation scale days assumed for red noise processes Logicals for stochastic optimals SOstate isUstr Surface u stress SOstate isVstr Surface v stress Lstate isTsur Surface tracer flux NT values expected Stochastic optimals surface standard deviations SO_sdev isUstr Surface u stress SO_sdev isVstr Surface u stress SO_ sdev isTsur Surface tracer flux NT values expected Logical switches to activate the writing of history averages fields Hout idUvel 3 D u velocity component Hout idVvel 3 D v velocity component Hout idWvel 3 D w velocity component 84 Hout idOvel 3 D Q vertical velocity Hout idUbar 2 D u velocity component Hout idVbar 2 D v velocity component Hout idFsur Free surface Hout idBath Time dependent bathymetry Hout idTvar Tracer type variables potential temperature sa
36. equation for the temperature T The surface flux to the air is Qai H LE es LW 1 a SW e 0 T3 273 198 The formulas for sensible heat latent heat and incoming longwave and shortwave radiations are the same as in Parkinson and Washington and are shown in Appendix E The sensible heat is A Qui T RETA 2 TITEL A Qi2 Fr hi N Ti Qio 7 Fr Figure 8 Diagram of internal ice temperatures and fluxes The hashed layer is the snow 40 Variable Value Description Qw 0 10 shortwave albedo of water Qi 0 60 shortwave albedo of wet ice Qi 0 65 shortwave albedo of dry ice Qs 0 72 shortwave albedo of wet snow Qs 0 85 shortwave albedo of dry snow Ck snow correction factor Cpi 2093 J kg K7 specific heat of ice Cpo 3990 J kg K7 specific heat of water Ew 0 97 longwave emissivity of water Ei 0 97 longwave emissivity of ice Es 0 99 longwave emissivity of snow E T r enthalpy of the ice brine system Fr heat flux from the ocean into the ice H sensible heat a fraction of the solar heating transmitted through a lead into the water below ki 2 04 W m t K71 thermal conductivity of ice ks 0 31 W m7 K71 thermal conductivity of snow Li 302 MJ m latent heat of fusion of ice Ls 110 MJ m latent heat of fusion of snow LE latent heat LW incoming longwave radiation m 0 054 C PSU coefficient in linear T S mS equation contribution to A equation from freezing water Q
37. fine compared to the first radius of deformation Arakawa and Lamb 2 A Vi j 1 bug Ph f Qij wisry O Vi J Figure 1 Placement of variables on an Arakawa C grid 4 1 2 Vertical grid The vertical discretization also uses a second order finite difference approximation Just as we use a staggered horizontal grid the model was found to be more well behaved with a staggered vertical grid The vertical grid is shown in Fig op wo Figure 2 Placement of variables on staggered vertical grid 13 A B C O O O O 0 J O O 0 x u points M N O xX O xX O xX O v points O p points 0 0 0 0 e Y points Figure 3 Masked region within the domain 4 2 Masking of land areas ROMS has the ability to work with interior land areas although the computations occur over the entire model domain One grid cell is shown in Fig 1 while several cells are shown in Fig including two land cells The process of defining which areas are to be masked is external to ROMS and is usually accomplished in Matlab this section describes how the masking affects the computation of the various terms in the equations of motion 4 2 1 Velocity At the end of every time step the values of many variables within the masked region are set to zero by multiplying by the mask for either the u v or p poi
38. gt VO Do Dt VC 2 ou ma Fo C 5 An equation of state is also required p p T S P 6 The variables are shown in Table 3 1 An overbar represents a time average and a prime represents a fluctuation about the mean These equations are closed by parameterizing the Reynolds stresses and turbulent tracer fluxes as 2 OC uw Ko vw Ku C w Kez 7 Equations and express the momentum balance in the x and y directions respectively The time evolution of all scalar concentration fields including those for T x y z t and S z y z t are governed by the advective diffusive equation 5 The equation of state is given by equation 6 In the Boussinesq approximation density variations are neglected in the momentum equations except in their contribution to the buoyancy force in the vertical momentum equation 3 Under the hydrostatic approximation it is further assumed that the vertical pressure gradient balances the buoyancy force Lastly equation expresses the continuity equation for an incompressible fluid For the moment the effects of forcing and horizontal dissipation will be represented by the schematic terms F and D respectively The horizontal and vertical mixing will be described more fully in 4 10 1 Variable Description C x y 2 t scalar quantity i e temperature salinity nutrient concentration Du Pv Do optional horizontal diffusive terms Fu Fu Fo forcing source terms f x y Corioli
39. idLrad Longwave radiation flux Hout idSrad Shortwave radiation flux Hout idEmPf E P flux Hout idevap Evaporation rate Hout idrain Precipitation rate Hout idDano Density anomaly Hout idVvis Vertical viscosity coefficient Hout idTdif Vertical diffusion coefficient for temperature Hout idSdif Vertical diffusion coefficient for salinity Hout idHsbl Depth of the surface boundary layer Hout idHbbl Depth of the bottom boundary layer Hout idMtke Turbulent kinetic energy Hout idMtls Turbulent length scale Ice fields Hout idUice Ice u velocity Hout idVice Ice v velocity Hout idAice Ice concentration Hout idHice Ice thickness Hout idTice Ice surface temperature Hout idHsno Snow thickness Hout idTimid Ice internal temperature Hout idSfwat Surface water melt ponds Hout idTauiw Tauiw Hout idChuiw Chuiw Hout idAgeice Ice age Hout idSig11 Ice internal stress 11 component Hout idSig12 Ice internal stress 12 component Hout idSig22 Ice internal stress 22 component Hout idSOmk Under ice salinity Hout idTOmk under ice temperature Hout idWfr Frazil ice growth rate Hout idWai Air ice melt rate Hout idWao Air ocean ice growth rate 86 Hout idWio Ice ocean ice melt growth rate Hout idWro Surface water runoff rate Inert passive tracers Hout inert Passive tracers NPT values Sediment tracers Hout idBott Sediment tracers MBOTP values User parameters
40. in the averages file names for handling one file per record one record per day for fifty years ifdef NEP5 20 FORMAT a _ i5 5 nc else 20 FORMAT a _ i4 4 nc endif Also in output F is code to ask for a new stations file for each run rather than appending to the old one ifdef NEP5 CALL def_station ng true else CALL def_station ng ldefout ng endif A better fix someday would be to allow the station files to be numbered just like the averages files Why not set Idefout to be true It is shared by many outputs including the floats and you can only restart floats with Idefout set to false One feature which got turned on after the two decade run 35 is SRELAXATION in which Pm nudging the sea surface salinity to a value from a forcing file If you check the nudging code you will find that the timescale for nudging is the same as that used for the open boundary conditions That isn t quite what I want so there s this change to set__vbc F elif defined SRELAXATION stflx i j isalt Tnudg isalt ng Hz i j N ng amp ifdef NEP5 60 days from Tnudg of 360 days amp 6 0_r8 amp endif amp t i j N mg nrhs isalt sss i j There is also a kludge in step__floats F for a yet untested float reuse strategy Having a handy tag like NEP5 allows me to search for these little hacks later The file Apps NEP ana_hmixcoef h tried to sneak by without any NEP5 but the locati
41. meant to be off set to an empty string the string off will test true on an ifdef test 138 G 3 Multiple Source Directories the ROMS Way There s more than one way to divide your sources into separate directories The choices we have made include nonrecursive make and putting the temporary files in their own SCRATCH_DIR directory These include the f90 files which have been through the C preprocessor object files module files and libraries G 3 1 Directory Structure The directory structure of the source code has the top directory a Master directory a ROMS directory with a number of subdirectories and several other directories Master contains the main program while the rest contain sources for libraries and other files Note that the bulk of the source code gets compiled into files that become libraries with the ar command one library per directory There is also a Compilers directory for the system and compiler specific Makefile components G 3 2 Conditionally Including Components The makefile will build the lists of libraries to create and source files to compile They start out empty at the top of the makefile sources i libraries That s simple enough but the list of directories to search for these sources will depend on the options chosen by the user not just in the make options 2 41 but inside the ROMS_HEADER file as well How does this happen Once make knows how to find the ROMS_ HEADER it is used
42. noclean Do not clean already compiled objects Note that the default is to compile serially and to issue a make clean before compiling It is left as an exercise for the user if they prefer different default behavior There are also a few variables which are not recognized by the ROMS makefile but are used locally inside the build script These are MY_PROJECT_DIR This is used in setting SCRATCH_ DIR and BINDIR MY_ROMS_ SRC Set the path to the user s local current ROMS source code This is used so that the script can be run from any directory not necessarily only from the top ROMS directory 2 5 Running ROMS ROMS expects to read a number of variables from an ASCII file details of the file are in 7 1 12 For serial or OpenMP execution the syntax is oceanS or ocean0 lt ocean in gt roms out amp while MPI execution requires oceanM ocean in gt roms out 4 so that each process can read the file Realistically you would only want to run relatively small applications such as UPWELLING interactively on the command line as shown here Also for either of the parallel options you will have to provide some information to ROMS and to the operating system about how many threads or processes to use Parallel computers may also have some sort of batch queuing system in place in which you would submit a job script I have easy access to two Linux clusters with differing details in the systems requiring different job scrip
43. p p T S z from Jackett and McDougall 34 and a linear p T S The linear form is p RO Tcoef T TO Scoef S S0 or p RO Tcoef T TO depending on whether or not SALINITY is defined Specify which equation of state you would like to use with the NONLIN_EOS C preprocessor flag in your application include file The linear coefficients RO TO Tcoef SO and Scoef are set in ocean in Note that we are computing in situ density from potential temperature and salinity Some of the vertical mixing schemes require potential density and some other fields which are computed by rho_eos as well 7 1 10 Boundary conditions The horizontal boundary conditions are provided by the subroutines in u8dbc_im v3dbc_im u2dbc_im v2dbc_im t3dbc_im and zetabc They are called every time step and provide the boundary values for the fields u v u v all tracers and respectively They are currently configured for a closed basin a periodic channel a doubly periodic domain or a domain with various open boundary conditions Each side is controlled independently with WEST being the i 1 boundary EAST being the i L boundary SOUTH being the j 1 boundary and NORTH being the j M boundary These choices are made via the cpp options in your include file Many of the choices for open boundaries require that the model have some boundary val ues for the field in question These can be specified in the appropriate ana_xxx h file for sa
44. r o Oz a h t Ox dy and 3 0 0 z z Zz See ae get te In the stretched coordinate system the vertical boundary conditions become top C 0 5p ou I x y t TE 35 78 2 9 t Ko OC _ Qe Hz Oo pocp and bottom o 1 Bee OP sa Note the simplification of the boundary conditions on vertical velocity that arises from the o coordinate transformation 3 5 Horizontal curvilinear coordinates In many applications of interest e g flow adjacent to a coastal boundary the fluid may be confined horizontally within an irregular region In such problems a horizontal coordinate system which conforms to the irregular lateral boundaries is advantageous It is often also true in many geo physical problems that the simulated flow fields have regions of enhanced structure e g boundary currents or fronts which occupy a relatively small fraction of the physical computational domain In these problems added efficiency can be gained by placing more computational resolution in such regions The requirement for a boundary following coordinate system and for a laterally variable grid resolution can both be met for suitably smooth domains by introducing an appropriate orthogonal 11 coordinate transformation in the horizontal Let the new coordinates be x y and n x y where the relationship of horizontal arc length to the differential distance is given by ds e 5 ae 14 C R
45. rain fall rate NA_ SEDIMENT Define for analytic sediment initial fields NA_SMFLUX Define for an analytic kinematic surface momentum stress NA_SPFLUX Define for analytic surface passive tracers fluxes NA_ SPINNING Define for an analytic time varying rotation force NA_SRFLUX Define for an analytic kinematic surface shortwave radiation NA_ SSFLUX Define for an analytic kinematic surface freshwater flux NA_ SSH Define for an analytic sea surface height NA_SSS Define for an analytic sea surface salinity NA_ SST Define for an analytic SST and 0Q OSST NA_STFLUX Define for an analytic kinematic surface heat flux NA_ TAIR Define for analytic surface air temperature NA_ TCLIMA Define for an analytic tracer climatology NA_ TOBC Define for analytic tracer open boundary conditions NA_VMIX Define for analytic vertical mixing coefficients NA__WIND Define for analytic surface winds ANA __WWAVE Define for an analytic wind induced wave field AAA AAA AAA AAA ee ee ee AAA Horizontal mixing of momentum MIX_GEO_ UV Define for viscosity along constant z geopotential surfaces MIX_S_UV Define for viscosity along constant s surfaces VISC_ GRID Define for horizontally variable viscosity coefficient Horizontal mixing of tracers DIFF_ GRID Define for horizontally variable diffusion coefficient MIX_GEO_ TS Define for diffusion along constant z geopotential surfaces MIX_ISO_TS Define for diffusion along constant potential density
46. requiring large values of v D 2 Horizontal Diffusion We have chosen anything from zero to the value of the horizontal viscosity for the horizontal diffusion coefficient One common choice is an order of magnitude smaller than the viscosity D 3 Vertical Viscosity and Diffusion ROMS stores the vertical mixing coefficients in arrays with three spatial dimensions called Akv and Akt Akt also has a fourth dimension specifying which tracer so that temperature and salt can have differing values Both Akt and Akv are stored at w points in the model horizontal averaging is done to obtain Akv at the horizontal u and v points The values for these coefficients can be set in a number of ways as described in 126 E Radiant heat fluxes As was seen in the model thermodynamics requires fluxes of latent and sensible heat and long wave and shortwave radiation We follow the lead of Parkinson and Washington in computing these terms E 1 Shortwave radiation The Zillman equation for radiation under cloudless skies is S cos Z cos Z 2 7 e x 1075 1 085 cos Z 0 10 Qo 258 where the variables are as in Table 8 The cosine of the zenith angle is computed using the formula cos Z sin sin cos cos cos HA 259 The declination is 23 44 x cos 172 day of year x 27 365 260 and the hour angle is HA 12 hours solar time x 7 12 261 The correction for cloudiness is given by SWl Qo 1 0 6c
47. step the expensive diffusion operations are not carried out during the predictor step The preceeding notes on tracer advection refer to all but the MPDATA option The MPDATA algorithm has its own predictor corrector with emphasis on not allowing values to exceed their original range it therefore gives up the constancy preservation This is most noticeable in shallow areas with large tides 4 8 Advection schemes The advection of a tracer C has an equation of the form H C f o 0 a E n o Ot mn ot j T yl 46 22 where we have introduced the advective fluxes H res Hue 47 n pi Hav A m H 9C F 49 mn 4 8 1 Second order Centered The simplest form of the advective fluxes is the centered second order Ts AG H uC Fe 2 50 T H vC Tl z p E 51 H Q0 52 mn This scheme is known to have some unfortunate properties in the presence of strong gradients such as large over and under shoots of tracers leading to the need for large amounts of horizontal smoothing ROMS provides alternative advection schemes with better behavior in many situations but retains this one for comparison purposes 4 8 2 Fourth order Centered The barotropic advection is centered fourth order unless you specifically pick centered second order as your horizontal advection scheme To get fourth order create gradient terms 53 54 55 The fluxes now become q H e 10G Ez Ps pol
48. theirs conflict s show all options Selecting p will behave as the old version If you get a conflict you need to do one of three things e Merge the conflicted text by hand by examining and editing the conflict markers within the file e Copy one of the temporary files on top of your working file e Run svn revert lt filename gt to throw away all of your local changes Once you ve resolved the conflict you need to let Subversion know by running svn resolved This removes the three temporary files and Subversion no longer considers the file to be in a state of conflict More on this below 1 5 1 Merging conflicts by hand To merge your changes with those from the latest revision in the repository open ana__hmixcoef h in your favorite editor and look for a string of lt characters You should see something like this lt lt lt lt lt lt lt mine Hifndef DISTRIBUTE IF Lanafile and tile eq 0 THEN else IF Lanafile THEN endif ifdef DISTRIBUTE IF Lanafile THEN else IF Lanafile and tile eq 0 THEN endif gt gt gt gt gt gt gt r291 9 After comparing the two code blocks separated by the symbols you decide that you prefer the logic from the repository so you remove the conflict markers and your code so the section now looks like this 151 ifdef DISTRIBUTE IF Lanafile THEN else IF Lanafile and tile eq 0 THEN endif Now you can save
49. version described here is that available through the myroms org svn site with modifications to include sea ice and other minor changes Contents 1 Introduction 2 Getting started 2k INYTOMS OTO sh e e e ee ok Gow ag a eS eM RE ie bane Sate Boe ye he Stray ge ed eae ee Br Ge aad es we ats foe dh e e aba eee Soe a ee eg oe 2 4 Compiling ROMS 2 4 1 Environment Variables for makel o e 2 4 2 Providing the Environment 0002 a a AAA A 2 0 Running ROMS e OTE A Ae Rear Gt ee Sv eee eee 3 Ocean Model Formulation 3 1 Equations of motion 3 2 Vertical boundary conditions o 3 3 Horizontal boundary conditions o a a e a 3 4 Terrain following coordinate system e 3 5 Horizontal curvilinear coordinates o a a a a a a 4 Numerical Solution Technique 4 1 Vertical and horizontal discretization a a a a o e a a a a 4 1 1 Horizontal eTidl s s e as asa diare dOP dar aa as a A 4 1 2 Vertical grid 4 2 Masking of land areas 4 2 1 Velocity a ss crime a a a a a ae Ed oan ecc 4 2 3 Wetting and drying o sos scoe pr a a ee ee EEE do e cr eR ET a ok Oe A A 4 4 Conservation properties ee OR E AN 4 6 Density in the mode coupling e A a8 Ppt bey Ad Aes 4 8 Advection schemes oo ooa e e a a e a 4 8 1 Second order Cemberedl e 48 2 Eourth order Centered sp ser tk esa aake aa a a a aaee e ee ee ee
50. we are simply performing one dimensional operations here Instead we use a centered fourth order scheme in the vertical when the third order upwind option is turned on F Oi jki i C bta Zo ike Cig ke 66 1 16 16 One advantage of UTOPIA over MPDATA is that it can be used on variables having both negative and positive values Therefore it can be used on velocity as well as scalars For the u velocity we have mn Za 1 Pu Au 0 H u z z Pu v 8 Hv n 2 z Pele ba 68 o Hw 1 9 9 1 F a quis 16 bk 16 Ui j k 1 16 ij h42 69 24 while for the v velocity we have y H u 82 Hu ma i 70 ell n v5 70 uN Hav H v EN 2 2 i ar m ar r 71 o H w 1 9 9 1 Fe ES TA 16 Mik l Tg i AZ 72 In all these terms the second derivatives are evaluated at an upstream location 4 9 Determination of the vertical velocity and density fields Having obtained a complete specification of the u v T and S fields at the next time level by the methods outlined above the vertical velocity and density fields can be calculated The vertical velocity is obtained by combining equations and to obtain O H u H v oO AQ O Du O Dv 73 aaa eae ea a a Solving for H Q mn and using the semi discrete notation of 4 4 we obtain S p aS al AQ uD uD uH vH a fja a 2 6 A Al 2 JE 74 The integral is actually computed as a sum from the
51. we will be looking at is nl_ocean h 6 2 3 ROMS_ initialize This is called at the beginning of the run and therefore starts off by finding out how many parallel processes are running and which one this is then calls the following ROMS routines 48 initialize__parallel is in the mod_parallel module and sets up a few variables including some for the built in profiling wclock__on is in timers F and initializes a timer for the built in profiling inp_ par reads in the ASCII input file s used by ROMS mod_arrays allocates and initializes the dynamically sized arrays in ROMS based on the grid sizes read in by inp_ par initial reads in the initial conditions from a NetCDF file or computes them analytically Likewise for the grid plus it sets up the vertical grid spacing to be used and many other details get_ data reads in the first record of time varying forcing fields boundary conditions etc 6 2 4 ROMS _ run This loops over all the steps from the starting iteration to the ending iteration in its argument list The loop consists of calls to get_ data reads in the second and subsequent records of time varying forcing fields boundary conditions etc main3d or main2d solves the full equations described in 4 main3d or the depth integrated version only main2d 6 2 5 ROMS _ finalize This is called at the end of the run whether it was otherwise successful or not The routines called are wrt_rst if the run had an er
52. well for a wave propagating normal to the boundary but has problems when waves approach the boundary at an angle Raymond and Kuo have modified the scheme to account for propagation in all three directions In ROMS only the two horizontal directions are accounted for with the recommended RADIATION_ 2D option Op _ OY 09 ne 6658 0170 173 where Fo de nE 174 ae 35 For bn y 175 MONO __ p 176 These terms are evaluated at the closest interior point in a manner consistent with the time stepping scheme used The phase velocities are limited so that the local CFL condition is satisfied They are then applied to the boundary point using equation 173 again using a consistent time stepping scheme Raymond and Kuo give the form used for centered differencing and a leapfrog time step while ROMS uses one sided differences The radiation approach is appropriate for waves leaving the domain A check is made to see which way the phase velocity is headed If it is entering the domain a zero gradient condition is applied unless the next option is also specified 4 13 7 Mixed radiation nudging boundary condition As described in Marchesiello et al 49 ROMS has an option for providing radiation conditions on outflow and nudging to a known exterior value on inflow This is implemented as a variation on the radiation condition requiring two timescales the inflow nudging timescale and the outflow nudging times
53. 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 Current Step 0333333333333333 0666936481246126 1001008081383191 1335675599008453 1670995344713988 2006952520463535 2343461315221331 2680365077488784 3017436598143192 3354378543044362 3690824080946679 4026337758325830 4360416678801077 4692492050905598 5021931174029188 5348039938429148 5670065921276978 5987202166780104 6298591744489550 6603333185976138 6900486906129490 7189082721405745 7468128583421599 7736620652363947 7993554840756058 8237939964192925 8468812641730052 8685254094681655 8886408998655913 9071506549726561 9239883911711814 9391012217603267 9524525304259083 9640251365547480 9738247715198163 9818838856691055 9882658063583314 9930692679747333 9964333355064079 9985427438187774 41 0 0217382098544504 0 0010909315881854 0 9890102622088882 0 9996336754069628 42 0 0109897377911118 0 0003663245930371 1 0000000000000000 0 9999999999999998 ndtfast nfast 30 42 nfast ndtfast 1 40000 Centers of gravity and integrals values must be 1 1 approx 1 2 1 1 1 000000000000 1 047601458608 0 523800729304 1 000000000000 1 000000000000 Power filter parameters Fgamma gamma 0 28400 0 18933 00000000E 00 km 00000000E 00 km 00000000E 00 km 00000000E 00 km 57503054E 00 m 30290891E 01 m Minimum X grid spacing DXmin Maximum X grid spacing DXmax Minimum Y grid spacing DYmin Maximum
54. 0 00 2 317054E 04 6 579497E 02 6 579499E 02 3 884376E 11 WRT_HIS wrote history fields Index 1 1 into time record 0000004 97 WRT_AVG wrote WRT_DIAGS wrote 288 1 00 00 00 WRT_HIS wrote WRT_AVG wrote WRT_DIAGS wrote WRT_RST wrote 1440 5 00 00 00 3 884376E 11 WRT_HIS wrote 0000021 WRT_AVG wrote 0000020 WRT_DIAGS wrote 0000020 WRT_RST wrote 0000001 averaged fields into time record 0000003 diagnostics fields into time record 0000003 5 382882E 04 6 579497E 02 6 579503E 02 3 884376E 11 history fields Index 1 1 into time record 0000005 averaged fields into time record 0000004 diagnostics fields into time record 0000004 re start fields Index 1 1 into time record 0000001 2 476731E 02 6 579533E 02 6 579781E 02 history fields Index 1 1 into time record averaged fields into time record diagnostics fields into time record re start fields Index 1 1 into time record Elapsed CPU time seconds Thread 0 CPU 557 021 Total 557 021 Nonlinear model elapsed time profile e A Baba D EEEE MERE 0 015 0 0027 Reading of input data cee eee eee eee 0 002 0 0004 Processing of input data oooooococooocmomoo o 0 127 0 0229 Processing of output time averaged data 55 303 9 9284 Computation of vertical boundary conditions 0 419 0 0752 Computation of global information integrals 2 316 0 4159 Writing of output data
55. 00 25 863 39 126 0 2605826 14 175 30 436 46 696 0 3178051 15 750 35 581 55 412 0 3862333 17 325 41 426 65 527 95 4 0 7500000 0 4682798 18 900 48 121 3 0 8125000 0 5668375 20 475 55 846 2 0 8750000 0 6853816 22 050 64 818 1 0 9375000 0 8280918 23 625 75 298 O 1 0000000 1 0000000 25 200 87 600 Time Splitting Weights ndtfast 30 nfast 42 Primary Secondary Accumulated to 1 0 0008094437383769 0 0333333333333333 0 0008094437383769 2 0 0014053566728197 0 0333603147912792 0 0022148004111966 3 0 0017877524645903 0 0334071600137066 0 0040025528757869 4 0 0019566842408176 0 0334667517625262 0 0059592371166046 5 0 0019122901320372 0 0335319745705535 0 0078715272486418 6 0 0016548570247459 0 0335957175749547 0 0095263842733877 7 0 0011849025289723 0 0336508794757796 0 0107112868023601 8 0 0005032751608632 0 0336903762267453 0 0112145619632232 9 0 0003887272597151 0 0337071520654408 0 0108258347035082 10 0 0014892209965583 0 0336941944901169 0 0093366137069498 11 0 0027955815694920 0 0336445537902317 0 0065410321374578 12 0 0043042707117221 0 0335513677379153 0 0022367614257357 13 0 0060106451121704 0 0334078920475245 0 0037738836864347 14 0 0079087469427945 0 0332075372104522 0 0116826306292293 15 0 0099910761708919 0 0329439123123590 0 0216737068001212 16 0 0122483446563884 0 0326108764399960 0 0339220514565096 17 0 0146692120341107 0 0322025982847830 0 0485912634906203 18 0 0172400033810439 0 031713624550
56. 1 while the light red function has weights to produce an estimate at time 1 SCRUM 3 0 Rutgers AGRIF UCLA Non hydrostatic Reference Barotropic LF TR LF AM3 with LF AM3 with Gen FB Gen FB mode FB feedback FB feedback AB3 AM4 AB3 AMA 2 D amex iter V2 2 1 85 2 185 3 1730 1 78 1 3 D momenta AB3 AB3 LF AM3 LF AM3 AB3 mod Tracers AB3 LF TR LF AM3 LF AM3 AB3 mod Internal AB3 Gen FB LF AM3 LF AM3 Gen FB waves AB3 TR FB feedback FB feedback AB3 AM4 Qmax advect 0 72 0 72 1 587 1 587 0 78 Qmax Cor 0 72 0 72 1 587 1 587 0 78 Omax int W 0 72 1 1 14 1 2 1 85 2 1 85 2 1 78 1 Table 3 The time stepping schemes used in the various ROMS versions w t is the Courant number and w ck is the frequency for a wave component with wavenumber k 4 4 Conservation properties From Shchepetkin and McWilliams 2005 73 we have a tracer concentration equation in advective form OC V C 0 24 Se u 24 and also a tracer concentration equation in conservation form e V uC 0 25 The continuity equation V u 0 26 can be used to get from one tracer equation to the other As a consequence of eq 24 if the tracer is spatially uniform it will remain so regardless of the velocity field constancy preservation On the other hand as a consequence of 25 the volume integral of the tracer concentration is conserved i
57. 125 where is a non dimensional coordinate ranging from 0 to 1 indicating depth within the surface boundary layer The x subscript stands for one of momentum temperature and salinity 30 Surface Boundary layer depth The boundary layer depth Asy is calculated as the minimum of the Ekman depth estimated as he 0 7u f 126 where u is the friction velocity uz 4 T2 72 p the Monin Obukhov depth L u3 By 127 where k 0 4 is von Karman s contant and By is the surface buoyancy flux and the shallowest depth at which a critical bulk Richardson number is reached The critical bulk Richardson number Rie is typically in the range 0 25 0 5 The bulk Richardson number Rip is calculated as Riv z A IV V d 12 V d where d is distance down from the surface B is the buoyancy B is the buoyancy at a near surface reference depth V is the mean horizontal velocity V the velocity at the near surface reference depth and V is an estimate of the turbulent velocity contribution to velocity shear The turbulent velocity shear term in this equation is given by LMD as 1 2 vaa LEDO yd 129 Rick where C is the ratio of interior N to N at the entrainment depth Gr is ratio of entrainment flux to surface buoyancy flux Cc and e are constants and ws is the turbulent velocity scale for scalars LMD derive based on the expected behavior in the pure convective limit The empirical rule of convection states that
58. 1977 G K Batchelor An introduction to fluid dynamics Cambridge University Press 1967 A Beckmann and D B Haidvogel Numerical simulation of flow around a tall isolated seamount part i Problem formulation and model accuracy J Phys Oceanogr 23 1736 1753 1993 W P Budgell Numerical simulation of ice ocean variability in the barents sea region Towards dynamical downscaling Ocean Dynamics 2005 doi 10 1007 s10236 005 0008 3 J A Carton B S Giese and S A Grodsky Sea level rise and the warming of the oceans in the soda ocean reanalysis J Geophys Res 110 2005 doi 10 1029 2004JC002817 F Chai R C Dugdale T H Peng F P Wilkerson and R T Barber One dimensional ecosystem model of the equatorial pacific upwelling system part i model development and silicon and nitrogen cycle Deep Sea Res I 49 2713 2745 2002 D C Chapman Numerical treatment of cross shelf open boundaries in a barotropic coastal ocean model J Phys Oceanogr 15 1060 1075 1985 Ben Collins Sussman Brian W Fitzpatrick and C Michael Pilato Version Control with Sub version O Reilly amp Associates Inc Sebastopol CA second edition 2008 http svnbook red bean com E E Ebert and J A Curry An intermediate one dimensional thermodynamic sea ice model for investigating ice atmosphere interactions J Geophys Res 98 10085 10109 1993 C W Fairall E F Bradley J E Hare A A Grachev and J B Edson
59. 3127 0 0658312668716642 19 0 0199444086685725 0 0311389577709445 0 0857756755402367 20 0 0227631639997064 0 0304741441486588 0 1085388395399431 21 0 0256737146312911 0 0297153720153352 0 1342125541712341 22 0 0286498597812016 0 0288595815276255 0 1628624139524358 23 0 0316613792205220 0 0279045862015855 0 1945237931729577 24 0 0346736416507075 0 0268492068942347 0 2291974348236652 25 0 0376471948657328 0 0256934188392112 0 2668446296893980 26 0 0405373376992232 0 0244385123436867 0 3073819673886213 27 0 0432936737565710 0 0230872677537126 0 3506756411451923 28 0 0458596469320356 0 0216441452951603 0 3965352880772279 29 0 0481720587108285 0 0201154903974257 0 4447073467880564 30 0 0501605672561820 0 0185097551070648 0 4948679140442383 31 0 0517471682814031 0 0168377361985254 0 5466150823256414 32 0 0528456577069106 0 0151128305891453 0 5994607400325521 33 0 0533610761022577 0 0133513086655816 0 6528218161348097 34 0 0531891349131380 0 0115726061288397 0 7060109510479476 35 0 0522156244733761 0 0097996349650684 0 7582265755213236 36 0 0503158038019030 0 0080591141492892 0 8085423793232266 37 0 0473537721847154 0 0063819206892258 0 8558961515079421 38 0 0431818225418188 0 0048034616164019 0 8990779740497608 39 0 0376397765791564 0 0033640675316746 0 9367177506289172 40 0 0305543017255206 0 0021094083123694 0 9672720523544379 96 77 341 91 216 107 586 126 971 150 000 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
60. BO Se ee A 4 8 4 Third order Upwind 000000 eee eee o pene Gabe 4 10 Horizontal Mime oaa ee i bee BR AE DE eR he be Pa 4 10 2 Transverse stress tensor o ooo a a e e aa 4 10 3 Rotated Transverse Stress Tensor oono o a e a o 410 4 Biarno nacio ia a ee ee le Ree as a e ee ek eed ib baie Boe ee kee eth Ee ee oe Ge ee ee Be a ES YR ERA Ra ed OR eee A eee ee i pate Rhee Gt an ene ara nod 4 13 Boundary Conditions or ed oa Gee ee pa A Oe ee oe eh eee bt Sy ioe win aration a ip pin od ee She Sip ete E ep eet weed a ete oe E ee ae da ido oho hy Bee ba ee a dea eee eed Peas i Gagan ees 5 Ice Model Formulation Ost Dynamics sar 4 4 ebb ee Ree OS Bae be RE AA E wR ee 5 2 Thermodynamics oo ee 5 2 1 Ocean surface boundary conditions 2 0 a ee 5 22 Frazil ice formation 5 2 3 Differences from Mellor and Kanthal 6_ Details of the Code 6 1 Directory Structurel e 6 2 Main subroutines a a aoa a a a a a a a 021 master E s toos ete ee he e apa A ee Be Ae Be ne E ee 6 22 Ocean Control E p socre cc Da a a Se eee ee E E AGA E E E EE da ee E E E bake da aw E oe a Ge E E E E a a A O A C20 mand gi ae brava aa a a aa a a ee aan be eA RR rara aa ad a 64 Modules vita beh a be ee a a e a e Ge Re aS 65 Fumctionals 2 2 a2 2 24 2 eA rad a a be ee a 6 6 Other subroutines and functions 1 ee ph gece Fe ee ene gee Geen ae ee ee g
61. Bulk parameterization of air sea fluxes Updates and verification for the coare algorithm J Climate 16 571 591 2003 K N Fedorov Layer thicknesses and effective diffusivities in the diffusive thermocline con vection in the ocean In J C J Nihoul and B M Jamart editors Small scale turbulence and mixing in the ocean pages 471 479 Elsevier New York 1988 K Fennel J Wilkin J Levin J Moisan J O Reilly and D Haidvogel Nitrogen cycling in the mid atlantic bight and implications for the north atlantic nitrogen budget Results from a three dimensional model Global Biogeochem Cycles 20 2006 doi 10 1029 2005GB002456 R A Flather A tidal model of the northwest european continental shelf Memoires de la Soci t Royale de Sciences de Li ge 6 141 164 1976 M G G Foreman Manual for tidal heights analysis and prediction Technical Report 77 10 Institute of Ocean Sciences Patricia Bay 1977 revised 1996 Pacific Marine Science Report M G G Foreman Manual for tidal currents analysis and prediction Technical Report 77 10 Institute of Ocean Sciences Patricia Bay 1978 revised 1996 Pacific Marine Science Report J A Francis and A Schweiger A new window opens on the arctic Trans Amer Geophys Un 81 77 83 2000 153 18 19 20 21 22 23 24 31 32 33 34 P J S Franks J S Wroblewski and G R Flierl Behavior of a simple plankton model with foo
62. C3P Buoyancy production coefficient plus GLS_SIGK Constant Schmidt number non dimensional for turbulent kinetic energy diffusivity GLS_SIGP Constant Schmidt number non dimensional for turbulent generic statisti cal field psi Constants used in surface TKE flux computation CHARNOK_ ALPHA Charnok surface roughness ZOS_HSIG_ ALPHA Roughness from wave amplitude SZ_ ALPHA Roughness from wave dissipation CRGBAN_CW Craig and Banner wave breaking Bottom drag coefficients RDRG Linear bottom drag coefficient RDRG2 Quadratic bottom drag coefficient Zob Bottom roughness Zos Surface roughness under ice shelves Height of atmospheric measurements for bulk flux parameterizations BLK_ZQ Height of air humidity values BLK_ZT Height of air temperature values BLK_ZW Height of wind values 82 Wetting and drying parameter DCRIT Minimum depth for dry cells Various parameters WTYPE Jerlov water type LEVSFRC Deepest level to apply surface momentum stresses as a body force Used when the C preprocessor option BODYFORCE is defined LEVBFRC Shallowest level to apply bottom momentum stresses as a body force Used when the C preprocessor option BODYFORCE is defined Vertical coordinate parameters Vtransform Transformation equation 1 for old style Vstretching Stretching function 1 for old style Vertical o coordinates parameters theta_s o coordinate surface control parameter 0 lt theta_s lt 20 theta_
63. CLM_NUDGING Define for nudging to 2 D momentum climatology M3CLM_NUDGING Define for nudging to 3 D momentum climatology TCLM_NUDGING Define for nudging to tracer climatology ZCLM_NUDGING Define for nudging to sea surface height climatology Ecosystem models BIO_FENNEL Define for Fennel 13 nitrogen based model 63 BIO_ SEDIMENT Define to restore fallen material to the nutrient pool CARBON Define to add carbon constituents DENITRIFICATION Define to add denitrification processes OXYGEN Define to add oxygen dynamics OCMIP_OXYGEN_SC Define if Schmidt number from Keeling et al 36 RIVER BIOLOGY Define for river biology point sources TALK _NONCONSERV Define for nonconservative computation of alkalinity BEST_NPZ Define for Gibson et al personal communication Bering Sea model STATIONARY Define for extra output BENTHIC Define for benthic components ICE_BIO Define for ice algae JELLY Define for jellyfish CLIM_ICE_1D Define if 1 D with ice BIO_UMAINE Define for Chai et al 7 model BIO_GOANPZ Define for Hinckley et al 80 Gulf of Alaska model NPZD_FRANKS Define for NPZD model of Franks et al 18 NPZD_ IRON Define for NPZD model with iron limitation NPZD_POWELL Define for NPZD model of Powell et al 64 IRON_LIMIT Define for iron limitation on phytoplankton growth IRON_RELAX Define for nudging to iron over the shelf ECOSIM Define for bio optical EcoSim model NEMURO Define for Nemuro ecosyst
64. Getting started 2 1 myroms org Starting off with ROMS is not the easiest thing to do and it just seems to be getting more complex as time goes by There are some resources however beginning with the electronic home for ROMS users at Go to register which gives you access to the subversion server for the code and to the discussion forum for all things ROMS There is also a wiki a bug tracking system and even a developer blog The wiki contains parts of this manual but the nature of wikis is that they can be more fluid with more authors than a static document such as this Dave Robertson robertsonQmarine rutgers edu is the one to talk to if you would like to contribute to the wiki 2 2 Prerequisites As mentioned in Chapter 1 ROMS has some external requirements These are e UNIX or UNIX like environment such as e A Fortran 90 compiler e The compiled with the above compiler including the Fortran 90 interface the subversion revision control software See Appendix I and the ROMS wiki e Gnu make version 3 81 or higher Appendix G contains more than you ever wanted to know about this software e A the one from gnu with the traditional flag works well See Appendix F e The Perl scripting language e Matlab is optional but it is a common tool for pre and post processing of ROMS files Make sure you ve got the right environment before attempting to download or compile ROMS 2 3 Acquiring the ROMS code The main R
65. LOATS Define for simulated Lagrangian drifters FLOATS VWALK Define if floats do vertical random walk VWALK_ FORWARD Define for forward time stepping of vertical random walk DEBUGGING Define to suppress timestamps for easier comparisons between files Analytic fields ANA_ BIOLOGY Define for analytic biology initial conditions NA_BPFLUX Define for an analytic bottom passive tracer flux NA_BSFLUX Define for an analytic bottom salt flux NA_BTFLUX Define for an analytic bottom heat flux NA_ CLOUD Define for an analytic cloud fraction NA_ DIAG Define for customized diagnostics NA_FSOBC Define for analytic free surface boundary conditions NA_ GRID Define for an analytic model grid set up NA_HUMIDITY Define for analytic surface air humidity NA_ ICE Define for analytic ice initial conditions NA_ INITIAL Define for analytic initial conditions NA_ M2CLIMA Define for an analytic 2D momentum climatology NA_ M20BC Define for analytic 2D momentum boundary conditions rrr rr re gt P amp P S amp S p 59 NA_ M3CLIMA Define for an analytic 3D momentum climatology NA_ M30OBC Define for analytic 3D momentum boundary conditions NA_ MASK Define for an analytic mask NA_ PAIR Define for an analytic surface air pressure NA_ PASSIVE Define for analytic initial conditions for inert tracers NA_ PERTURB Define for analytic perturbation of initial conditions NA_ PSOURCE Define for analytic point sources NA_RAIN Define for analytic
66. NUSER Number of user parameters USER Values of user parameters NUSER values NetCDF 4 HDF5 parameters NC_ SHUFFLE If non zero turn on shuffle filter NC_DEFLATE If non zero turn on deflate filter NC_DLEVEL Deflate level 0 9 Input NetCDF file names GRDNAME Grid file ININAME Initial conditions file ITLNAME Initial tangent linear file IRPNAME Initial representer file IADNAME Initial adjoint file CLMNAME Climatology file BRYNAME Boundary condition file FWDNAME Forward model file ADSNAME Adjoint sensitivity functional file Forcing NetCDF files NFFILES Number of forcing files FRCNAME Forcing files Output NetCDF file names GSTNAME GST analysis restart file RSTNAME Restart file HISNAME History file TLMNAME Tangent linear file TLFNAME Impulse forcing file for tangent linear model ADJNAME Adjoint file AVGNAME Averages file STANAME Stations file FLTNAME Lagrangian floats file 87 ASCII input file names APARNAM Assimilation parameters SPOSNAM Stations positions FPOSNAM Initial drifter positions IPARNAM Ice parameters BPARNAM Biology parameters SPARNAM Sediment parameters USRNAME User s generic input The bottom of the sample files contain comments describing some of these in greater detail 7 1 13 User variables and subroutines It is possible for the user to add new variables and functionality though it is discouraged The design goal is to isolate the most common feature
67. OCS Study MMS 2009 062 Technical Manual for a Coupled Sea Ice Ocean Circulation Model Version 3 Katherine S Hedstr m Arctic Region Supercomputing Center University of Alaska Fairbanks U S Department of the Interior Minerals Management Service Anchorage Alaska Contract No M0O7PC13368 OCS Study MMS 2009 062 Technical Manual for a Coupled Sea Ice Ocean Circulation Model Version 3 Katherine S Hedstr m Arctic Region Supercomputing Center University of Alaska Fairbanks Nov 2009 This study was funded by the Alaska Outer Continental Shelf Region of the Minerals Manage ment Service U S Department of the Interior Anchorage Alaska through Contract MO7PC13368 with Rutgers University Institute of Marine and Coastal Sciences The opinions findings conclusions or recommendations expressed in this report or product are those of the authors and do not necessarily reflect the views of the U S Department of the Interior nor does mention of trade names or commercial products constitute endorsement or recommenda tion for use by the Federal Government This document was prepared with IAT E X xfig and inkscape Acknowledgments The ROMS model is descended from the SPEM and SCRUM models but has been entirely rewritten by Sasha Shchepetkin Hernan Arango and John Warner with many many other con tributors I am indebted to every one of them for their hard work Bill Hibler first came up with the viscous plast
68. OMS code is available for download via svn at https www myroms org svn src The version of the model described in this document is a merger between ROMS 3 2 and a sea ice model The sea ice code is a branch off a different repository and requires special access contact Dave Robertson as above for more information ROMS comes with several cases all ready to go at the flip of a switch Try these out first and get to understand how they are set up e describes how to pick the cases and set up the build environment g6 7 lists all the ROMS options that can be added to your case e 6 5 lists the fields which can be provided to ROMS via analytic expressions 7 1 19 lists the input parameters ROMS reads from a text file at run time e Chapters 6 and 7 are meant to be informative for the simple and not so simple cases If that isn t the case please let me know In addition to this manual there are some other ROMS resources e You may be best served by going to the ROMS wikil which includes sections called Getting Started and Tutorials e Don t be afraid to use the forum It has everything from employment opportunities to debug ging help Posting there can get you help from one of several people improving your odds of success over private emails Registered users get an email once a day about new postings so you might have to wait a day or more for a reply e There have been ROMS meetings and classes in which a tutorial s
69. R RKO PPPOE PI YY DOOR OO XX XXX 22 NN DOOR XXX XX 222 NN BOXKKXXKKKKKKXKKKKXY RAIN ROXKKXXKKKKKKKKKKKKY RRA DARRO RX XX NN KAXKKXXKKKKKKKKKKXXKXK PONY OOOO OOO OOOO PON 199090000000 RRR ARAN RAK Y 9 INN RI KY 9 RXR ORR RARER C9 DI RRRA RRR KRY 9 RXR ARRAY RRR RRR Y 9 ICAA CONNY 9 RXR YER RRA BRK RRR KKY 9 RXR Les RARA OONN mN 5 RA A Y CS RXR RRNA RI v RXR ORR RORY DI RRR AN RORY v RXR RRR RRR RORY v RARER RORY v RXR KY CCAA RR AAN OOO AAA RORY OOO AAA AA RORY OOO AAA RRR YY OOO AAA RORY OOO ORR OOOO ONO 0000000000 ACA OOOO BRO OO000000 A ORAM RA OOO OOOOOOOOOOOO0000 AAA OOO OOOOOOOOONOOO000 AAA OOO OOOOOOOOOOOO0000 AAA OO OOOOOOOOOOOOOO000 RORY AA OO OOOOOOOOOOOO0000 AXA ANN RRA POLIS ChunkSizeJ in Y seus RAYNOR IRR Prin RAN ANNO AAN NN Y A NN Y A III RAN PONIA OA NN III OOO NN III OA NN III AAN 4 MMR KERR A RIOS O II na III IO III III II O OR ONO AAN AAA Y II tae III III hy III O RAN PONIAN ON O NN toe AAN PONIAN II III aul II III 900000000 NA DI LAI II DAA AA IA ORINAR j d PANY IAH A ORIO ARENA RY d O DI RI II II II II O RN A ONIS A OO DI AIDA d e WN DORR RRNA II III RRR YORK E KRY ROR RRXXY IA EY ERR RY XR XE IA XR RODD RY YY ROR RRR RAY YAK EY EPR rr XK EY A RA KYRA NI RX XX IIA A A RRR RRR II III OA NN O NN OOO IU NOOO PORRA A A Aol 2S SeSe SSSSSSeSeS SST OSOS OSOS SS Jti
70. Subversion places three extra unversioned not under Subversion control files in your working copy filename mine This is your file as it existed in your working copy local copy before you updated your working copy This file has only your latest changes in it If Subversion considers the file to be unmergeable then the mine file isn t created since it would be identical to the working file filename rOLDREV This is the file that was the BASE revision before you updated your working copy That is the file that you checked out before you made your latest edits filename rNEWREV This is the file that your Subversion client just received from the server when you updated your working copy This file corresponds to the HEAD latest revision of the repository For example let s say you checked out revision 280 and made some changes to a file for instance User Functionals ana_hmixcoef h Now you want to update your ROMS source code to take advantage of a new algorithm but when you run svn update your ana_hmixcoef h is now in conflict with the new version in the repository 150 gt svn update U Version U ROMS Modules mod_mixing F U ROMS Functionals ana_hmixcoef h C User Functionals ana_hmixcoef h Updated to revision 291 The above is with an older version of svn The latest greatest does this Conflict discovered in ROMS Utility inp_par F Select p postpone df diff full e edit mc mine conflict tc
71. W m K Stefan Boltzmann constant La air temperature Ta dew point temperature Tsfe surface temperature of the water ice snow Vwg geostrophic wind speed Z solar zenith angle Table 8 Variables used in computing the incoming radiation and latent and sensible heat The latent heat is also given by a standard aerodynamic formula LE paLCEVwg d10m 4s 270 Note that these need to be computed independently for the ice covered and ice free portions of each gridbox since the empirical factors a and b and the factor L differ depending on the surface type 128 F The C preprocessor The C preprocessor cpp is a standalone program which comes with most C compilers On older UNIX systems it was not in the default path but in lib or in usr lib Most modern systems have some version of gec which comes with a quite capable cpp Just be sure to give it the traditional flag This Appendix describes the C preprocessor as used in ROMS with the Fortran language A more complete description is given by Kernighan and Ritchie 37 More practical advice on using cpp is given by Hazard 27 F 1 File inclusion Placing common blocks in smaller files which are then included in each subroutine is the easiest way to make sure that the common blocks are declared consistently Many Fortran compilers support an include statement where the compiler replaces the line include file h with the contents of file h file h is ass
72. Y grid spacing DYmax Minimum Z grid spacing DZmin Maximum Z grid spacing DZmax NOrPrRPrRPRP re Minimum barotropic Courant Number 2 22358627E 01 Maximum barotropic Courant Number 5 42494240E 01 Maximum Coriolis Courant Number 2 47800000E 02 Maximum grid stiffness ratios rx0 6 931666E 02 Beckmann and Haidvogel rxi 1 188435E 00 Haney Initial basin volumes TotVolume 3 8843755884E 11 m3 MinVolume 1 6168290365E 06 m3 MaxVolume 2 3029089059E 07 m3 Max Min 1 4243366824E 01 NL ROMS TOMS started time stepping Grid 01 TimeSteps 00000001 00001440 STEP Day HH MM SS KINETIC_ENRG POTEN_ENRG TOTAL_ENRG NET_VOLUME 0 O 00 00 00 0 000000E 00 6 579497E 02 6 579497E 02 3 884376E 11 DEF_HIS creating history file ocean_his nc WRT_HIS wrote history fields Index 1 1 into time record 0000001 DEF_AVG creating average file ocean_avg nc DEF_DIAGS creating diagnostics file ocean_dia nc 72 O 06 00 00 8 366194E 06 6 579497E 02 6 579497E 02 3 884376E 11 WRT_HIS wrote history fields Index 1 1 into time record 0000002 WRT_AVG wrote averaged fields into time record 0000001 WRT_DIAGS wrote diagnostics fields into time record 0000001 144 O 12 00 00 7 416156E 05 6 579497E 02 6 579498E 02 3 884376E 11 WRT_HIS wrote history fields Index 1 1 into time record 0000003 WRT_AVG wrote averaged fields into time record 0000002 WRT_DIAGS wrote diagnostics fields into time record 0000002 216 O 18 0
73. a ana__psource 7 1 12 ocean in ROMS expects to read a number of variables from an ASCII file as described in Example input files are in ROMS External with names like ocean_grav__adj in where grav_adj refers to the name of the application Lines beginning with are comments and will be ignored by ROMS on reading them The input is organized as key value pairs separated by one or two equals signs It is possible for ROMS to run on more than one grid simultaneously with the number of grids being set at compile time via the build script and or the makefile If there is one equals sign ROMS will use the corresponding value for all grids If there are two ROMS will read a value for each grid Thus far our domains have used just one grid since the inter grid coupling has not been released ROMS will ignore the parameters not needed by the current simulation e g the GLS param eters will not be read if you are not using that mixing scheme However I believe all the example files contain all possible parameters except the ice branch ones The input parameters are in groups as follows Header TITLE A text string to put in the output files MyAppCPP The shorthand name for this application VARNAME The location of the varinfo dat file containing information about fields to read write from to NetCDF files Grid dimension parameters Lm Number of i direction INTERIOR RHO points Mm Number of j direction INTERIOR RHO
74. a shallower depth above which we wish to have more resolution C c is defined as sinh 90 tamh 9 0 5 tanh 50 sinh9 2tanh 30 Cl o 1 244 where 0 and b are surface and bottom control parameters Their ranges are 0 lt 0 lt 20 and 0 lt b lt 1 respectively Equation 243 leads to z for 0 and z h for o 1 Some features of this coordinate system e It is a generalization of the traditional o coordinate system Letting 0 go to zero and using L Hopital s rule we get z h 1l o0 h 245 which is the classic o coordinate e It is infinitely differentiable in o e The larger the value of 0 the more resolution is kept above he e For b 0 the resolution all goes to the surface as 0 is increased e For b 1 the resolution goes to both the surface and the bottom equally as 0 is increased e Some problems turn out to be sensitive to the value of 0 used Figure 25 shows the o surfaces for several values of 0 and b for one of our domains It was produced by a Matlab tool written by Hernan Arango which is available from our web site see 92 1 We find it convenient to define Oz OC o z h he 24 H 5 T Cth hh 0 246 The derivative of C o can be computed analytically oC o cosh 00 coth 50 b 6 247 Oo sinh 0 2 cosh o 3 dd However we choose to compute H discretely as Az Ao since this leads to the vertical sum of H being exactly the tot
75. able second_long_variable becomes double precision really_long_variable second_long_variable and you run the risk of creating lines which are longer than 72 characters in length Also make sure that your macros will not be found anywhere else in the code I used to use define DOUBLE for double precision until it was pointed out to me that DOUBLE PRECISION is perfectly valid Fortran Since a string that is simply defined becomes 1 the macro processor would turn this into 1 PRECISION F 7 Modern Fortran I started working on these ocean models before 1990 much less before Fortran 90 compilers were generally available Fortran 90 s kind feature as used in ROMS is a better way to handle the BIGREAL type declarations On the other hand Fortran 90 does not include conditional compi lation However it is deemed useful enough that the Fortran 2003 standard has the coco means of conditional compilation We might start using this in about ten years 132 G Makefiles One of the software development tools which comes with the UNIX operating system is called make Make has many uses but is most commonly used to keep track of how a large program should be compiled ROMS now requires the gnu version of make sometimes known as gmake Mecklenburg 52 This appendix describes generic make the gnu enhancements to make as well as describing the makefile structure used in ROMS This material first appeared in the ARSC HPC Newsletter
76. acian or biharmonic diffusion D to the turbulent kinetic energy equation The form of this equation in the model coordinates becomes O H 42 Huq a O Hwg E HQ O Ky N Ot mn oE n On m Os mn Os mnH Os 2H K du Ov 7 2H K 2H H HOROR A sa Pide 113 mn Oz Oz mn mnByl mn The vertical boundary conditions are top z a y t mm O vie i 5 02 mnH og 8 62 2_ H K N 2 and bottom z h x y q 0 K dq _ B E 2 maA E Po K y u 0 1 Km 33e p 672 H K N 0 The equation is timestepped much like the model tracer equations including an implicit solve for the vertical operations and an option for using the third order upwind advection There is also an equation for the turbulent length scale l D o alg 2 Dt ta Oz Oz 2 114 TW 114 K 1E P Py 29 where W is the wall proximity function X LY 1 1 aT 116 The form of this equation in the model coordinates becomes 0 Hl g H ugl E Huq a2 H Ql o K 3l B Ot mn ot n an m Os mn Os mnH Os 3 Hz LE P Pp Hed W 4 Hz T mn mnB mn Da 17 where D is the horizontal diffusion of the quantity gel Given these solutions for q and l the vertical viscosity and diffusivity coefficients are Km alSm K miackgrouna 118 Ks qlSn iaa 119 Ky qlSq pd 120 and the stability coefficients Sm Sp and Sy are found b
77. ad eh as Sons aie py eee a ees eau e bas eto hase aes ck oe ia es oe kG ae doh eee de tre oe oe Bee AA AS Kah Gat A e dd pe ae pa e Gel ano dies dia ads atada ss A A asa ere FA ae eee om se a ok oh are nad oes a SS aU eran ee eee ea id nea H sfmakedepend a a sn R ad arar casa ad asa e psa a ed caso bah I 4 Code changes L5 Contlicts 4 6 2 22h seb bbe e eik Ba eae ee A a a 1 5 1 Merging conflicts by hand 20 20000 eee 1 5 2 Copying a file onto your working file 2 aaa 1 5 3 Punting Using svn revert 2 2 ea 133 133 134 134 135 135 136 136 137 138 139 139 139 140 142 142 143 146 147 List of Figures 1 Placement of variables on an Arakawa C grid o 13 2 Placement of variables on staggered vertical grid o 13 3 Masked region within the domaid e o 14 4 Diagrams of the time stepping and mode coupling used in various ROMS versions a Rutgers University ROMS from myroms org b ROMS AGRIF c UCLA ROMS described in 73 d non hydrostatic ROMS 35 In all the curved arrows update 5 The split time stepping used in the model o e 19 6 Weights for the barotropic time stepping The upper panel shows the primary weights centered at time n 5 7 Diagram of the different locations where ice melting and freezing can occur 40 8 Diagram of internal ice tem
78. adourny and K Maynard Formulations of lateral diffusion in geophysical fluid dynamics models In C A Lin R Laprise and H Ritchie editors Numerical Methods of Atmospheric and Oceanic Modelling pages 547 556 NRC Research Press 1997 A J Semtner Jr A model for the thermodynamic growth of sea ice in numerical investigations of climate J Phys Oceanogr 6 379 389 1976 A F Shchepetkin and J C McWilliams Quasi monotone advection schemes based on explicit locally adaptive dissipation Mon Wea Rev 126 1541 1580 1998 A F Shchepetkin and J C McWilliams A method for computing horizontal pressure gradient force in an oceanic model with a nonaligned vertical coordinate J Geophys Res 108 1 34 2003 A F Shchepetkin and J C McWilliams The regional ocean modeling system roms A split explicit free surface topography following coordinates oceanic model Ocean Modeling 9 347 404 2005 A F Shchepetkin and J C McWilliams A correction note for ocean forecasting in terrain following coordinates formulation and skill assessment of the regional ocean modeling system Ocean Modeling 2008 submitted A F Shchepetkin and J C McWilliams Computational kernel algorithms for fine scale multi process long time oceanic simulations In R Temam and J Tribbia editors Handbook of Numerical Analysis Computational Methods for the Ocean and the Atmosphere Elsevier Science 2009 in press J Smagorinsky
79. ai heat flux out of the snow ice surface Qao heat flux out of the ocean surface Qi2 heat flux up out of the ice Qio heat flux up into the ice Qs heat flux up through the snow r brine fraction in ice Pi 910 m3 kg density of ice Si 5 PSU salinity of the ice SW incoming shortwave radiation o 5 67 x 1078 W m K Stefan Boltzmann constant To temperature of the bottom of the ice T temperature of the interior of the ice T temperature at the upper surface of the ice T temperature at the upper surface of the snow T freezing temperature Tmelt i mS melting temperature of ice Tmelt_s 0 C melting temperature of snow Wai melt rate on the upper ice snow surface Wao freeze rate at the air water interface Wp rate of frazil ice growth Wio freeze rate at the ice water interface Wo Wai rate of run off of surface melt water Table 5 Variables used in the ice thermodynamics 41 a function of 73 as is the heat flux through the snow Qs Setting Qa Qs we can solve for T3 by setting quee T AT and linearizing in AT The temperature 73 is found by an iterative solution of the surface heat flux balance using the previous value of T in equation 206 As in Parkinson and Washington if 73 is found to be above the melting temperature it is set to Tmelt and the extra energy goes into melting the snow or ice Qu Qu A 199 Lg E T3 1 E T1 R1 200 Note that L3 1 r L plus a small sensible heat correction We a
80. al water depth D Note that though we have used this form of o coordinate ROMS is written in such a way as to work with a variety of vertical mappings There is one feature which is critical however If the free surface is at rest 0 you get one solution for the level depths 20 k In the case of nonzero the displacements must be proportional to and to the original distance from the bottom z k 2 k 248 or h This ensures that the vertical mass fluxes generated by a purely barotropic motion will vanish at Az k Az k 1 249 every interface 123 Cad MM pains Section m MI 200 22PPPOETEIS e ES Figure 25 The o surfaces for the North Atlantic with a 0 0 0001 and b 0 b 0 8 and b 0 c 0 8 and b 1 d The actual values used in this domain were 0 5 and b 0 4 124 C Horizontal curvilinear coordinates The requirement for a boundary following coordinate system and for a laterally variable grid res olution can both be met for suitably smooth domains by introducing an appropriate orthogonal coordinate transformation in the horizontal Let the new coordinates be x y and n x y where the relationship of horizontal arc length to the differential distance is given by ds d 250 ds dn 251 Here m 7 and n 1 are the scale factors which relate the differential dista
81. algorithm using the largest number of neighbor points In ROMS the halo area would be two grids points wide unless the MPDATA advection scheme is used in which case it needs three The variable GHOST __ POINTS is set accordingly if defined TS_MPDATA defined UV_VIS4 define GHOST_POINTS 3 if defined DISTRIBUTE defined EW_PERIODIC defined NS_PERIODIC define THREE_GHOST endif else define GHOST_POINTS 2 endif The number of tiles is set in the input file as NtileI and NtileJ For an MPI job the product of the two must equal the number of MPI processes For an OpenMP job the number of threads must be a multiple of the number of tiles For instance for NtileI 4 and NtileJ 6 you must have 24 MPI processes while 2 3 4 6 8 12 and 24 are all valid numbers of OpenMP threads Also a serial run could have 24 tiles and would just compute them sequentially Once the input file has been read we can compute the tile sizes ChunkSizeI Lm NtileI 1 NtileI ChunkSizeJ Mm NtileJ 1 NtileJ MarginI NtileI ChunkSizeI Lm 2 MarginJ NtileJ ChunkSizeJ Mm 2 66 Mo x o x O X O X O X O X O X O x u points O v points O p points Figure 12 The whole grid Note that there are Lm by Mm interior computational points The points on the thick outer line and those outside it are provided by the boundary conditions 67 ChunkSizel PIROPO PIRATA RI LI OI NNN RI LI A RO
82. ance for the grid and the wind forcing 7 1 1 Case Name First you need to decide on a name for your particular application or configuration This name is provided via the ROMS_ APPLICATION in either the makefile or the build script This name should be reasonably short all uppercase with spaces converted to underscores For example let s say we pick the name WIKI_TEST This name gets defined during the build so you can add code protected by ifdef WIKI_ TEST as needed This would be a good time to either copy the makefile or the build Script to create one specific to this case prior to editing it 75 7 1 2 Case specific Include File Each application has its own include file included by cppdefs h The name of this file is the name of your application WIKI__TEST here turned into lower case with h appended wiki_test h The location of this file is set by MY_HEADER_ DIR pointing to User Include or some other location of your choosing The complete list of options to be set prior to compilation are listed in Place those you need in the wiki_test h file These include algorithm choices e g advection and turbulence closure schemes boundary conditions output options averages diagnostics stations floats and application modules biology sediments Each line should be of the form define SOME_VAR Note that any undefined variable need not be mentioned Also note that if you copy a predefined application from ROMS
83. and write classic NetCDF 3 files Care must be taken though in the event of an error ROMS has been cleaned up so that the master process will broadcast its return state to the other processes and they can all die gracefully together when there is a problem An example of a routine which reads from disk is get__grid called from initial Each MPI process calls get__ grid CALL get_grid ng iNLM ifdef DISTRIBUTE CALL mp_bcasti ng iNLM exit_flag endif if exit_flag ne NoError RETURN If any one of the processes has trouble it will enter into the exit_ flag which is then shared by all To read in an array variable all processes in get__grid uses nf_fread2d and friends status nf_fread2d ng model ncname ncGRDid ng var_name it var_id it 0 gtype Vsize LBi UBi LBj UBj Fscl Fmin Fmax GRID ng rmask GRID ng rmask IF status ne nf90_noerr THEN exit_flag 2 ioerror status EXIT END IF amp amp amp amp amp amp 2 2 2 F amp F e Within nf fread2d we get to a call to the NetCDF library from just the master process IF InpThread THEN status nf90_get_var ncid ncvarid wrk start total END IF ifdef DISTRIBUTE CALL mp_bcasti ng model status endif IF status ne nf90_noerr THEN exit_flag 2 ioerror status nf_fread2d status RETURN END IF At this point the master process has the entire 2 D array stored in wrk This then needs to be divvied out to the various tiles t
84. andard make feature for providing a list of directories for make to search for files of different types Here we are saying that F files can be found in the directories provided in the modules list and so on for the others 2 For each directory in the modules list make will include the file Module mk that is found there More on these later 3 For each directory in the includes list add that directory to the list searched by cpp with the I flag G 3 3 User defined make Functions The Module mk fragments mentioned before call some user defined functions It s time to show these functions and talk about how they work They get defined in the top Makefile 140 call source dir to binary dir directory list source dir to binary dir addprefix SCRATCH_DIR notdir 1 call source to object source file list source to object call source dir to binary dir subst F 0 filter F 1 call 90 source source file list 90 source call source dir to binary dir subst F 90 1 call make library library name source file list define make library libraries SCRATCH_DIR 1 sources 2 SCRATCH_DIR 1 call source dir to binary dir subst F 0 2 AR ARFLAGS 7 RANLIB endef call one compile rule binary file 90 file source files define one compile rule 1 2 3 cd SCRATCH_DIR FC c FFLAGS notdir 2
85. b o coordinate bottom control parameter 0 lt theta_b lt 1 Tcline Width of the surface or bottom boundary layer in which higher vertical resolution is required during stretching Mean Density and Brunt Vaisala frequency RHOO Mean density used in the Boussinesq approximation BVF_BAK Background Brunt Vaisala frequency squared Time parameters DSTART Time stamp assigned to model initialization days TIDE_START Time of tidal origin relative to model origin TIME_ REF Reference time in the format yyyymmdd dd or else a special value for specific calendars as documented in the ocean in files Nudging time scales TNUDG Time scale days of nudging towards tracer climatology at the interior and at the boundaries A value is expected for each active tracer ZNUDG Time scale days of nudging towards free surface climatology at the interior and at the boundaries M2NUDG Time scale days of nudging towards 2 D momentum climatology at the interior and at the boundaries M3NUDG Time scale days of nudging towards 3 D momentum climatology at the interior and at the boundaries Open boundary factor OBCFAC Ratio of inflow and outflow nudging time scales Linear equation of state parameters 83 RO Background density value used in the linear equation of state TO Background potential temperature constant SO Background salinity constant TCOEF Thermal expansion coefficient in the linear equation of state SCOEF Saline contraction
86. both values are positive the vectors are drawn at interior RHO points xxx LGRID Longitude latitude grid spacing The grid is drawn at LGRID spacing starting from specified map lower corner If LGRID is negative the latitude labels in the map are rotated 90 degrees to avoid label congestion if any NPAGE Number of plots per page Set this parameter to a negative value 1 2 or 4 to activate preservation of the plot aspect ratio Plotting Fields classification derived fields 1 IDUTOT total velocity component in the XI direction cm s 2 IDVTOT total velocity component in the ETA direction cm s 3 IDTVEC total velocity vectors cm s 4 IDTMAG total velocity vector magnitude cm s As you can see there are comments describing what needs to be done Please see the variable ID file for the complete list of fields which can be plotted this list changes as Hernan adds the ability to plot new fields Also check your default cnt file for other vector and contour parameters The palette file includes two number systems one in the scale from 0 to 255 and the other from 0 to 1 The ROMS plotting uses the first while the SEOM plotting uses the second 119 A Model Time stepping Schemes Numerical time stepping uses a discrete approximation to w Ot where represents one of u v C or and F t represents all the right hand side terms In ROMS the goal is to find time stepping schemes which are accu
87. bottom upwards and also as a sum from the top downwards The value used is a linear combination of the two weighted so that the surface down value is used near the surface while the other is used near the bottom The density is obtained from temperature and salinity via an equation of state ROMS provides a choice of a nonlinear equation of state p p T S z or a linear equation of state p p T The nonlinear equation of state has been modified and now corresponds to the UNESCO equation of state as derived by Jackett and McDougall 34 It computes in situ density as a function of potential temperature salinity and pressure Warning although we have used it quite extensively McDougall personal communication claims that the single variable 9 p T equation of state is not dynamically appropriate as is He has worked out the extra source and sink terms required arising from vertical motions and the compressibility of water They are quite complicated and we have not implemented them to see if they alter the flow 4 10 Horizontal mixing In Chapter 3 the diffusive terms were written simply as Du Du Dr and Dg The vertical compo nent of these terms was described in 4 12 Here we describe the ROMS options for representing the horizontal component of these terms 4 10 1 Deviatory stress tensor Note this material was copied from the wiki where it was contributed by Hernan Arango The hor izontal components of the div
88. by cpp to generate an include file telling make about these other options MAKE_MACROS Compilers make_macros mk MACROS shell cpp P ROMS_CPPFLAGS Compilers make_macros h gt MAKE_MACROS CLEAN MAKE_MACROS The make_macros h file contains blocks such as ifdef SWAN_COUPLING USE_SWAN on else USE_SWAN endif The resulting make_macros mk file will simply end up with either USE_SWAN on or USE_SWAN This file can then be included by the makefile and the variable USE_SWAN will have the correct state for this particular compilation We can now use it and all the similar flags to build a list of directories We initialize two lists modules includes ROMS Include 139 Add the adjoint bits ifdef USE_ADJOINT modules ROMS Adjoint endif ifdef USE_ADJOINT includes ROMS Adjoint endif Add the bits we ll always need modules ROMS Nonlinear ROMS Functionals ROMS Utility ROMS Modules includes ROMS Nonlinear ROMS Utility ROMS Drivers Then we add in some more ifdef USE_SWAN modules Waves SWAN Src includes Waves SWAN Src endif modules Master includes Master Compilers Now that our lists are complete let s put them to use vpath F modules vpath h includes vpath 90 SCRATCH_DIR vpath o SCRATCH_DIR include addsuffix Module mk modules CPPFLAGS patsubst 1 includes 1 vpath is a st
89. cale These timescales are provided in the input to ROMS 7 1 12 36 5 Ice Model Formulation The sea ice component of ROMS is a combination of the elastic viscous plastic EVP rheology Hunke and Dukowicz 33 Hunke 32 and simple one layer ice and snow thermodynamics with a molecular sublayer under the ice Mellor and Kantha 55 It is tightly coupled having the same grid Arakawa C and timestep as the ocean and sharing the same parallel coding structure for use with MPI or OpenMP Budgell 5 5 1 Dynamics The momentum equations describe the change in ice snow velocity due to the combined effects of the Coriolis force surface ocean tilt air and water stress and internal ice stress 177 and 178 Mo mjo Mge 43 Fa 177 Ov E ij 1 Ma M fu Moz HT Tyg EY 178 In this model we neglect the nonlinear advection terms as well as the curvilinear terms in the internal ice stress Nonlinear formulas are used for both the ocean ice and air ice surface stress Ta paCalViol Vio 179 1 Ca Ca 1 cos 2a min h 1 5 180 Tw PwCwl Uw Ul Uw U 181 The force due to the internal ice stress is given by the divergence of the stress tensor o The rheology is given by the stress strain relation of the medium We would like to emulate the viscous plastic rheology of Hibler 1979 29 f P Tij 2m61j n kk ij 5 Sis 182 A 1 Ou du e 183 Ei 2 T P P
90. call to output or anything else that would be needing the surface boundary conditions Not in the trunk code 50 ana__vmix is called if there s an analytic profile for the vertical mixing coefficient lmd_ vmix is called when using the K profile parameterization of vertical mixing 41 and 40 bvf_mix computes the vertical mixing as a function of the Brunt Vaisala frequency hmixing computes time dependent horizontal mixing coefficients 76 BI and 23 omega computes the 2 vertical velocity from the horizontal divergences wvelocity computes the physical vertical velocity for the model output set_ zeta sets the surface elevation to the time mean over the last baroclinic time step set_diags accumulates the time average of the diagnostics fields set_ filter accumulates a weighted sum using a Lanczos filter for detiding the most important of the output fields Not in the trunk code set_ avg accumulates time averaged fields for the averages output set_ avg2 accumulates the time averaged surface fields for the second averages output Not in the trunk code output writes to various output NetCDF files rhs3d computes right hand sides of the three dimensional velocity and tracer fields my25_ prestep computes the predictor step for turbulent kinetic energy prognostic variables tke and gls gls_ prestep computes the predictor step for turbulent kinetic energy prognostic variables tke and gls step2d computes the depth in
91. cean in The forcing is applied over the levels from levsfre to N The above caution about vertical resolution also applies to the surface fluxes of T and S although BODYFORCE only refers to wind stress not the surface tracer fluxes b Climatology One way to force the model is via a nudging to the tracer and or momentum climatologies Nudging to tracers was used in the North Atlantic simulations in sponge layers along the northern 78 and southern boundaries Set the climatologies in ana_tclima h or in a file read by get_ data set TCLM_NUDGING in wiki_test h and also set the array Tnudgcof in ana_ nudgcoef h c Tides There is also more than one way to force with tides One way is to provide boundary conditions with enough temporal resolution to resolve the tides Another is to provide ROMS with the tidal constituents at all grid points and to have ROMS reconstruct the tidal currents UV__TIDES and or elevations SSH_TIDES for any given time An example of such a tidal forcing file is in Data ROMS Forcing test_head_ frc nc The non trunk code with ice etc includes the TIDES_ ASTRO option to add on the long period tides from Foreman 15 and 16 There is also an option to include the tidal potential forcing term POT_ TIDES requiring the tidal potential to be included in the tides forcing file d Rivers Point sources can be used to provide river inflow to the model These too can be specified in a forcing file if not provided vi
92. cess to NAG f90 and IBM xlf produce a private my__module mod file if you define module My_ Module in file mod f90 This file is used by the compiler when you use the module as a consistency check type safe programming If foo f90 uses that module you will need the following dependency information foo o my_module mod my_module mod mod o This says that before compiling foo f90 we need to have the file my__module mod This file in turn depends on mod o so that mod f90 must be compiled before foo f90 The sgi is similar except that it uses the file MY _MODULE kmo to store the private module information Use sfmakedepend g on the SGI Rather than creating extra module files the Cray and Parasoft compilers store the module information in the object file and then files which use the modules need to be compiled with extra flags pointing to the module object files For instance if foo f90 uses My __Module which was defined in mod f90 then you will need to compile mod f90 first and provide the Cray compiler with the extra option p mod o when compiling foo f90 When using the Cray use sfmakedepend c to get the dependency information 147 foo o mod o CFT FFLAGS c p mod o foo f90 CFT and FFLAGS are assumed to be previously defined as the name of the compiler and the compiler options respectively Note These f90 module dependencies can confuse some versions of make especially of the System V variety We use gnu make becau
93. ch do not depend on the surface elevation and therefore remain constant in time ana_hmixcoef Computes the horizontal mixing coefficients ana_nudgcoef Computes the nudging time scales ana_initial Analytic initial conditions for momentum and active tracers ana_ passive Analytic initial conditions for passive tracers ana_ biology Analytic initial conditions for ecosystem tracers ana_sediment Analytic initial conditions for sediment tracers ana_ice Analytic initial conditions for ice variables get_ state Reads initial fields from disk either restart or from some other source which has been converted to the appropriate format of NetCDF file set_ depth Compute time evolving depths set_ massflux Compute initial horizontal mass fluxes get_idata Read in time invariant forcing data stiffness Compute grid stiffness grid_ coords Convert initial float and station locations to fractional grid coordinates 6 4 Modules Now that we are using Fortran 90 the method of choice for managing data structures is modules The ROMS Modules directory contains all of the ROMS modules that contain globally used variables The complete list is mod_arrays F This actually has no data structures but has the routine that calls the allocate and initialize routines for all the others mod_average F If AVERAGES is defined this will provide the storage for the running means of the fields you are averaging 52 mod_ average2 F If AVERAGES2 is defined
94. chesiello and J C McWilliams Evaluation and application of the roms 1 way embedding procedure to the central california upwelling system Ocean Modeling 12 157 187 2006 H Peters M C Gregg and J M Toole On the parameterization of equatorial turbulence J Geophys Res 93 1199 1218 1988 N A Phillips A coordinate system having some special advantages for numerical forecasting J Meteorology 14 2 184 185 1957 T M Powell C V W Lewis E N Curchitser D B Haidvogel A J Hermann and E L Dobbins Results from a three dimensional nested biological physical model of the califor nia current system and comparisons with statistics from satellite imagery J Geophys Res 111 C07018 2006 doi 10 1029 2004JC002506 W H Press B P Flannery S A Teukolsky and W T Vetterling Numerical Recipes The Art of Scientific Computing Cambridge University Press 1986 P J Rasch Conservative shape preserving two dimensional transport on a spherical reduced grid Mon Wea Rev 122 1337 1350 1994 W H Raymond and H L Kuo A radiation boundary condition for multi dimensional flows Quart J R Met Soc 110 535 551 1984 R Rew G Davis S Emmerson and H Davies NetCDF User s Guide Unidata University Corporation for Atmospheric Research Boulder Colorado 1996 Version 2 4 156 69 70 71 72 73 74 76 TT 78 79 80 81 82 83 84 85 R S
95. covered many things so far but there s still a few bits which might be confusing 1 There can be rare cases where you might have special code for some systems You can check which system you are on in the F file with ifdef X86_64 special stuff endif To be sure this is defined on each X86__64 system it has to be passed to cpp CPPFLAGS D shell echo 0S tr _ tr a z A Z CPPFLAGS D shell echo CPU tr _ tr a z A Z CPPFLAGS D shell echo FORT tr _ tr a z A Z CPPFLAGS D ROOT_DIR ROOTDIR ifdef ROMS_APPLICATION CPPFLAGS ROMS_CPPFLAGS CPPFLAGS DNestedGrids NestedGrids MDEPFLAGS DROMS_HEADER HEADER endif This guarantees that CPPFLAGS will have terms in it such as DLINUX DX86_64 DPGI D ROOT_DIR export staffdata kate roms trunk DSHOREFACE D HEADER shoreface h D ROMS_HEADER shoreface h DNestedGrids 1 143 2 For mod_ strings F there is a special case SCRATCH_DIR mod_strings f90 CPPFLAGS DMY_OS 0S DMY_CPU CPU DMY_FORT FORT DMY_FC FC DMY_FFLAGS FFLAGS allowing ROMS to report in its output Operating system Linux CPU hardware x86_64 Compiler system pgi Compiler command pgf90 Compiler flags 03 tp k8 64 Mfree Local Root export staffdata kate roms trunk Header Dir ROMS Include Header file shoreface h
96. d this will provide the storage for the float tracking variables mod_ forces F This provides the storage for the surface and bottom forcing fields mod_ fourdvar F If either FOUR_DVAR or VERIFICATION is defined this will set up the variational data assimilation variables mod_ grid F This provides the storage for the model grid fields mod_ ice F If ICE_MODEL is defined this will provide storage for the ice fields mod_iounits F This contains a number of variables used by the I O including file names and file IDs mod_ kinds F This contains the integers associated with the various integer and real Fortran types If you find more systems supporting 128 bit reals let us know mod_mixing F This contains the arrays for the various optional horizontal and vertical mixing parameterizations mod_ncparam F This contains all sorts of parameters relating to the NetCDF I O files includ ing that read from the varinfo dat file The parameters MV and NV are set here giving the maximum number of variables that can be read this is a change from the trunk code 53 mod_nesting F If NESTING is defined this module defines generic structures used for nesting composed and mosaic grids Not yet functional mod_netcdf F This brings in netcdf mod and defines a few type variables based upon it mod_ obs F If either ASSIMILATION or NUDGING is defined this contains variables for the observed fields mod_ocean F This contains the 2 D and 3 D field
97. d level acclimation by herbivores Mar Biol 91 121 129 1986 N G Freeman A M Hale and M B Danard A modified sigma equations approach to the numerical modeling of great lake hydrodynamics J Geophys Res 77 6 1050 1060 1972 B Galperin L H Kantha S Hassid and A Rosati A quasi equilibrium turbulent energy model for geophysical flows J Atmos Sci 45 55 62 1988 J M N T Gray and P D Killworth Sea ice ridging schemes J Phys Oceanogr 26 2420 2428 1996 S M Griffies A Gnanadesikan R C Pacanowski V Larichev J K Dukowicz and R D Smith Isoneutral diffusion in a z coordinate ocean model J Phys Oceanogr 28 805 830 1998 S M Griffies and R W Hallberg Biharmonic friction with a smagorinsky like viscosity for use in large scale eddy permitting ocean models Mon Wea Rev 128 2935 2946 2000 D B Haidvogel H G Arango W P Budgell B D Cornuelle E Curchitser E Di Lorenzo K Fennel W R Geyer A J Hermann L Lanerolle J Levin J C McWilliams A J Miller A M Moore T M Powell A F Shchepetkin C R Sherwood R P Signell J C Warner and J Wilkin Ocean forecasting in terrain following coordinates formulation and skill assessment of the regional ocean modeling system J Comp Phys 227 3429 3430 2007 doi 10 1016 j jcp 2007 01 016 D B Haidvogel H G Arango K Hedstrom A Beckmann P Malanotte Rizzoli and A F Shchepetkin Model evaluation
98. e dn 15 Here m amp 7 and n amp 7 are the scale factors which relate the differential distances A An to the actual physical arc lengths Appendix C contains the curvilinear version of several common vector quantities Denoting the velocity components in the new coordinate system by v u 16 and v 17 the equations of motion 8 13 can be re written see e g Arakawa and Lamb 2 as A 25 Hua ot mn DEN n dix m i mn fN e 1 Um rr a al H 0p gpdz OC 1 O Kym Ou H n 5 Po OE 07 mn Oo E mn Fut Pa 18 Hv H uv H v H v0 z n EE f 1 1 B To 005 o Ca 5 5 e aude F D 19 H mn mn o a 22 HyvC 9 HQC _ On m d00X mn oO H C 09 H uC dt mn OE n 1 Kc dC H mn Os E mn or Bal 20 p p T S P 21 Oo E gHzp do Po O HN 9 Ha Hoy d HO _ 5 5 5 E GS o e All boundary conditions remain unchanged 12 4 Numerical Solution Technique 4 1 Vertical and horizontal discretization 4 1 1 Horizontal grid In the horizontal 7 a traditional centered second order finite difference approximation is adopted In particular the horizontal arrangement of variables is as shown in Fig This is equivalent to the well known Arakawa C grid which is well suited for problems with horizontal resolution that is
99. e gece aves e e e a haemo a eres Benes ees a Sane he eee oe Se E eee Bee a Oe oe NN 6 9 1 ROMS internal numbers NO O O O Guo 69 3 Code syntax s s sasso sreap sannast ka sa ma haa pao AE s se es ye OE we oO ee e ork ee e HE aw 7 Configuring ROMS for a Specific Application 7 1 Configuring ROMS 7 1 1 Case Namel 2 ee 7 1 2 Case specific Include File 2 2 o eee 1 3 Functional cria SSA SE oR RAKES ER aak Ea DG be bees BORG eed Reed eas hee ee eS ee ER eae She ee oe Se bee aoe ee ee ee A ST ca og ak A KS Rae RE are ee ee eae a he AEREO ee ea A ale e a eo he A ae oe Be ek a ee ee we ahs 7 1 8 Initial conditions 2 ee Be dee oe BG He ete oe ee eae Sh ee 7 1 10 Boundary conditions a a 000 eee ee ee ll Model forcing vu baa es ha BRAG EO ee a A Oe ER Ae RR eee T112 eanin circadiano 7 1 13 User variables and subroutines a ao a a a a a a a 37 37 38 43 44 45 46 46 48 48 48 48 49 49 49 52 52 54 56 56 65 65 66 68 69 72 7 2 Upwelling Downwelling Example 4 24 0 24 4 era as 7 2 Cppdets bl rea ee ees w BREN A ee Re ey ee eS Oh ea koe ob Oy Bie wedi a on oe eee ii oe ee 123 ana CTI oso ico PRL RR a RE a Oe Ea EG Oe ee 7 2 4 Initial conditions and the equation of state 2 Tee eesti A 7 2 6 Model forcing 2 maaa e a a aaa ee ee 1 2 1 ocean in 7 2 8 Output 7 3 Northeast Pacificexample 2 a 02 ee 7 3 1 nep5 h
100. e SPONGE ifdef SOLVE3D define TS_DIF2 define MIX_GEO_TS define DIFF_GRID endif vertical mixing ifdef SOLVE3D define SOLAR_SOURCE define LMD_MIXING ifdef LMD_MIXING define LMD_RIMIX define LMD_CONVEC define LMD_SKPP undef LMD_BKPP define LMD_NONLOCAL define LMD_SHAPIRO undef LMD_DDMIX endif undef GLS_MIXING undef MY25_MIXING 106 if defined GLS_MIXING defined MY25_MIXING define KANTHA_CLAYSON define N252_HORAVG endif endif surface forcing ifdef SOLVE3D define CORE_FORCING define BULK_FLUXES define CCSM_FLUXES if defined BULK_FLUXES defined CCSM_FLUXES define LONGWAVE_OUT define DIURNAL_SRFLUX define EMINUSP undef ANA_SRFLUX undef ALBEDO undef LONGWAVE endif endif surface and side corrections ifdef SOLVE3D define SRELAXATION undef QCORRECTION endif point sources rivers line sources Using Runoff instead now ifdef SOLVE3D define RUNOFF define UV_PSOURCE define ANA_PSOURCE undef TS_PSOURCE endif tides define LTIDES ifdef LTIDES define FILTERED define SSH_TIDES define UV_TIDES define ADD_FSOBC define ADD_M20BC undef RAMP_TIDES define TIDES_ASTRO HH HH H HF H OF 107 define POT_TIDES define UV_LDRAG define RDRG_GRID define DRAG_LIMITER undef UV_QDRAG else define UV_QDRAG endif Boundary conditions careful wit
101. e a different extension Each subdirectory is resetting the local_src variable That s OK because we re saving the values in the global sources variable inside the make library function which also adds the local library to the libraries list The compile rules function uses this local__src variable to generate the rules for compiling each file placing the resulting files in the Build directory G 3 5 Main Program The main program is in a directory called Master and its Module mk is similar to the library one 142 local_sub Master local_src wildcard local_sub F local_objs subst F o local_src local_objs addprefix SCRATCH_DIR notdir local_objs sources local_src ifdef LD_WINDOWS BIN libraries local_objs LD FFLAGS local_objs o libraries LIBS_WIN32 LDFLAGS else BIN libraries local_objs LD FFLAGS LDFLAGS local_objs o libraries LIBS endif eval compile rules Instead of a rule for building a library we have a rule for building a binary BIN In this case the name of the ROMS binary is defined elsewhere The binary depends on the libraries getting compiled first as well as the local sources During the link the libraries are compiled from the sources in the other directories while LIBS are exteral libraries such as NetCDF G 3 6 Top Level Makefile Now we get to the glue that holds it all together We ve
102. econd order and fourth order pressure gradient algo rithms in a o coordinate ocean model Int J Num Meth Fluids 18 361 383 1994 155 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 R Mechlenburg Managing Projects with GNU Make O Reilly amp Associates Inc Sebastopol CA 2005 G L Mellor The three dimensional current and surface wave equations J Phys Oceanogr 33 1978 1989 2003 G L Mellor Some consequences of the three dimensional currents and surface wave equation J Phys Oceanogr 35 2291 2298 2005 G L Mellor and L Kantha An ice ocean coupled model J Geophys Res 94 10 937 10 954 1989 G L Mellor and T Yamada A hierarchy of turbulence closure models for planetary boundary layers J Atmos Sci 31 1791 1806 1974 G L Mellor and T Yamada Development of a turbulence closure model for geophysical fluid problems Rev Geophys Space Phys 20 851 875 1982 W A Oost G J Komen C M J Jacobs and C van Oort New evidence for a relation between wind stress and wave age from measurements during asgamage Bound Layer Meteor 103 409 438 2002 I Orlanski A simple boundary condition for unbounded hyperbolic flows J Comp Phys 21 3 251 269 July 1976 C L Parkinson and W M Washington A large scale numerical model of sea ice J Geophys Res 84 6565 6575 1979 P Penven L Debreu P Mar
103. ee surface EAST_FSGRADIENT Define for a gradient condition on the free surface EAST _FSRADIATION Define for a radiation condition on the free surface EAST_FSNUDGING Define for an active passive nudging term on the free surface EAST_FSCLAMPED Define for a clamped free surface EAST _M2FLATHER Define for a Flather condition on the 2 D momentum EAST_M2GRADIENT Define for a gradient condition on the 2 D momentum EAST M2RADIATION Define for a radiation condition on the 2 D momentum EAST_M2REDUCED Define for a reduced physics condition on the 2 D momentum EAST_M2NUDGING Define for an active passive nudging term on the 2 D momen tum EAST_M2CLAMPED Define for clamped 2 D momentum EAST_M3GRADIENT Define for a gradient condition on the 3 D momentum 62 EAST M3RADIATION Define for a radiation condition on the 3 D momentum EAST_M3NUDGING Define for an active passive nudging term on the 3 D momen tum EAST_M3CLAMPED Define for clamped 3 D momentum EAST_KGRADIENT Define for a gradient condition on the TKE fields EAST KRADIATION Define for a radiation condition on the TKE fields EAST_TGRADIENT Define for a gradient condition on the tracers EAST TRADIATION Define for a radiation condition on the tracers EAST_TNUDGING Define for an active passive nudging term on the tracers EAST_TCLAMPED Define for clamped tracers EAST_VOLCONS Define for Eastern edge mass conservation enforcement Tides The tidal data is processed in terms of t
104. elling 0 km 25 150 0 143 8 137 5 131 3 125 0 118 8 112 6 106 3 100 1 93 8 87 6 81 4 75 1 68 9 62 6 96 4 50 2 43 9 37 7 31 4 25 2 0 km Min 2 5200E 01 Max 1 5000E 02 Bathymetry at RHO points m Figure 17 The upwelling downwelling bathymetry 100 ROMS 3 2 Upwelling 1 00 Day Min 7 4098E 00 Max 1 5913E 01 Total Velocity Vectors cm s at Level 16 Figure 18 Surface velocities after one day showing the flow to the left of the wind southern hemisphere 101 ROMS 3 2 Upwelling 1 00 Day 100 0 km 50 50 Min 4 3757E 00 Max 1 8683E 00 Total V velocity cm s Min 1 5903E 01 Max 9 9623E 01 Total U velocity cm s 0 0 m m 100 100 0 km 50 0 km 50 7 Min 1 4606E 01 Max 2 1898E 01 Min 8 5761E 02 Max 9 28642E 02 Potential Temperature C Vertical Velocity cm day Figure 19 Constant slices of the u v T and w fields at day 1 102 ROMS 3 2 Upwelling 5 00 Day 0 0 100 100 0 km 50 0 km 50 Min 6 4684E 01 Max 1 1970E O1 Min 6 8943E 00 Max 6 3215E 00 Total U velocity cm s Total V velocity cm s 0 0 m m bg E 100 100 E 521E d E 0 km 50 Min 1 4498E 01 Max 2 2011E 01 Min 6 3530E 03 Max 4 9777E 03 Potential Temperature C Vertical Velocity em day Figure 20 Constant slices of the u v T and w fields at day 5 103 Bottom Topography MIN DEPTH
105. em model 38 Need to choose a zooplankton grazing option HOLLING_GRAZING or IVLEV_EXPLICIT The de fault implicit IVLEV algorithm does not work yet BIO_SEDIMENT Define to restore fallen material to the nutrient pool HOLLING_GRAZING Define for Holling type s shaped curve grazing im plicit IVLEV_EXPLICIT Define for Ivlev explicit grazing algorithm Sediment transport model SEDIMENT Define to activate sediment transport model 87 BEDLOAD_MPM Define to activate Meyer Peter Mueller bed load BEDLOAD_SOULSBY Define to activate Soulsby wave current bed load RIVER_SEDIMENT Define to activate river sediment point sources SED_DENS Define to allow sediment to affect equation of state SED_MORPH Define to allow bottom model elevation to evolve SUSPLOAD Define to activate suspended load transport Nearshore options WET_DRY Define to allow wetting and drying of cells NEARSHORE MELLOR Define for radiation stress terms from waves NetCDF input output 64 DEFLATE Define to set compression of NetCDF 4 HDF5 format files NETCDF4 Define to create NetCDF 4 HDF5 format files PARALLEL_IO Define to create NetCDF 4 HDF5 format files with MPI I O NO_READ_GHOST Define to not include ghost points during read scatter NO_WRITE_GRID Define to omit writing grid arrays PERFECT_RESTART Define to include perfect restart variables READ_WATER Define to only read water points WRITE_WATER Define to only write water points RST_SINGLE Define to w
106. en interior point computations The routines mp__exchange2d mp_ exchange3d and mp exchange4d can be used to update the neighbors to figure halo points of up to four arrays at a time Each of these routines call tile out which tiles are neighboring and whether or The mp not there really is a neighboring tile on each side exchangexd routines then call _send wait mpi_irecv mpi mpi 10n sav ng irect the north south di then in 1011 The exchanges happen first in the east west direct Y E O a O ost a Po o H a0 yo g as o T ES Ho yo o A ie Q o Y A ie a ti 2 H o hal a 5 o H a0 d lt a o 00 g a a o ES O ost So z yo H amp yo Q Q o a E 4 ie co n o H Y z The updated halo po 68 as 1g needing an update is shown in F 123 45 6 7 8 9 101112131415 1 2 3 4 5 gt EAS 2 3 4 5 a 12 3 4 5 6 7 8 9 10 11121314 15 Figure 14 A choice of numbering schemes a each tile is numbered the same and b each tile retains the numbering of the parent domain 6 9 3 Code syntax In main3d many function calls are surrounded by OpenMP parallel code such as OMP PARALLEL DO PRIVATE thread subs tile SHARED ng numthreads DO thread 0 numthreads 1
107. ergence of the stress tensor Wajsowicz 1993 86 in nondimesional orthogonal curvilinear coordinates 7 s with dimensional spatially varying metric factors 2 t H and velocity components u v wHz are given by 25 i E _ mn O Az0 O HO n O Ogs POSEN al n al m Ue mn gt 1 1 1 OH Hoe 5 Homa 2 088 DE 76 aa n mn 0 H One Oe 0 Ons SM all n al m Fas mn E re 1 1 1 OH Heong ae 7 Hevea Oss On 78 where Te Ay 1 egg v Am en Om v Am egg Am 1 Cm Oss 2V ess Ten Ong 2 Am egn Ogs 2 Km gs Ons 2Ku Ens and the strain field is ege mo mao 85 Cnn 1 s H mn A 86 C 7 UA me oe 87 2egy BO pa 88 Degg q a HS 89 2 ems a ii En He 90 Here Am 7 and Ky y s are the spatially varying horizontal and vertical viscosity coef ficients respectively and v is another very small often neglected horizontal viscosity coefficient Notice that because of the generalized terrain following vertical coordinates of ROMS we need to transform the horizontal partial derivatives from constant z surfaces to constant s surfaces And the vertical metric or level thickness is the Jacobian of the transformation H oz Also in these models the vertical velocity is computed as wile and has units of m s 26 4 10 2 Transverse stress tensor
108. ertical integral of equation is o Du o Duu A o Duv Df Ot mn OE n On m mn O af _A 1 _ _ Dd oc ae Oan 7 2 n 97 Da es 1 Fut Dr 718 33 mn mn where 2 includes the E term D is the horizontal viscosity and the vertical viscosity only contributes through the upper and lower boundary conditions The corresponding vertical integral 18 of equation is Do 0 Pmj 2 Dvvy Dfu timn dE n On m mn _dafl ee oh a o D Opa OC ae a u ae 35 D 1 Fo Dh 73 1 34 We also need the vertical integral of equation 23 shown above as eq 30 The presence of a free surface introduces waves which propagate at a speed of gh These waves usually impose a more severe time step limit than any of the internal processes We have therefore chosen to solve the full equations by means of a split time step In other words the depth integrated equations 33 34 and are integrated using a short time step and the values of u and Y are used to replace those found by integrating the full equations on a longer time step A diagram of the barotropic time stepping is shown in Fig Barotropic steps m 0 m M m M n n 1 Figure 5 The split time stepping used in the model Some of the terms in equations and are updated on the short time step while others are not The contributions from the slow terms are computed once per long time step and stored If we call
109. ession is included as part of the program e There are various resources from these online I ve heard good things about the from Manu Di Lorenzo 2 4 Compiling ROMS 2 4 1 Environment Variables for make ROMS has a growing list of choices the user must make about the compilation before starting the compile process set in user defined variables Since we now use gnu make it is possible to set the value of these variables in the Unix environment rather than necessarily inside the Makefile see G The user definable variables understood by the ROMS makefile are ROMS_ APPLICATION Set the cpp option defining the particular application This is used for setting up options inside the code specific to this application and also determines the name of the h header file for it This can be either a predefined case such as BENCHMARK or one of your own such as NEP5 MY_ HEADER_ DIR Sets the path to the user s header file if any It can be left empty for the standard cases where benchmark h and the like are found in ROMS Include which is already in the search path In the case of NEP5 this is set to Apps NEP where nep5 h resides MY_ ANALYTICAL_ DIR Sets the path to the user s analytic files described in if any This can be User Functionals or some other location I tend to place both the header file and the functionals in the same directory one directory per application MY_CPP_FLAGS Set tunable cpp options Sometimes it is desirable to
110. experiments in the north atlantic basin Simulations in non linear terrain following coordinates Dyn Atmos Ocean 32 239 281 2000 D B Haidvogel and A Beckmann Numerical Ocean Circulation Modeling Imperial College Press 1999 W P Hazard Using cpp to aid portability Computer Language 8 11 49 54 1991 K S Hedstrom Technical manual for a coupled sea ice ocean circulation model version 2 Technical report Institute of Marine and Coastal Sciences New Brunswick NJ June 2000 OCS Study MMS 2000 047 W D Hibler III A dynamic thermodynamic sea ice model J Phys Oceanogr 9 815 846 1979 S Hinckley K O Coyle G Gibson A J Hermann and E L Dobbins A biophysical npz model with iron for the gulf of alaska Reproducing the differences between an oceanic hnlc ecosystem and a classical northern temperate shelf ecosystem Deep Sea Res IT 2009 in press W R Holland J C Chow and F O Bryan Application of a third order upwind scheme in the ncar ocean model J Climate 11 1487 1493 1998 E C Hunke Viscous plastic sea ice dynamics with the evp model linearization issues J Comp Phys 170 18 38 2001 E C Hunke and J K Dukowicz An elastic viscous plastic model for sea ice dynamics J Phys Oceanogr 27 1849 1868 1997 D R Jackett and T J McDougall Stabilization of hydrographic data J Atmos Ocean Tech 12 381 389 1995 154 35 36 41 42 43 44 45 46
111. for shallow vertically well resolved domains TS_DIF2 Define to compute horizontal Laplacian diffusion TS_DIF4 Define to compute horizontal biharmonic diffusion TS_SMAGORINSKY Define for Smagorinsky like diffusion TS_FIXED Define for a diagnostic calculation in which the tracer fields do not change in time T_PASSIVE Define for passive tracers SALINITY Define if salinity is used as one of the active tracers NONLIN_EOS Define to use the nonlinear equation of state 57 QCORRECTION Define to use the net heat flux correction SCORRECTION Define to use freshwater flux correction SOLAR SOURCE Define to use solar radiation source term SRELAXATION Define to use salinity relaxation as a freshwater flux TS_PSOURCE Define for point sources sinks Pressure gradient options If no option is selected the pressure gradient term is computed using standard density Jacobian algorithm Notice that there are two quartic pressure Jacobian options They differ on how the WENO reconciliation step is done and in the monotonicity constraining algorithms DJ_GRADPS Define for splines density Jacobian 72 PJ_GRADP Define for finite volume Pressure Jacobian 45 PJ_GRADPSQ2 Define for quartic 2 Pressure Jacobian 72 PJ_GRADPSQ4 Define for quartic 4 Pressure Jacobian 72 DJ_GRADPS Define for weighted density Jacobian 78 ATM_PRESS Define to impose atmospheric sea level pressure onto the sea surface Atmospheric boundary
112. fusion of each tracer A value is expected for each of the NAT NPT tracers TNU4 Constant mixing coefficient for the horizontal biharmonic diffusion of each tracer A value is expected for each of the NAT NPT tracers Horizontal viscosity coefficients VISC2 Constant mixing coefficient for the horizontal Laplacian viscosity VISC4 Constant mixing coefficient for the horizontal biharmonic viscosity Vertical mixing coefficients for tracers AKT_BAK Background vertical mixing coefficient for the tracers NAT NPT val ues Vertical mixing coefficient for momentum 81 AKV_BAK Background vertical mixing coefficient for momentum Turbulent closure parameters AKK_BAK Background vertical mixing coefficient for turbulent kinetic energy AKP_BAK Background vertical mixing coefficient for turbulent kinetic energy TKENU2 Constant mixing coefficient for the horizontal Laplacian diffusion of turbulent kinetic energy TKENU4 Constant mixing coefficient for the horizontal biharmonic diffusion of turbu lent kinetic energy Generic length scale turbulence closure parameters GLS_P Stability exponent GLS_M Turbulent kinetic energy exponent GLS_N Turbulent length scale exponent GLS_Kmin Minimum value of specific turbulent kinetic energy GLS_Pmin Minimum value of dissipation GLS_CMUO Stability coefficient GLS_ C1 Shear production coefficient GLS_ C2 Dissipation coefficient GLS_C3M Buoyancy production coefficient minus GLS_
113. ght want to average over the last day of a thirty day run NAVG Number of time steps between writing time averaged data into the averages file NDEFAVG Number of time steps between starting new averages files NTSDIA Starting time step for the accumulation of output diagnostics data For in stance you might want to write diagnostics for the last day of a thirty day run NDIA Number of time steps between writing diagnostics data into the diagnostics file NDEFDIA Number of time steps between starting new diagnostics files Tangent linear and adjoint output parameters LeycleTLM Logical true to cycle between two records NTLM Number of time steps between writing tangent linear data NDEFTLM Number of time steps between starting new tangent linear files LcycleADJ Logical true to cycle between two records of the restart file NADJ Number of time steps between writing adjoint data NDEFADJ Number of time steps between starting new adjoint files NSFF Number of time steps between 4DVAR adjustment of surface forcing fluxes NOBC Number of time steps between 4DVAR adjustment of open boundary fields Check pointing GST restart parameters LrstGST GST restart switch MaxIterGST Maximum number of iterations NGST Check pointing interval Ritz GST parameter Ritz_tol Relative accuracy of the Ritz values computed in the GST analysis Horizontal mixing of tracers TNU2 Constant mixing coefficient for the horizontal Laplacian dif
114. gnu make and Perl All input and output is done in NetCDF Network Common Data Format requires the NetCDF library e Options include serial parallel with MPI and parallel with OpenMP The above list hasn t changed so very much in the past ten to fifteen years but many of the numerical details have changed a great deal Examples include consistent temporal averaging of the barotropic mode to guarantee both exact conservation and constancy preservation properties for tracers redefined barotropic pressure gradient terms to account for local variations in the density field vertical interpolation performed using conservative parabolic splines and higher order quasi monotone advection algorithms ROMS now comes with a full suite of advanced data assimilation routines these options are beyond the scope of this document Chapter 2 has some information on getting started with ROMS Chapters 3 and 4 describe the model physics and numerical techniques and contain information from Shchepetkin and McWilliams and Haidvogel et al 24 Chapter 5 describes the ice equations and Chapter 6 lists the model subroutines and functions As distributed ROMS is ready to run with a number of example problems The process of configuring ROMS for a particular application and running it is described in Chapter 7 including a discussion of a few example applications Chapter 8 describes Hernan Arango s plotting programs cnt cent sec and csec 2
115. h grid orientation define EASTERN_WALL define NORTHERN_WALL undef WESTERN_WALL undef SOUTHERN_WALL define RADIATION_2D ifndef NORTHERN_WALL define NORTH_FSCHAPMAN define NORTH_M2FLATHER ifdef SOLVE3D define NORTH_M3RADIATION define NORTH_M3NUDGING define NORTH_TRADIATION define NORTH_TNUDGING define NORTH_MIGRADIENT endif Hendif ifndef WESTERN_WALL define WEST_FSCHAPMAN define WEST_M2FLATHER ifdef SOLVE3D define WEST_M3RADIATION define WEST_M3NUDGING define WEST_TRADIATION define WEST_TNUDGING define WEST_MIGRADIENT endif Hendif ifndef SOUTHERN_WALL ifdef SOLVE3D define SOUTH_FSCHAPMAN define SOUTH_M2FLATHER define SOUTH_M3RADIATION define SOUTH_M3NUDGING define SOUTH_TRADIATION 108 define SOUTH_TNUDGING define SOUTH_MIGRADIENT endif endif ifndef EASTERN_WALL define EAST_FSCHAPMAN define EAST_M2FLATHER ifdef SOLVE3D define EAST_M3RADIATION define EAST_M3NUDGING define EAST_TRADIATION define EAST_TNUDGING define EAST_MIGRADIENT endif endif roms quirks ifdef SOLVE3D define ANA_BSFLUX define ANA_BTFLUX else define ANA_SMFLUX endif Biological model options define NEMURO define LIMIT_BIO_AKT undef BIO_GOANPZ Sarah Hinckley s 11 box model undef BEST_NPZ Georgina Gibsons BEST NPZ model if defined BEST_NPZ
116. hat ana_ sbflux ana_ stflux etc are set correctly in this case taking the default of zero rather than explicitly setting anything for UPWELLING However we do set ana_ srflux to be non zero in case we opt to turn on one of the biology models 7 2 7 ocean in The model has been set up to run for five days with an internal time step of 300 s and an external time step of 10 s NTIMES 1440 DT 300 0d0 NDTFAST 30 We will write history averages and diagnostics records every 1 4 day restart records once a day NRREC 0 LcycleRST T NRST 288 LDEFOUT T NHIS 72 NDEFHIS 0 NTSAVG 1 NAVG 72 NDEFAVG 0 NTSDIA 1 NDIA 72 NDEFDIA 0 The value of the linear bottom friction coefficient rdrg is set to 3 0 x 1074 and the channel walls are set to be free slip RDRG 3 0d 04 m s GAMMA2 1 0d0 The vertical stretching is set to a modest value of theta__s 3 91 Vtransform 1 transformation equation stretching function surface stretching parameter bottom stretching parameter critical depth m Vstretching 1 THETA_S 3 0d0 THETA_B 0 0d0 TCLINE 50 0d0 7 2 8 Output The model writes some information to standard out after setting ninfo to 72 Process Information Thread 0 pid 32022 is active Model Input Parameters ROMS TOMS version 3 2 Friday October 30 2009 9 11 20 AM Wind Driven Upwelling Downwelling over a Periodic Channel Operating system Li
117. he Functionals directory contains analytical F which conditionally includes code bits for com puting analytic values for a wide variety of fields Many are alternates for reading from NetCDF files especially for idealized problems ana_aiobc Computes open boundary conditions for the ice concentration ana_biology Computes analytic initial conditions for the biology tracers ana_bmflux Computes analytic kinematic bottom momentum flux ana_btflux Computes analytic kinematic bottom flux of tracer type variables ana_cloud Computes analytic cloud fraction ana_diag Computes customized diagnostics 54 ana_fsobc Computes analytic open boundary conditions for the free surface ana_ grid Sets up an analytic grid ana_hiobc Computes open boundary conditions for the ice thickness ana_hmixcoef Computes spatially variable horizontal mixing coefficients ana_hsnobc Computes open boundary conditions for the snow thickness ana_humid Computes analytic atmospheric humidity ana_ice Computes analytic initial conditions for the sea ice ana_ initial Sets up analytic initial conditions for the ocean ana_m2clima Sets up an analytic climatology for the two dimensional momentum ana_m2obc Computes open boundary conditions for the two dimensional momentum ana_m3clima Sets up an analytic climatology for the three dimensional momentum ana_m3obc Computes open boundary conditions for the three dimensional momentum ana_mask Sets up an analytic mask ana_ncep
118. how many grid points to use and can be afforded There are three parameters in ocean in which specify the grid size and one parameter for the number of active tracers 76 Lm Number of finite difference points in Mm Number of finite difference points in n N Number of finite difference points in the vertical NAT Number of active tracers The number of biological tracers is set in the biology in file There are no constraints on these except Lm gt 2 Mm gt 2 N gt 2 and NAT gt 1 Lm and Mm should be at least 3 if the domain is periodic in that direction 7 1 6 x y grid The subroutine get__grid or ana_ grid is called by initial to set the grid arrays the bathymetry and the Coriolis parameter Most of the simple test problems have their grid information specified in ana_ grid h in the directory ROMS Functionals More realistic problems require a NetCDF grid file produced by the grid generation programs described in Wilkin and Hedstrom 89 by the Matlab SeaGrid or by some other method The variables which are read by get__grid are xl el spherical f h pm pn x_rho y_rho lon_rho lat_rho angle If the grid is curved get__grid will also read dndx dmde Likewise if MASKING is defined it will read mask_rho mask_u mask_v mask_ psi 7 1 7 n grid Before providing initial conditions and boundary conditions the user must understand the model grid The fields are laid out on an Arakawa C grid as in Fig 1
119. howing tides This is from a snapshot in a history file the averages files have been detided 114 ROMS 3 2 NEP5 467 00 Day Min 8 1842E 05 Max 2 9995E D1 Fraction of Cell Covered by Ice Figure 23 Ice concentration averaged over the month of April 1959 115 ROMS 3 0 NEP5 487 00 Day I 092 o 424 J m 100 200 200 Min 9 2Y44E 01 Max 64270E 00 Potential Temperature C 400 600 Figure 24 Vertical slice if temperature averaged over the month of April 1959 The slice is across the Bering Sea shelf showing the transition from vertically mixed at the coast a two layer system at mid shelf then a thermocline over the shelf break 116 8 Plotting Programs for Postprocessing Hernan Arango has provided ROMS with some Fortran programs for creating plots from the NetCDF history and restart files Some example plots are shown in Other plotting options are available via Matlab Python and NCL Here we describe only the Fortran programs available via svn from myroms org There are four plotting programs cnt ccnt sec csec creates black and white plots of the horizontal fields including constant depth plots of the 3 D fields creates color plots of the horizontal fields including constant depth plots of the 3 D fields creates black and white plots of vertical slices through the 3 D fields It includes on option of finding equal spaced points along a straight line through the curvi
120. i F at u e 3 OE 56 H m AG rie z Es c a Di 57 a Ay no gt o 1067 aon mn e 3 00 ce 4 8 3 Fourth order Akima An alternate fourth order algorithm is that by Akima OC OC OC OC G 2 I F 59 BE BE iy DET E a4 2 gn 32C OC Z je 60 On on j l On j On j l pS 32C OC OC OC ee Oo k of Z i dera en With the fluxes as in 4 8 4 Third order Upwind There is a class of third order upwind advection schemes both one dimensional Leonard 44 and two dimensional Rasch and Shchepetkin and McWilliams 711 This scheme is known as UTOPIA Uniformly Third Order Polynomial Interpolation Algorithm Applying flux limiters to UTOPIA is explored in Thuburn 83 although it is not implemented in ROMS The two dimensional formulation in Rasch contains terms of order u2C and uC including cross terms uvC The terms which are nonlinear in velocity have been dropped in ROMS leaving one extra upwind term in the computation of the advective fluxes H u RC Fe c vaz 63 H v 3C P ma c y a 64 The second derivative terms are centered on a p point in the grid but are needed at a u or v point in the flux The upstream value is used He zZ Ps E max 0 ui j k Cons min 0 ui jk Ci j k 65 The value of y in the model is t while that in Rasch is Because the third order upwind scheme is designed to be two dimensional it is not used in the vertical though one might argue that
121. ibrary gets included during the final link stage However we are now using the Fortran 90 version of it which requires its module information as well We just copy the mod files into the Build directory NETCDF_MODFILE netcdf mod TYPESIZES_MODFILE typesizes mod SCRATCH_DIR NETCDF_MODFILE SCRATCH_DIR cp f NETCDF_INCDIR NETCDF_MODFILE SCRATCH_DIR SCRATCH_DIR TYPESIZES_MODFILE SCRATCH_DIR cp f NETCDF_INCDIR TYPESIZES_MODFILE SCRATCH_DIR Old versions of NetCDF do not have the typesizes mod file in which case it has to be removed from the following dependency list SCRATCH_DIR MakeDepend makefile SCRATCH_DIR NETCDF_MODFILE SCRATCH_DIR TYPESIZES_MODFILE SCRATCH_DIR Then there is MakeDepend itself This file is created by the Perl script sfmakedepend MakeDepend gets created by make depend and included on make s second pass through the makefile depend SCRATCH_DIR SFMAKEDEPEND MDEPFLAGS sources gt SCRATCH_DIR MakeDepend ifneq MAKECMDGOALS clean include SCRATCH_DIR MakeDepend endif The dash before the include tells make to ignore errors so that make depend will succeed before the file exists The MakeDepend file will contain the include and module dependen cies for each source file such as Build mod_diags o tile h cppdefs h globaldefs h shoreface h Build mod_diags f90 tile h cppdefs h globaldefs h shoreface h Bu
122. ic rheology we are using Paul Budgell has rewritten the dynamic sea ice model improving the solution procedure and making the water stress term implicit in time then changing it again to use the elastic viscous plastic rheology of Hunke and Dukowicz I am very grateful that he is allowing us to use his version of the code The sea ice thermodynamics is derived from Sirpa Hakkinen s implementation of the Mellor Kantha scheme She was kind enough to allow Paul and I to start with her code Thanks to the internet community for providing great tools like Perl patch cpp svn and gmake to aid in software development and to make it more fun This work was supported in part by a grant of HPC resources from the Arctic Region Super computing Center and the DoD High Performance Computing Modernization Program Development and testing of the ROMS model has been funded by many including the USGS Coastal and Marine Program the Office of Naval Research the National Ocean Partnership Pro gram UNIX is a registered trademark of the Open Group Cygwin is a registered trademark of Red Hat Inc Abstract The Regional Ocean Modeling System ROMS authored by many most notably Sasha Shchepetkin is one approach to regional and basin scale ocean modeling This user s manual for ROMS describes the model equations and algorithms as well as additional user configurations necessary for specific applications ROMS itself has now branched out as well the
123. idal constituents classified by period The tidal forcing is computed for the full horizontal grid If requested the tidal forcing is added to the processed open boundary data Both tidal elevation and tidal currents are required to force the model properly However if only the tidal elevation is available the tidal currents at the open boundary can be esti mated by reduced physics Only the pressure gradient Coriolis and surface and bottom stresses terms are considered at the open boundary See u2dbc_im F or v2dbc_im F for details Notice that there is an additional option FSOBC_REDUCED for the computation of the pressure gradient term in both Flather or reduced physics conditions _M2FLATHER _M2REDUCED SSH_TIDES Define if imposing tidal elevation UV_TIDES Define if imposing tidal currents RAMP_ TIDES Define if ramping over one day tidal forcing from zero ADD_FSOBC Define to tidal elevation to processed OBC data ADD_M2OBC Define to tidal currents to processed OBC data TIDES_ASTRO Define to add contributions from the long period tides as done by Foreman Climatology M2CLIMATOLOGY Define for processing the 2 D momentum climatology arrays M38CLIMATOLOGY Define for processing the 3 D momentum climatology arrays OCLIMATOLOGY Define for processing the vertical momentum climatology arrays TCLIMATOLOGY Define for processing the tracer climatology arrays ZCLIMATOLOGY Define for processing the sea surface height climatology arrays M2
124. ild mod_diags o Build mod_kinds o Build mod_param o Build mod_diags f90 Note that the h files are included by cpp so that both the f90 and o files become out of date when an include file is modified Without the module dependencies make would try to build the sources in the wrong order and the compiler would fail with a complaint about not finding mod_param mod for instance 145 G 4 Final warnings The cost of this nifty make stuff is 1 We re a little closer to the gnu make bugs here and we need a newer version of gnu make than before version 3 81 3 80 if you re lucky Hence this stuff at the top of the makefile NEED_VERSION 3 80 3 81 if filter MAKE_VERSION NEED_VERSION error This makefile requires one of GNU make version NEED_VERSION 2 The Makefile dependencies get just a little trickier every change we make Note that F90 has potentially both include and module use dependencies The book example uses the C compiler to produce its own dependencies for each source file into a corresponding d file to be included by make Our Fortran compilers are not so smart For these hairy compiles it s critical to have accurate dependency information unless we re willing to make clean between compiles 146 H sfmakedepend The other Perl script I use with Fortran modifies the Makefile to include dependency information much like the X11 program makedepend l originally wrote fmakedepend which was used
125. in several segments and has been updated and restructured here G 1 Introduction to Portable make Make gets its instructions from a description file by default named makefile or Makefile This file is called the Makefile but other files can be used by invoking make with the f option make f Makefile yukon The ancester to ROMS originally had a Makefile that looked something like model main o init o plot o lt TAB gt 77 o model main o init o plot o main o main f lt TAB gt 77 c O main f init o init f lt TAB gt 77 c O init f plot o plot f lt TAB gt 77 c 00 plot f clean lt TAB gt rm o core The default thing to build is model the first target The syntax is target dependencies lt TAB gt command lt TAB gt command The target model depends on the object files main o and friends They have to exist and be up to date before model s link command can be run If the target is out of date according to the timestamp on the file then the commands will be run Note that the tab is required on the command lines The other targets tell make how to create the object files and that they in turn have dependen cies To compile model simply type make Make will look for the file makefile read it and do whatever is necessary to make model up to date If you edit init f that file will be newer than init o Make would see that init o is out of date and run the f77 c O init f command Now
126. init o is newer than model so the link command f77 o model main o init o plot o must be executed To clean up type make clean and the clean target will be brought up to date The clean target has no dependencies What happens in that case is that the command rm o core will always be executed The original version of this Makefile turned off optimization on plot f due to a compiler bug but hopefully you won t ever have to worry about that 133 G 1 1 Macros Make supports a simple string substitution macro Set it with MY_MACRO nothing today and refer to it with MY_MACRO The convention is to put the macros near the top of your Makefile and to use upper case Also use separate macros for the name of your compiler and the flags it needs F90 f90 F9OFLAGS 03 LIBDIR usr local lib LIBS L LIBDIR Imylib Let s rewrite our Makefile using macros IBM version F90 x1f90_r F90FLAGS 03 qstrict LDFLAGS bmaxdata 0x40000000 model main o init o plot o F90 LDFLAGS o model main o init o plot o main o main f F90 c F9OFLAGS main f init o init f F90 c F9OFLAGS init f plot o plot f F90 c F9OFLAGS plot f clean rm o core Now when we change computers we only have to change the compiler name in one place Likewise if we want to try a different optimization level we only specify that in one place Notice that you can use comments b
127. interior mixing parameteriza tion represents double diffusive mixing From limited sources of laboratory and field data LMD parameterize the salt fingering case R gt 1 0 4 Rp 1 2 3 O oe 1x 10 4 1 GE for 1 0 lt Ry lt Ry 19 ass 0 otherwise va R 0 74 149 For diffusive convection 0 lt Rp lt 1 0 LMD suggest several formulations from the literature and choose the one with the most significant impact on mixing Fedorov Je Vy 1 579 0 909 exp 4 6 exp 0 54 R 1 150 for temperature For other scalars r e 0 85R5 R for 0 5 lt Rp lt 1 0 AR j v0 15 R otherwise Internal wave generated mixing Internal wave generated mixing serves as the background mixing in the LMD scheme It is specified as a constant for both scalars and momentum Eddy diffusivity is estimated based on the data of Ledwell et al 43 While Peters et al 62 suggest eddy viscosity should be 7 to 10 times larger than diffusivity for gradient Richardson numbers below approximately 0 7 Therefore LMD use y 10x10 m s gt 152 y 1 0 x 10 m s 153 4 12 Timestepping vertical viscosity and diffusion The Da Dy and Dg terms in equations 18 20 represent both horizontal and vertical mixing processes The horizontal options will be covered in 4 10 1 The model has several options for computing the vertical coefficients these will be described in 4 11 The vertical viscosity and diffusion terms ha
128. ion of ice blocks in which no ridging occurs as long as there is any open water An optional ridging term can be added Gray and Killworth 21 OA O uA O vA AQ A V vH V 0 SA DA 0 lt AK lt 1 192 Ot Ox Oy where a A is an arbitrary function such that a 0 0 a 1 1 and 0 lt a A lt 1 The ridging term leads to an increase in h under convergent flow as would be produced by ridging The function a A should be chosen so that it is near zero until the ice concentration is large enough that ridging is expected to occur then should increase smoothly to one The symbols used in these equations along with the values for the constants are listed in Table a Note that Hibler s hy variable is equivalent to our Ah combination his Ay is the average thickness over the whole gridbox while our h is the average thickness over the ice covered fraction of the gridbox 5 2 Thermodynamics The thermodynamics is based on calculating how much ice grows and melts on each of the surface bottom and sides of the ice floes as well as frazil ice formation Mellor and Kantha 55 Once the ice tracers are advected the ice concentration and thickness are timestepped according to the terms on the right 189 and is DAh _ Po Dt Pi A Wio Wai 1 gt A Wao Wir 193 DA PA ao r lt aL 194 Di eli O 1 A W 1 A Wy O lt A 9 The term Ah is the effective thickness a measure of the ice volume Its evo
129. itialized with dense water at one end and light water at the other At time zero the water is released and it generates two propagating fronts as the light water rushes to fill the top and the dense water rushes to fill the bottom This configuration was used to test various advection schemes OVERFLOW This configuration is similar to the GRAV_ADJ problem but is initialized with dense water in the shallow part of a domain with a sloping bottom SEAMOUNT The seamount test was used to test the pressure gradient errors It has an idealized seamount in a periodic channel See Beckmann and Haidvogel 4 and McCalpin for more information UPWELLING The upwelling downwelling example was contributed by Anthony Macks and Ja son Middleton 47 and consists of a periodic channel with shelves on each side There is along channel wind forcing and the Coriolis term leads to upwelling on one side and downwelling on the other side If you run it for several days without vertical mixing you end up with dense water over light water Some NetCDF input files for the ROMS examples can be found under Data ROMS in the ROMS distribution The ASCII input files are under ROMS External 7 1 Configuring ROMS The four main sections you need to change in ROMS are the makefile or build bash an include file with cpp options any analytic functions and the ASCII input file If more realistic fields are desired you will have to provide other input files as well for inst
130. ity te out tracer 01 temp te out tracer 02 salt ocean_rst nc ocean_his nc ocean_avg nc ocean_dia nc ROMS External varinfo dat Grid 01 0041x0080x0016 tiling 001x001 Jstr Jend Npts 1 80 52480 Tile minimum and maximum fractional grid coordinates interior points only tile Xmin Xmax 0 2 50 43 50 0 2 50 43 50 0 2 50 43 50 Ymin Ymax grid 0 50 82 50 RHO points 0 50 82 50 U points 0 50 82 50 V points Activated C preprocessing Options UPWELLING Wind Driven Upwelling Downwelling over a Periodic Channel ANA_BSFLUX Analytical kinematic bottom salinity flux ANA_BTFLUX Analytical kinematic bottom temperature flux ANA_GRID Analytical grid set up ANA_INITIAL Analytical initial conditions ANA_SMFLUX Analytical kinematic surface momentum flux 94 ANA_SSFLUX ANA_STFLUX ANA_VMIX ASSUMED_SHAPE AVERAGES DIAGNOSTICS_TS DIAGNOSTICS_UV DJ_GRADPS DOUBLE_PRECISION EW_PERIODIC MIX_S_TS MIX_S_UV NONLINEAR INONLIN_EOS POWER_LAW PROFILE IRST_SINGLE SALINITY SOLVE3D SPLINES TS_U3HADVECTION TS_C4VADVECTION TS_DIF2 UV_ADV UV_COR UV_U3HADVECTION UV_C4VADVECTION UV_LDRAG UV_VIS2 VAR_RHO_2D Analytical kinematic surface salinity flux Analytical kinematic surface temperature flux Analytical vertical mixing coefficients Using assumed shape arrays Writing out time averaged fields Computing and writing tracer diagnostic terms Computing and writing momentum diagnostic terms Parabolic Spli
131. ke supports an assignment which avoids the need for FFLAGS2 entirely 136 lmd_skpp o FFLAGS 0 inlinefrom lmd_wscale f90 What this means is that for the target lmd__skpp o only append the inlining directive to FFLAGS I think this is pretty cool One last kind of assignment is to set the value only if there is no value from somewhere else the environment for instance FC mpxlf90_r If we used or we would override the value from the environment G 2 3 Include and a Few Functions One reasonably portable make feature is the include directive It can be used to clean up the Makefile by putting bits in an include file The syntax is simply include file and we use it liberally to keep the project information neat We can also include a file with the system compiler information in it assuming we have some way of deciding which file to include We can use uname s to find out which operating system we re using We also need to know which compiler we re using One user defined variable is called FORT the name of the Fortran compiler This value is combined with the result of uname s to provide a machine and compiler combination For instance ftn on Linux is the Cray cross compiler This would link to a different copy of the NetCDF library and use different compiler flags than the Intel compiler which might be on the same system The user sets FORT FORT ftn OS shell uname s sed s g
132. l terms including v terms in the u equation and vice versa These extra terms were found to be small in a test problem and have been left out of the model 4 11 Vertical mixing schemes ROMS contains a variety of methods for setting the vertical viscous and diffusive coefficients The choices range from simply choosing fixed values to the KPP generic lengthscale GLS and Mellor Yamada turbulence closure schemes See Large 40 for a review of surface ocean mixing schemes Many schemes have a background molecular value which is used when the turbulent processes are assumed to be small such as in the interior 28 4 11 1 Mellor Yamada One of the more popular closure schemes is that of Mellor and Yamada 56 57 They actually present a hierarchy of closures of increasing complexity ROMS provides only the Level 2 5 closure with the Galperin et al modifications as described in Allen et al i This closure scheme adds two prognostic equations one for the turbulent kinetic energy 5d and one for the turbulent kinetic energy times a length scale q 1 The turbulent kinetic energy equation is D o q 7 5 Oz Ks 5 F Py a 103 where P is the shear production P is the buoyant production and q is the dissipation of turbulent kinetic energy These terms are given by Bo 163 TES 110 P KN 111 3 q gt 112 bd Bil 112 where B is a constant One can also add a traditional horizontal Lapl
133. l mixing when operating along constant s surfaces The transverse stress tensor rotated along geopotentials constant depth is then given by where 27 o o ronal e Mn A Blt 2 988 06m 260 5 nau nfs pO 26 2 ae nD qo ie E nv mA Lait oe noe O m no 1 O mu n Am n No ot Tag H Os On H s i nee a a ofA je 02 1 Oz 1 O nv O nv 1 Oz 1 O mu O mu AN Ga a m IO MIE Of an 105 Oz 1 Oz 1 O mv O mv 1 Oz 1 O nu O nu Mr s on n mz ds aE or Notice that transverse stress tensor remains invariant under coordinate transformation The rotated tensor and retains the same properties as the unrotated tensor and 92 The additional terms that arise from the slopes of s surfaces along geopotentials are discretized using a modified version of the triad approach of Griffies et al 1998 22 4 10 4 Biharmonic Biharmonic mixing has less of a physical justification and is used for damping of numerical noise at the 2Ax scale The biharmonic operator is V4 V2V the corresponding term is computed using a temporary variable Y _ mn 1v4H m ow O v H n Ow Oe sel n de j A m a 107 oe O H m OY O Hyn OY zm on sel n 0 oad m m we where w is any of u v T and S Note that u and v are treated as independent scalar quantities rather than as a vector The complete Laplacian operator on a vector quantity u contains additiona
134. layer There are three ways to provide longwave radiation in the at mospheric boundary layer 1 Compute the net longwave radiation internally using the Berliand 1952 equation LONGWAVE as function of air temperature sea surface temperature relative humidity and cloud fraction 2 provide read longwave down welling radiation only and then add outgoing longwave radiation LONGWAVE_ OUT as a function of the model sea surface temperature 3 provide net longwave radiation default The shortwave radiation can be computed using the global albedo equation with a cloud correction Alternatively input shortwave radiation data computed from averaged data with snapshots greater or equal than 24 hours can be modulated by the local diurnal cycle which is a function longitude latitude and day of year BULK_ FLUXES Define for bulk flux computation CCSM_FLUXES Define for CCSM version of bulk flux computation NCEP_FLUXES Define if NCEP forcing files are used NL_BULK_FLUXES Define to use bulk fluxes computed by nonlinear forward model COOL_ SKIN Define for cool skin correction LONGWAVE Define to compute net longwave radiation LONGWAVE_ OUT Define to compute outgoing longwave radiation EMINUSP Define to compute evaporation minus precipitation COARE_TAYLOR_YELLAND Define to use Taylor and Yelland wave roughness 182 COARE_OOST Define to use Oost et al wave roughness 58 DEEPWATER_ WAVES Define to use deep water wave
135. le A tiled grid with some ROMS tile variables Figure 13 a o H o l H n Y A Z O m o a D E lan yo g lan mod_param F MarginI and MarginJ are zero if the numbers work out perfectly i e Lm NtileI and Mm NtileJ are integers The tile numbers match the MPI process numbers Some internal ROMS numbers are shown in Fig Each tile can be numbered from 1 to Chunksizel or it can retain the numbering it n o Pa 45 fan a0 z a0 a0 l O o yo g o a o Nn o H amp g ie del a o yo g Q oO o Nn 2 a me o Nn ie a o S as a gt E a0 o a ie a o a ares a S TS a g m ae 5 B F Fan 5 Ey Jo g 2 O yo Q o ed AQ z a7 2 pi as a a0 o aq a O n g O E ia o o r o E o Q n O p gt Q A fav a G a ES n Aa E sy Oo A pin n av a o n 3 S O F 3 Pa ee o z D 8 2 g 3 S S 3 In picking a numbering scheme for indices within a tile there are two common choices as shown in F indices for each tile Some of ing d end mning an With the tile sizes known we can assign beg the details depend on whether or not the doma ap z as shown in F irection in that di ic iod is per n 6 9 2 MPI exchange For MPI jobs the ghost points need to be updated betwe
136. linear grid creates color plots of vertical slices through the 3 D fields It includes on option of finding equal spaced points along a straight line through the curvilinear grid All of these program come with example input files For instance the input file for cent is called ccnt in and is as follows 1996 1 year and starting year day use yearday lt 0 for no time label ROMS 3 2 Coarse Arctic ocean with Budgell ice dynamics ice thermodynamics 8 NFIELDS number of fields to plot Line below field s types 42 45 46 47 48 49 50 121 122 123 124 field identification FLDID 1 NFIELDS 1 20 2 0 0 2 2 1 1 1 0 0 0 0 0 0 0 1 0 1 2 0 1 1 0 25 1 0 0 0 0 NLEVELS number of levels and or depths to plot 0 for all levels levels gt 0 or depths lt 0 to plot FLDLEV 1 NLEVELS NFIELDS number of fields to plot Line below field s types field identification FLDID 1 NFIELDS NLEVELS number of levels and or depths to plot 0 for all levels 2 3 4 5 levels gt 0 or depths lt 0 to plot FLDLEV 1 NLEVELS FRSTD first day to plot LASTD last day to plot DSKIP plot every other DSKIP days 0 0 plot at its own time frequency VINTRP vertical interpolation scheme O linear 1 cubic splines PMIN field minimum value for color palette 0 0 for default PMAX field maximum value for color palette 0 0 for default ICNT draw contours between color bands O no 1
137. linity etc Hout idUsms Surface u stress Hout idVsms Surface v stress Hout idUbms Bottom u stress Hout idVbms Bottom v stress Hout idUbrs Bottom u current stress Hout idVbrs Bottom v current stress Hout idUbws Bottom u wave stress Hout idVbws Bottom v wave stress Hout idUbcs Bottom max wave current u stress Hout idVbcs Bottom max wave current v stress Hout idUbot Bed wave orbital u velocity Hout idVbot Bed wave orbital v velocity Hout idUbur Bottom max wave current u stress Hout idVbur Bottom max wave current v stress Hout idW2xx 2D radiation stress Sss component Hout idW2xy 2D radiation stress Sy component Hout idW2yy 2D radiation stress Syy component Hout idU2rs 2D u radiation stress Hout idV2rs 2D v radiation stress Hout idU2Sd 2D u Stokes velocity Hout idV2Sd 2D v Stokes velocity Hout idW3xx 3D radiation stress Sss component Hout idW3xy 3D radiation stress Szy component Hout idW3yy 3D radiation stress Syy component Hout idW3zx 3D radiation stress Sz component Hout idW3zy 3D radiation stress Szy component Hout idU3rs 3D u radiation stress Hout idV3rs 3D v radiation stress Hout idU3Sd 3D u Stokes velocity Hout idV3Sd 3D v Stokes velocity Hout idWamp Wave height 85 Hout idWlen Wave length Hout idWdir Wave direction Hout idTsur Surface net heat and salt flux NAT values Hout idLhea Latent heat flux Hout idShea Sensible heat flux Hout
138. logy cases none of which have been defined Likewise we are using the default of ANA__VMIX as distributed We are asking for many other analytic functions too including the grid We are asking for diagnostic output with the DIAGNOSTICS_ TS and DIAGNOSTICS_UV 7 2 2 Model domain The flow does not vary in x so Lm can be small Set the values for Lm Mm N and NT in the input file Lm 41 Number of I direction INTERIOR RHO points Mm 80 Number of J direction INTERIOR RHO points N 16 Number of vertical levels NAT 2 Number of active tracers usually 2 7 2 3 ana_ grid For this geometry one has a choice of using one of the external grid generation programs or of using ana_ grid to create the grid analytically The code in ana_ grid h was modified to produce a bathymetry with a shelf on both walls of the channel when UPWELLING is defined The fluid depth ranges from 27 m on the shelves to 150 m in the center of the channel The horizontal grid spacing is uniform at 1 km and the Coriolis parameter f is set to a constant value suitable for Sydney Australia 7 2 4 Initial conditions and the equation of state We would like the initial conditions to be a motionless fluid with an exponential stratification The UPWELLING section of ana_initial h is configured accordingly The stratification can be provided by either T or S or by both T and S For simplicity we will only have an active temperature field and we will use the linea
139. lso has a major responsibility for American Indian reservation communities and for people who live in island territories under U S administration The Minerals Management Service Mission As a bureau of the Department of the Interior the Minerals Management Service s MMS primary responsibilities are to manage the mineral resources located on the Nation s Outer Continental Shelf OCS collect revenue from the Federal OCS and onshore Federal and Indian lands and distribute those revenues Moreover in working to meet its responsibilities the Offshore Minerals Management Program administers the OCS competitive leasing program and oversees the safe and environmentally sound exploration and production of our Nation s offshore natural gas oil and other mineral resources The MMS Minerals Revenue Management Program meets its responsibilities by ensuring the efficient timely and accurate collection and disbursement of revenue from mineral leasing and production due to Indian tribes and allottees States and the U S Treasury The MMS strives to fulfill its responsibilities through the general guiding principles of 1 being responsive to the public s concerns and interests by maintaining a dialogue with all potentially affected parties and 2 carrying out its programs with an emphasis on working to enhance the quality of life for all Americans by lending MMS assistance and expertise to economic development and environmental protection
140. lution equation is simply quantifying the change in the amount of ice The ice concentration equation is more interesting in that it provides the partitioning between ice melt growth on the sides vs on the top and bottom The parameter controls this and has differing values for ice melt and retreat In principle most of the ice growth is assumed to happen at the base of the ice while rather more of the melt happens on the sides of the ice due to warming of the water in the leads The heat fluxes through the ice are based on a simple one layer Semtner 70 type model with snow on top The temperature is assumed to be linear within the snow and within the ice 38 Variable Value Description A x y t ice concentration a A ridging function Ca nonlinear air drag coefficient Ca 2 2 x 1078 air drag coefficient Cw 10 x 1073 water drag coefficient Dr Ds DA diffusion terms bi Kronecker delta function E Young s modulus e 2 eccentricity of the elliptical yield curve eij x y t strain rate tensor n x y t nonlinear shear viscosity f x y Coriolis parameter Esa internal ice stress g 9 8 m s7 acceleration of gravity H Heaviside function hi x y t ice thickness of ice covered fraction ho 1m ice cutoff thickness h x 054 snow thickness on ice covered fraction M x y t ice mass density times thickness P x y t ice pressure or strength P C 2 75 x 104 20 ice strength parameters Sh Ss SA thermodynamic terms Cle y
141. lux profiles are expressed as analytical fits to atmospheric surface boundary layer data In stable conditions they vary linearly with the stability parameter as bc 14 5 131 In near neutral unstable conditions common Businger Dyer forms are used which match with the formulation for stable conditions at 0 Near neutral conditions are defined as 0 2 lt lt 0 132 31 for momentum and 10 lt lt 0 133 for scalars T he non dimensional flux profiles in this regime are m 1 164 134 gs 1 164 135 In more unstable conditions is chosen to match the Businger Dyer forms and with the free convective limit Here the flux profiles are m 1 26 8 38 136 bs 28 86 98 96 1 3 137 The shape function The non dimensional shape function G c is a third order polynomial with coefficients chosen to match the interior viscosity at the bottom of the boundary layer and Monin Obukhov similarity theory approaching the surface This function is defined as a 3rd order polynomial G 0 ao ayo ago ago 138 with the coefficients specified to match surface boundary conditions and to smoothly blend with the interior de 0 139 ay 1 140 O Ve Rsot OxValh Vvalh swze l MA O a hw2 1 eu Ve Asot OxValh vye h Oswx 1 A a hw2 1 Soe where v h is the viscosity calculated by the interior parameterization at the boundary layer depth Countergradient flux
142. mic interpolation of Ur Vr SSW_ FORM_DRAG_ COR Define to activate form drag coefficient SSW_ZOBIO Define for biogenic bedform roughness for ripples SSW _ZOBL Define for bedload roughness for ripples SSW _ZORIP Define for bedform roughness for ripples ICE_MODEL Define to use ice component of the model see 5 ICE_THERMO Define for ice thermodynamics ICE_MK Define for Mellor Kantha 55 ice thermodynamics this is the only choice ICE_SMOOTH Define to smooth some ice fields ICE_ALB_EC92 Define for albedo computations from Ebert and Curry 10 ICE_MOMENTUM Define for momentum component of the ice ICE_MOM_ BULK Define for alternate ice water stress computation ICE_EVP Define for elastic viscous plastic rheology 33 and 32 ICE_ADVECT Define for advection of ice tracers ICE_SMOLAR Define to use MPDATA for ice tracers no other option ICE_UPWIND Define for upwind advection not available ICE_BULK_FLUXES Define for ice part of bulk flux computation Boundary conditions EW_PERIODIC Define for periodic boundaries in the direction SPONGE Define to use SSH climatology as 2D inflow data NS_PERIODIC Define for periodic boundaries in the 7 direction RADIATION_ 2D Define for tangential phase speed in radiation conditions Detailed eastern open boundary conditions Other sides have similar If none of these are defined it is assumed to be a closed wall EAST_FSCHAPMAN Define for a Chapman condition on the fr
143. momentum advection need to be activated except for UV_ADV The 3rd order upstream split advection UV_U3ADV_ SPLIT can be used to correct for the spurious mixing of the advection operator in terrain following coordinates If this is the case the advection operator is split in advective and viscosity components and several internal flags are activated in globaldefs h Notice that horizontal and vertical advection of momentum is 4th order centered plus biharmonic viscosity to correct for spurious mixing UV_ADV Define to compute the momentum advection terms CURVGRID Define to compute the extra non linear advection terms which arise when using curvilinear coordinates UV_COR Define to compute the Coriolis term UV_U3ADV_ SPLIT Define for 3rd order upstream split momentum advection UV_C2ADVECTION Define for 2nd order centered advection 56 UV_C4ADVECTION Define for 4rd order centered advection UV_SADVECTION Define for splines vertical advection for shallow vertically well resolved domains UV_VIS2 Define to compute the horizontal Laplacian viscosity UV_VIS4 Define to compute the horizontal biharmonic viscosity UV_SMAGORINSKY Define for Smagorinsky like viscosity UV_LOGDRAG Define for logarithmic bottom friction UV_LDRAG Define for linear bottom friction UV_QDRAG Define for quadratic bottom friction RDRG_GRID Define for spatially variable bottom drag DRAG_LIMITER Define for bottom drag limiter UV_PSOURCE Define for point
144. mp assigned to model initialization days Reference time for units attribute yyyymmdd dd Nudging relaxation time scale days for tracer 01 temp Nudging relaxation time scale days for tracer 02 salt Nudging relaxation time scale days for free surface Nudging relaxation time scale days for 2D momentum Nudging relaxation time scale days for 3D momentum Factor between passive and active open boundary conditions Background potential temperature C constant Background salinity PSU constant 93 1027 000 RO 1 7000E 04 Tcoef 0 0000E 00 Scoef 1 000 gamma2 Hout idFsur Hout idUbar Hout idVbar Hout idUvel Hout idVvel Hout idWvel Hout idOvel Hout idTvar Hout idTvar HHA HAHAHA AH Output Input Files Output Restart Output History Output Averages Output Diagnostics IO Variable Information Tile partition information tile Istr Tend Number of tracers 0 1 41 Bac o The Sal Sli Wri Wri Wri Wri Wri Wri Wri Wri Wri File File File File File for kground density kg m3 used in linear Equation f State rmal expansion coefficient 1 Celsius ine contraction coefficient 1 PSU pperiness variable free slip 1 0 or no slip 1 0 te out free surface te out 2D U momentum component te out 2D V momentum component te out 3D U momentum component te out 3D V momentum component te out W momentum component te out omega vertical veloc
145. n KPP MY25_ MIXING Define to activate Mellor Yamada Level 2 5 closure KANTHA_ CLAYSON Define for Kantha and Clayson stability function K C2ADVECTION Define for 2th order centered advection K C4ADVECTION Define for 4th order centered advection N2S2_ HORAVG Define for horizontal smoothing of buoyancy shear PP_MIXING Define to activate Pacanowski Philander closure RI_HORAVG Define for horizontal Richardson number smoothing RI_VERAVG Define for vertical Richardson number smoothing L L L L L L Bottom boundary layer The Options MB__ZOBL and MB__ZORIP should be activated con currently MB_BBL Define to activate Meinte Blaas BBL closure MB_CALC_ZNOT Define to compute bottom roughness internally MB_CALC_ UB Define to compute bottom orbital velocity internally MB_ZOBIO Define for biogenic bedform roughness for ripples MB_ZOBL Define for bedload roughness for ripples MB_ZORIP Define for bedform roughness for ripples SG_BBL Define to activate Styles Glenn bottom boundary layer formulation SG_CALC_ZNOT Define to compute bottom roughness internally SG_CALC_ UB Define to compute bottom orbital velocity internally SG_LOGINT Define for logarithmic interpolation of Ur Vr SSW_BBL Define to activate Sherwood Signell Warner bottom boundary layer closure SSW_CALC_ZNOT Define to compute bottom roughness internally 61 Sea ice SSW_CALC_ UB Define to compute bottom orbital velocity internally SSW_LOGINT Define for logarith
146. n the absence of internal sources and fluxes through the boundary Both properties are valuable and should be retained when constructing numerical ocean models The semi discrete form of the tracer equation is HE FO O H C uH C vH C H A 1 0 Km0C t mn ngs Js m J e mn mn do 0 tPotFe ey Here d y and 6 denote simple centered finite difference approximations to 0 0 0 01 and 0 00 with the differences taken over the distances A An and Ag respectively Az is the vertical distance from one p point to another CY CY and CY represent averages taken over the distances A An and Ao The finite volume version of the same equation is no different except that a quantity C is defined as the volume averaged concentration over the grid box AV mn H C 5 n do 28 C A Jay mn att The quantity EE is the flux through an interface between adjacent grid boxes n 17 This method of averaging was chosen because it internally conserves first moments in the model domain although it is still possible to exchange mass and energy through the open boundaries The method is similar to that used in Arakawa and Lamb 2 though their scheme also conserves enstrophy Instead we will focus on nearly retaining constancy preservation while coupling the barotropic depth integrated equations and the baroclinic equations The time step in eq is assumed to be from time n to time n 1 with the
147. nces A An to the actual physical arc lengths It is helpful to write the equations in vector notation and to use the formulas for div grad and curl in curvilinear coordinates see Batchelor Appendix 2 3 _ 09 00 Vo Emo Ma 252 n o a o fb V d mn FOO 253 Vxa mn e 5 2 254 a b e AA mdd 0 n0 Vo v Vo mn Te demon a where is a scalar and d is a vector with components a b and c 125 D Viscosity and Diffusion D 1 Horizontal viscosity The horizontal viscosity and diffusion coefficients are scalars which are read in from ocean in Several factors to consider when choosing these values are spindown time The spindown time on wavenumber k is Fn for the Laplacian operator and g for the biharmonic operator The smallest wavenumber corresponds to the length 2A and is k x leading to Az Azt At lt tdamp 23 7 or q TV Ttv 256 This time should be short enough to damp out the numerical noise which is being gen erated but long enough on the larger scales to retain the features you are interested in This time should also be resolved by the model timestep boundary layer thickness The western boundary layer has a thickness proportional to Agg Laz 2 d f 257 for the Laplacian and biharmonic viscosity respectively We have found that the model typically requires the boundary layer to be resolved with at least one grid cell This leads to coarse grids
148. nd 2 4 2 Local File Options BINDIR Directory in which to place the binary executable The default is gt the current top directory SCRATCH_ DIR Put the f90 and the temporary binary files in a build directory to avoid clutter The default is Build under the top directory It can also point to differing places if you want to keep these files for multiple projects at the same time each in their own directory 2 4 2 Providing the Environment Before compiling you will need to find out some background information e What is the name of your compiler e What is returned by uname s on your system e Is there a working NetCDF library e Where is it e Was it built with the above compiler e Do you have access to MPI or OpenMP As described more fully in the makefile will be looking for a file in the Compilers directory with the combination of your operating system and your compiler For instance using Linux and the Pathscale compiler the file would be called Linux path mk Is the corresponding file for your system and compiler in the Compilers directory If not you will have to create it following the existing examples there Next there are two ways to provide the location for the NetCDF files and optional HDF5 library One is by editing the corresponding lines in your system compiler file Another way is through the Unix environment variables If you are always going to be using the same compiler on each system you can edit yo
149. nd the compiler A note on directories make uses vpath to find the source file where it resides It would be possible to compile from the top directory and put the o file in Build with the appropriate arguments but I don t know how to get the mod file into Build short of a mv command Likewise if we compile in the top directory we need to know the compiler option to tell it to look in Build for the mod files it uses Doing a cd to Build before compiling is just simpler 6 The last compile rules is given a list of sources then calls one compile rule once per source file Again you can invoke make p to see how make internally transforms all this into actual targets and rules G 3 4 Library Module mk In each library directory there is a file called Module mk which gets included by the top level makefile These Module mk bits build onto the list of sources and libraries to be compiled and built respectively These Module mk files look something like local_sub ROMS Nonlinear local_lib libNLM a local_src wildcard local_sub x F eval call make library local_lib local_src eval compile rules First we provide the name of the current directory and the library to be built from the resident sources Next we use the wildcard function to search the subdirectory for these sources Note that every F file found will be compiled If you have half baked files that you don t want used make sure they hav
150. nes density Jacobian Shchepetkin 2002 Double precision arithmetic East West periodic boundaries Mixing of tracers along constant S surfaces Mixing of momentum along constant S surfaces Nonlinear Model Linear Equation of State for seawater Power law shape time averaging barotropic filter Time profiling activated Double precision fields in restart NetCDF file Using salinity Solving 3D Primitive Equations Conservative parabolic spline reconstruction Third order upstream horizontal advection of tracers Fourth order centered vertical advection of tracers Harmonic mixing of tracers Advection of momentun Coriolis term Third order upstream horizontal advection of 3D momentum Fourth order centered vertical advection of momentum Linear bottom stress Harmonic mixing of momentum Variable density barotropic mode INITIAL Configuring and initializing forward nonlinear model Vertical S coordinate System level S coord 16 0 0000000 15 0 0625000 14 0 1250000 13 0 1875000 12 0 2500000 11 0 3125000 10 0 3750000 9 0 4375000 8 0 5000000 7 0 5625000 6 0 6250000 5 0 6875000 Cs curve at_hmin over_slope at_hmax 0 0000000 0 000 0 000 0 000 0 0188264 1 575 2 750 3 925 0 0383166 3 150 5 541 7 932 0 0591578 4 725 8 417 12 108 0 0820849 6 300 11 422 16 544 0 1079063 7 875 14 608 21 342 0 1375324 9 450 18 032 26 614 0 1720078 11 025 21 758 32 492 0 2125480 12 6
151. ng schemes are also listed in Table and described in detail in and 75 the relevant ones are described in Appendix A 15 a Barotropic Mode m M m M f AB3 type forward extrapolation ATA AA AA a 5 PR s 5 0 V rhs u v U V Baroclinic Mode AB3 step for u v Tu v with saving T r h s terms n 1 2 1 U V T S n 1 2 ny2 1 b Barotropic Mode m M m M fe id ee b Be 7 S 0 V ths u v T Baroclinic Mode gt 28 U U np 1 u v gt 20 g TS n 1 2 nly 2 1 c a Barotropic Barotropic Mode Mode rhs2D u v rhs2D u u dd Jrhs u v T 5 UV Baroclinic Baroclinic Mode Mode gt 26 u v AB3 type forward u v extrapolation ni2 1 for u v T gt 26 e T S T S 2 AM4 type 3 y Tinterpolation n 1 2 n2 1 for u v Figure 4 Diagrams of the time stepping and mode coupling used in various ROMS versions a Rutgers University ROMS from myroms org b ROMS AGRIF c UCLA ROMS described in 73 d non hydrostatic ROMS 35 In all the curved arrows update the 3 D fields those with pillars are leapfrog in nature with the pillar representing the r h s terms Straight arrows indicate exchange between the barotropic and baroclinic modes The shape functions for the fast time steps show just one option out of many possibilities The grey function has weights to produce an estimate at time n
152. ns only two records with the older record being overwritten during the next write The history file can contain a subset of the restart fields for instance just the surface elevation and the surface temperature The averages file contains time averages of the model fields for instance daily or monthly means depending on NAVG The station file contains timeseries for specified points possibly quite frequently since each record is small For some machinery is in place to write multiple files numbering them _ 0001 _ 0002 etc NRREC Record number of the restart file to read as the initial conditions Set to 0 at the beginning of the run 1 to read the latest record LcycleRST Logical true to cycle between two records of the restart file NRST Number of time steps between writing of restart fields NSTA Number of time steps between writing fields into the stations file NFLT Number of time steps between writing fields into the floats file NINFO Number of time steps between calling diag to write some global information and check for NaN values History average diagnostic output parameters LDEFOUT True for creating new output files for stations history floats etc If false output is appended to these files 80 NHIS Number of time steps between writing history records NDEFHIS Number of time steps between starting new history files NTSAVG Starting time step for the accumulation of output time averaged data For instance you mi
153. nsistent formulation of the anisotropic stress tensor for use in models of the large scale ocean circulation J Comp Phys 105 333 338 1993 87 J C Warner C R Sherwood R P Signell C K Harris and H G Arango Development of a three dimensional regional coupled wave current and sediment transport model Computers amp Geosciences 34 1284 1306 2008 88 D J Webb B A De Cuevas and C S Richmond Improved advection schemes for ocean models J Atmos Ocean Tech 15 1171 1187 1998 89 J Wilkin and K Hedstrom User s manual for an orthogonal curvilinear grid generation package Institute for Naval Oceanography 1991 158 U S Department of the Interior MMS Minerals Management Service Alaska OCS Region The Department of the Interior Mission As the Nation s principal conservation agency the Department of the Interior has responsibility for most of our nationally owned public lands and natural resources This includes fostering sound use of our land and water resources protecting our fish wildlife and biological diversity preserving the environmental and cultural values of our national parks and historical places and providing for the enjoyment of life through outdoor recreation The Department assesses our energy and mineral resources and works to ensure that their development is in the best interests of all our people by encouraging stewardship and citizen participation in their care The Department a
154. nt in freezing equation n 7 59 x 1074 constant in freezing equation S1 salinity before freezing S2 salinity after freezing T temperature before freezing To temperature after freezing Table 7 Frazil ice variables Once we have a the value for Fr we can use it to find the ice growth rates 1 Wio mE 2 Fr 215 Wao Qu Fr 216 217 where Lo E To 1 E T r1 218 The ocean model receives the following heat and salt fluxes Fr AQio 1 A Qao _ WoLo 219 Fg Wo AWro Si So 1 A S0 P E 220 Wo AWio 1 AJWao 221 5 2 2 Frazil ice formation Following Steele et al 80 we check to see if any of the ocean temperatures are below freezing at the end of each timestep If so frazil ice is assumed to form changing the local temperature and salinity The ice that forms is assumed to instantly float up to the surface and add to the ice layer there We assume balances in the mass heat and salt before and after the ice is formed May May Mi 222 Man Cpw Ti L Maia Cow La L miCpiLa 223 Mwy S1 Miro 92 224 The variables are defined in Table 7 Defining y m my and dropping terms of order y leads to L Cpi pw pw So S1 1 y 226 44 We also want the final temperature and salinity to be on the freezing line which we approximate as Tp MS nz 227 We can then solve for y Ti mS nz q Cpi LCpw T
155. nternal ocean temperature ur friction velocity Tio 7p ya 20 roughness parameter Table 6 Ocean surface variables The difference between Qro and Qro goes into the enthalpy of the ice OE pihi 0 VE Qr Qr 208 We can use the chain rule to obtain an equation for timestepping Ti E OT Gye SS a y y o 209 Pi ES U Qr Qr 209 where To T Qro Qra i o T1 PG 2k T3 2 T a EE 2 C Ti hi 1 Ck 5 2 1 Ocean surface boundary conditions The ocean receives surface stresses from both the atmosphere and the ice according to the ice concentration OUw A 1 A K Ta 210 m Oz Po Po ao OVy A 1 A K iY TY 211 m Oz Po io r Po ao where the relevant variables are in table 6 The surface ocean is assumed to be at the freezing temperature for the surface salinity To mS where we use the salinity from the uppermost model point at z 5Az From this we can obtain a vertical temperature gradient for the upper ocean to use in the heat flux formula Fr z gt 0 212 og ort T 212 where C is 213 T2 Pak ln z z0 Br y 1 2 Br b Pr2 3 214 v 43 Variable Value Definition Chi 1994 J kg7 K specific heat of ice Cow 3987 J kg KT specific heat of water y Mi Megs fraction of water that froze L 3 16e5 J kg latent heat of fusion mi mass of ice formed Mar mass of water before freezing Mwz mass of water after freezing m 0 0543 consta
156. nts This is appropriate for the v points E and L in Fig 3 since the flow in and out of the land should be zero It is likewise appropriate for the u point at I but is not necessarily correct for point G The only term in the u equation that requires the u value at point G is the horizontal viscosity which has a term of the form eee Since point G is used in this term by both points A and M it is not sufficient to replace its value with that of the image point for A Instead the term oe is computed and the values at points D and K are replaced with the values appropriate for either free slip or no slip boundary conditions Likewise the term fv in the v equation must be corrected at the mask boundaries This is accomplished by having a fourth mask array defined at the y points in which the values are set to be no slip in metrics For no slip boundaries we count on the values inside the land point G having been zeroed out For point D the image point at G should contain minus the value of u at point A The desired value of 2u is therefore 2ua while instead we have simply ua In order to achieve the correct result we multiply by a mask which contains the value 2 at point D It also contains a 2 at point K so that ou there will acquire the desired value of 2um The corner point F is set to have a value of 1 4 2 2 Temperature salinity and surface elevation The handling of masks by the temperature salinity and surface elevation equations i
157. nux CPU hardware x86_64 Compiler system pgi Compiler command usr local pkg pgi current linux86 64 6 1 bin pgf90 Compiler flags 03 tp k8 64 Mfree SVN Root URL https www myroms org svn omlab branches kate trunk SVN Revision Local Root export staffdata kate roms kate_svn Header Dir export staffdata kate roms kate_svn ROMS Include Header file upwelling h Analytical Dir export staffdata kate roms kate_svn ROMS Functionals Resolution Grid 01 0041x0080x016 Parallel Threads 1 Tiling 001x001 Physical Parameters Grid 01 1440 ntimes Number of timesteps for 3 D equations 300 000 dt Timestep size s for 3 D equations 30 ndtfast Number of timesteps for 2 D equations between each 3D timestep 1 ERstr Starting ensemble perturbation run number 1 ERend Ending ensemble perturbation run number O mnrrec Number of restart records to read from disk T LcycleRST Switch to recycle time records in restart file 288 nRST Number of timesteps between the writing of data into restart fields 1 ninfo Number of timesteps between print of information 92 72 72 72 OOOOE 00 OOOOE 00 OOOOE 00 OOOOE 06 OOOOE 06 OOOOE 05 OOOOE 04 OOOOE 03 OOOOE 02 1 1 OOOOE 00 OOOOE 00 50 000 1025 000 0 000 0 00 OOOOE 00 OOOOE 00 OOOOE 00 OOOOE 00 OOOOE 00 OOOOE 00 14 000 35 000 ldefout nHIS ntsAVG nAVG ntsDIA nDIA tnu2 01
158. o come up with an estimate for time n 1 For reasons of efficiency we choose to use a split explicit time step integrating the depth integrated equations with a shorter time step than the full 3 D equations There is an integer ratio M between the time steps The exact details of how the time stepping is done vary from one version of ROMS to the next with the east coast ROMS described here being older than other branches Still all versions have these steps 1 Take a predictor step for at least the 3 D tracers to time n 5 2 Compute p and p for use in the depth integrated time steps from the density either at time n or time n z 3 Depth integrate the 3 D momentum right hand side terms at time n gt for use in the depth integrated time steps or extrapolate to obtain an estimate of those terms 4 Take all the depth integrated steps Store weighted time means of the u v fields centered at both time n 5 and time n 1 plus at time n 1 The latter requires this time stepping to extend past time n 1 using M steps rather than just M 5 Use the weighted time means from depth integrated fields to complete the corrector step for the 3 D fields to time n 1 Great care is taken to avoid the introduction of a mode splitting instability due to the use of shorter time steps for the depth integrated computations The mode coupling has evolved through the various ROMS versions as shown in Fig 4 from 74 The time steppi
159. o their copy of the array in question stored in the A argument to nf fread2d 73 ifdef DISTRIBUTE CALL mp_scatter2d ng model LBi UBi LBj UBj amp Nghost MyType Amin Amax if defined READ_WATER amp amp defined MASKING amp NWpts SCALARS ng lJwater wtype endif amp Npts wrk A Something similar happens when writing to output files 74 7 Configuring ROMS for a Specific Application This chapter describes the parts of ROMS for which the user is responsible when configuring it for a given application Section describes the process in a generic fashion while and step through the application of ROMS to upwelling downwelling and wind driven Northeast Pacific problems respectively As distributed ROMS is ready to run quite a few examples where the C preprocessor flags determine which is to be executed Some of these examples are described in Haidvogel and Beckmann 26 some are listed here BASIN This is a rectangular flat bottomed basin with double gyre wind forcing When run it produces a western boundary current flowing into a central Gulf Stream which goes unstable and generates eddies The goal is to run adiabatically to study the homogeniza tion of potential vorticity It earned its nickname of Big Bad Basin by taking a long time to run and causing difficulties for the spectral versions of SPEM GRAV_ADJ The gravitational adjustment problem takes place in a long narrow domain which is in
160. ommon newbie mistake is to forget that the commands after a target have to start with a tab G 2 gnu make Over the years the community has moved from the stance of writing portable Makefiles to a stance of just using a powerful portable make The previous section described a portable subset of make features We now delve into some of the more powerful tools available in gnu make 135 G 2 1 Make rules The core of make hasn t changed in decades but concentrating on gmake allows one to make use of its nifty little extras designed by real programmers to help with real projects The first change that matters to my Makefiles is the change from suffix rules to pattern rules I ve always found the SUFFIXES list to be odd especially since f90 is not in the default list Good riddance to all of that For a concrete example the old way to provide a rule for going from file f90 to file o is SUFFIXES o f90 F F90 90 0 lt TAB gt FC c FFLAGS lt while the new way is 4 0 90 lt TAB gt FC c FFLAGS lt In fact going to a uniform make means that we can figure out what symbols are defined and use their standard values in this case FC and FFLAGS are the built in default names for the compiler and its options If you have any questions about this you can always run make with the p or print data base option This prints out all the rules make knows about such as default COMPILE f FC
161. on t work as expected All the ROMS TOMS files are stored in a SVN repository on www myroms org with access controlled by requiring au thentication with the same ROMS Username Password combination assigned to registered users of the ROMS Forum This svn repository is the official version of the code which only the developers are allowed to change Users should download the ROMS code to their local machines using an svn client Don t attempt to use a regular web browser to browse or download files from the svn repository there are much better tools for interacting with the code repository We strongly recommend that users always check out the current trunk version since this has the most recent updates and bug fixes The tags version is kept largely as an historical record of stable releases at the conclusion of major code upgrades Below is a general description of how subversion works Please look at the svn book 9 for more detailed information and the wiki for more on some available GUI clients 1 2 Checking out the code WARNING It is strongly suggested that you checkout the ROMS source code using the same operating system you wish to compile and run ROMS on If you download the code on a Windows machine and wish to run it on a non Windows machine you will need convert the line endings with a utility like dos2unix or recode Even with these utilities you may still have problems compiling ROMS In order download source code from a Subve
162. on means that it will be used to set up the sponge layers for this domain 111 Northeast Pacific sponge areas l Iwrk 10 if defined UV_VIS2 DO i IstrR IendR DO j JstrR MIN Iwrk JendR cff 250 0 5_r8 1 0_r8 C0S pi REAL j r8 REAL Iwrk r8 visc2_r i j max cff visc2_r i j visc2_p i j max cff visc2_p i j END DO END DO DO i IstrR MIN Iwrk IendR DO j MAX JstrR i JendR cff 250 0 5_r8 1 0_r8 COS pi REAL i r8 REAL Iwrk r8 visc2_r i j max cff visc2_r i j visc2_p i j max cff visc2_p i j END DO END DO if defined TS_DIF2 DO itrc 1 NT ng DO j JstrR MIN Iwrk JendR cff 100 1 0_r8 C0S pix REAL j r8 REAL Iwrk r8 DO i IstrR lendR diff2 i j itrc max cff diff2 i j itrc END DO END DO DO i IstrR MIN Iwrk IendR DO j MAX JstrR i JendR cff 100 1 0_r8 COS pi REAL i r8 REAL Iwrk r8 diff2 i j itrc max cff diff2 i j itrc END DO END DO 7 3 3 Model domain A number of horizontal grid points was chosen to resolve the domain at roughly 10 km Values for Lm Mm and N are Lm 224 Number of I direction INTERIOR RHO points Mm 640 Number of J direction INTERIOR RHO points N 60 Number of vertical levels The vertical structure is the same as we ve been using for a good long time now though the minimum depth keeps creeping shallower Vtransform 1 transformation equation stretching function Vstretching 1 THETA_S 5 0d0 THETA_B 0 4d0
163. ontal velocity has a prescribed bottom stress which is a choice between linear quadratic or logarithmic terms The vertical concentration flux may also be prescribed at the bottom although it is usually set to zero 3 3 Horizontal boundary conditions As distributed the model can easily be configured for a periodic channel a doubly periodic domain or a closed basin Code is also included for open boundaries which may or may not work for your particular application Appropriate boundary conditions are provided for u uv T S and The model domain is logically rectangular but it is possible to mask out land areas on the boundary and in the interior Boundary conditions on these masked regions are straightforward with a choice of no slip or free slip walls If biharmonic friction is used a higher order boundary condition must also be provided The model currently has this built into the code where the biharmonic terms are calculated The high order boundary conditions used for u are 2 ves 0 on the eastern and western boundaries and a ves 0 on the northern and southern boundaries The boundary conditions for v and C are similar These boundary conditions were chosen because they preserve the property of no gain or loss of volume integrated momentum or scalar concentration 3 4 Terrain following coordinate system From the point of view of the computational model it is highly convenient to introduce a stretched vertical coordina
164. opic steps and 1 fS OP F dz 42 por Jn OF eo 20 m 0 m M m M Figure 6 Weights for the barotropic time stepping The upper panel shows the primary weights centered at time n 1 while the lower panel shows the secondary weights weights centered at time 1 is the vertically integrated pressure gradient The latter is a function of the bathymetry free surface gradient and the free surface itself as well as the vertical distribution of density The disadvantage of this approach is that after the barotropic time stepping is complete and the new free surface is substituted into the full baroclinic pressure gradient its vertical integral will no longer be equal to the sum of the new surface slope term and the original coupling term based on the old free surface This is one form of mode splitting error which can lead to trouble because the vertically integrated pressure gradient is not in balance with the barotropic mass flux Instead let us define the following E gr af me 4 a 3 p 5 one p oe foe pee 43 Changing the vertical coordinate to o yields 0 0 0 p f pdo gt 2 f f pdo do 44 1 o which implies that p and p are actually independent of as long as the density profile p p o does not change The vertically integrated pressure gradient becomes 1 9 0 pros __0h lg ol DOp _ 0h D D z i 45 el 2 PE Vota OP 45 In the case of uniform density po we obtain p P po but we
165. options I used to give guidelines on reasonable values for THETA_S but I no longer know what s reasonable o coordinates have long had a bad reputation because errors in the pressure gradient terms can lead to spurious currents These errors are must less troublesome than in the past due to code improvements and can also be controlled with some smoothing of the bathymetry This in turn changes the shape of the basin and leads to its own set of problems such as altered sill depths Also the currents will react to the change in shelf slope you are now solving a different problem You may want to explore a matlab tool for minimally smoothing the bathymetry found at http www liga ens fr dutour Bathymetry index html There remain bugs in ROMS If you find any please report them on the forum and or the bug tracking system at myroms org 3 Ocean Model Formulation 3 1 Equations of motion ROMS is a member of a general class of three dimensional free surface terrain following numer ical models that solve the Reynolds averaged Navier Stokes equations using the hydrostatic and Boussinesq assumptions The governing equations in Cartesian coordinates can be written u _ 0 8 f u aTa fu Te eu Fu Dy 1 v 0 5 Ov gt v vt Dy 2 an Vv fu Jy 2 vu V Fy 2 Od _ pg ao 3 with the continuity equation Ou v Ow ee es 4 Ox E Oy a Oz a and scalar transport given by OC 0 53 OC
166. or macro and expects one or more parameters to be passed by the caller We will show some user defined functions in G 3 3 G 2 4 Conditionals We used to have way too many Makefiles a separate one for each of the serial MPI OpenMP versions on each system if supported For instance the name of the IBM compiler changes when using MPI the options change for OpenMP The compiler options also change when using 64 bit addressing or for debugging A better way to organize this is to have the user select 64 bit or not MPI or not etc then use conditionals The complete list of user definable make variables is given in 2 4 1 Gnu make supports two kinds of if test ifdef and ifeq plus the negative versions ifndef ifneq The example from the book is ifdef COMSPEC We are running Windows else We are not on Windows endif An example from the IBM include file is ifdef USE_DEBUG FFLAGS g qfullpath else FFLAGS 03 qstrict endif To test for equality an example is ifeq USE_MPI on Do MPI things endif or ifeq USE_MPI on Do MPI things endif The user has to set values for the USE_MPI USE_DEBUG and USE_LARGE switches in the Makefile or the build script before the compiler dependent piece is included USE_MPI on USE_DEBUG USE_LARGE The Makefile uses the conditional assign in case a build script is used to set them in the environment Be sure to leave the switches
167. other terms being evaluated at time n gt for second order accuracy Setting C to 1 everywhere reduces eq to zy q H uH vH AQ o 0 29 TORRE il mn aa If this equation holds true for the step from time n to time n 1 then our constancy preservation will hold In a hydrostatic model such as ROMS the discrete continuity equation is needed to compute vertical velocity rather than grid box volume Ha the latter is controlled by changes in in the H 0 mn barotropic mode computations Here is the finite volume flux across the moving grid box interface vertically on the w grid The vertical integral of the continuity eq 23 using the vertical boundary conditions on Q is 756 p Of g uD UD e n 0 30 E Be a0 where is the surface elevation D h is the total depth and u v are the depth integrated horizontal velocities This equation and the corresponding 2 D momentum equations are time stepped on a shorter time step than eq and the other 3 D equations Due to the details in the mode coupling it is only possible to maintain constancy preservation to the accuracy of the barotropic time steps 4 5 Depth integrated equations The depth average of a quantity A is given by 7 A gt tado 31 where the overbar indicates a vertically averaged quantity and D C 9 t A n 32 is the total depth of the water column The v
168. otherwise have two new terms The accuracy of these terms depends on an accurate vertical integration of the density as described in Shchepetkin and McWilliams 2005 73 4 7 Time stepping internal velocity modes and tracers The momentum equations and 19 are advanced before the tracer equation by computing all the terms except the vertical viscosity and then using the implicit scheme described in 4 12 to find the new values for u and v The depth averaged component is then removed and replaced by the u and v computed as in A third order Adams Bashforth AB3 time stepping is used requiring multiple right hand side time levels see Appendix A These stored up r h s values can be used to extrapolate to a value at time n gt for use in the barotropic steps as shown in Fig The tracer concentration equation is advanced in a predictor corrector leapfrog trapezoidal step with great care taken to optimize both the conservation and constancy preserving properties of the continuous equations The corrector step can maintain both as long as it uses velocities and column depths which satisfy eq 37 This also requires tracer values centered at time n 5 obtained from the predictor step The vertical diffusion is computed as in The predictor step cannot be both constancy preserving and conservative it was therefore decided to make it constancy preserving Also since it is only being used to compute the advection for the corrector
169. pdate Subversion will remember where you checked out from before and see if a newer revision exists If so it will download and apply all the relevant changes 1 4 Code changes As you use ROMS you may find yourself adding new files and new chunks of code to existing files Unless you are a developer with write access to the repository at www myroms org you have no easy way to save your changes within the svn framework since any one directory can only point to one repository To keep getting updates from the trunk you must keep using the svn server at myroms org At the very least saving a tarball before fetching major updates is a prudent step I 5 Conflicts Sometimes when you make changes to your copy of the ROMS code those changes will conflict with changes made to the repository This means that the changes from the server overlapped with your own and now you have to manually choose between them Whenever a conflict occurs three things typically occur to assist you in resolving that conflict e Subversion halts during the update offering you several choices and remembers that the file is in a state of conflict if you don t clear it right then e If Subversion considers the file to be mergeable it places conflict markers special strings of text which delimit the sides of the conflict usually lt and gt characters into the file to visibly demonstrate the overlapping areas e For every conflicted file
170. peratures and fluxes The hashed layer is the snow 40 9 ROMS directory structure 2 0 a 47 10 ROMS main structure 48 11 Flow chart of the model main program e e 50 12 The whole grid Note that there are Lm by Mm interior computational points The points on the thick outer line and those outside it are provided by the boundary Leah e Els ake oe es Oe a Oe ee 67 13 A tiled grid with some ROMS tile variables 0 0 0 0 0 0 0 00 000048 68 14 A choice of numbering schemes a each tile is numbered the same and b each tile retains the numbering of the parent domain 2 4 69 15 Some ROMS variables for tiles for both a periodic and non periodic case Shown are the variables in the i direction the j direction is similar 70 16 A tiled grid with out of date halo regions shown in grey and the interior points color coded by tile a before an exchange and b after an exchange 71 17 The upwelling downwelling bathymetry 100 18 Surface velocities after one day showing the flow to the left of the wind southern A OA ey eae eae 19 Constant slices of the u v T and w fields at day Ll 102 20 Constant slices of the u v T and w fields at day 5 103 21 Bathymetry of the Northeast Pacific domain NEP5 104 22 Surface elevation after 200 days showing tides This is from a
171. points N Number of vertical levels Nbed Number of sediment bed layers NAT Number of active tracers usually 2 NPT Number of passive tracers NCS Number of cohesive mud sediment tracers NNS Number of non cohesive sand sediment tracers 79 Domain decomposition parameters NtileI Number of direction partitions NtileJ Number of j direction partitions Time stepping parameters NTIMES Number of time steps to evolve the 3 D equations in the current run This is actually the total number including any previous segments of the same run For instance if you already did a three month run and wish to continue for another three months set NTIMES to the number of steps needed for six months DT Time step in seconds for the 3 D equations NDTFAST Number of time steps for the 2 D equations to be executed each dt Model iteration loops parameters ERstr Starting ensemble run number ERend Ending ensemble run number Nouter Maximum number of 4DVAR outer loop iterations Ninner Maximum number of 4DVAR inner loop iterations Nintervals Number of time interval divisions for stochastic optimals computations Generalized Stability Theory GST parameters NEV Number of eigenvalues NCV Number of eigenvectors Input Output parameters ROMS has several possible output files The output files can include a restart file a history file an averages file and a station file for instance The restart file often contai
172. quations with potential temperature salinity and an equation of state Hydrostatic and Boussinesq approximations Optional third order upwind advection scheme Optional Smolarkiewicz advection scheme for tracers potential temperature salin ity etc Optional Lagrangian floats Option for point sources and sinks Orthogonal curvilinear coordinates Arakawa C grid Closed basin periodic prescribed radiation and gradient open boundary conditions Masking of land areas o terrain following coordinate Free surface Tridiagonal solve with implicit treatment of vertical viscosity and diffusivity Hunke and Dukowicz elastic viscous plastic dynamics Mellor Kantha thermodynamics Orthogonal curvilinear coordinates Arakawa C grid Smolarkiewicz advection of tracers Mixing options Horizontal Laplacian and biharmonic diffusion along constant s z or density surfaces Horizontal Laplacian and biharmonic viscosity along constant s or z surfaces Optional Smagorinsky horizontal viscosity and diffusion but not recommended for diffusion Horizontal free slip or no slip boundaries Vertical harmonic viscosity and diffusion with a spatially variable coefficient with options to compute the coefficients with Large et al 41 Mellor Yamada 56 or generic length scale GLS mixing schemes Implementation e Dimensional in meter kilogram second MKS units Fortran 90 Runs under UNIX requires the C preprocessor
173. r equation of state by setting NONLIN_ EOS to undef We want the density to be 26 35 at the bottom and 24 22 at the top with an e folding scale of 50 meters The initial temperature is set to TO 8e7 in ana_ initial The linear equation of state parameters are set in ocean__upwelling in RO 1027 0d0 kg m3 TO 14 0d0 Celsius SO 35 0d0 PSU TCOEF 1 7d 4 1 Celsius SCOEF 0 0d0 1 PSU 90 Since density does not depend on salinity we have a choice of how to handle the second tracer The salinity is set to a uniform value of SO though it could be left out entirely if we undefine SALINITY and set NAT to 1 7 2 5 Boundary conditions The option EW_PERIODIC has been chosen for the eastern and western edges None of the open boundary options have been selected for the Northern and Southern edges they are both then given the default wall conditions and no boundary values are required 7 2 6 Model forcing In this problem we want to resolve the surface Ekman layer and to use a surface wind stress rather than a body force We want the amplitude of the wind to ramp up with time so we modify ana_smflux h accordingly The wind will build to an amplitude of 0 1 Pascals po or 10 4m s We need to edit ana__vmix h to make sure that the vertical viscosity Akv is set to the value we want This must be large at the surface 10 m s to create a thick Ekman layer but has been chosen to decrease with depth We also need to check t
174. r parameterization Rev Geophys 32 363 403 1994 W G Large and S G Yeager The global climatology of an interannually varying air sea flux data set Clim Dyn 33 341 364 2008 DOI 10 1007 s00382 0008 0441 3 J R Ledwell A J Wilson and C S Low Evidence for slow mixing across the pycnocline from an open ocean tracer release experiment Nature 364 701 703 1993 B P Leonard A stable and accurate convective modelling procedure based on quadratic upstream interpolation Comput Method Appl Mech Eng 19 59 98 1979 S J Lin A finite volume integration method for computing pressure gradient force in general vertical coordinates Quart J R Met Soc 123 1749 1762 1997 Jon Loeliger Version Control with Git O Reilly amp Associates Inc Sebastopol CA first edition 2009 A Macks and J Middleton Numerical modelling of wind driven upwelling and downwelling University of New South Wales 1993 J Mailhot and R Benoit A finite element model of the atmospheric boundary layer suitable for use with numerical weather prediction models J Atmos Sci 39 2249 2266 1982 P Marchesiello J C McWilliams and A Shchepetkin Open boundary conditions for long term integration of regional oceanic models Ocean Modelling 3 1 20 2001 L Margolin and P K Smolarkievicz Antidiffusive velocities for multipass donor cell advection SIAM J Sci Comput pages 907 929 1998 John D McCalpin A comparison of s
175. rate where they are valid and damping on unresolved signals 75 Also the preference is for time stepping schemes requiring only one set of the right hand side terms so that different time stepping schemes can be used for different terms in the equations Finally as mentioned in Table 4 3 not all versions of ROMS use the same time stepping algorithm We list some time stepping schemes here which are used or have been used by the ROMS SCRUM family of models plus a few to help explain some of the more esoteric ones F t 230 A 1 Euler The simplest approximation is the Euler time step p t At o t At where you predict the next value based only on the current fields This method is accurate to first order in At however it is unconditionally unstable with respect to advection F t 231 A 2 Leapfrog The leapfrog time step is accurate to O At o t At o t At 2At This time step is more accurate but it is unconditionally unstable with respect to diffusion Also the even and odd time steps tend to diverge in a computational mode This computational mode can be damped by taking correction steps SCRUM s time step on the depth integrated equations was a leapfrog step with a trapezoidal correction LF TR on every step which uses a leapfrog step to obtain an initial guess of t At We will call the right hand side terms calculated from this initial guess F t At F t 232 p t A
176. re Obsolete Long unused versions of the boundary conditions are stored here Programs Not all computer architectures or compilers are the same The types F pro gram checks your compiler for the sizes of the Fortran floating point types Representer This is the representer of the forward model for data assimilation Sealce The sea ice model described in glis here Tangent This is the tangent linear of the forward model for data assimilation Utility Here are utility functions used by the various ROMS routines many dealing with I O Version A file containing the time and date of this svn revision also the svn URL User Some might choose to use this directory rather than the Apps directory It serves the same purpose but is arranged by file type rather than by application Waves The SWAN wave model is here 6 2 Main subroutines 6 2 1 master F The main program is in master F It is simply a shell including one of mct__coupler h esmf_coupler h or ocean h In our case ocean h contains the actual main program which initializes MPI if needed calls ROMS_ initialize calls ROMS_ run with arguments for how many steps to take then ROMS_ finalize and finally wraps up the MPI See Fig 6 2 2 ocean control F This is again a shell which includes one of many other files to do the actual work In this case the worker files all contain ROMS_ initialize ROMS_ run and ROMS_ finalize and live in the ROMS Drivers directory The driver file
177. re not storing water on the surface in melt pools so everything melted at the surface is assumed to flow into the ocean Wro a Inside the ice there are brine pockets in which there is salt water at the in situ freezing tem perature It is assumed that the ice has a uniform overall salinity of S and that the freezing temperature is a linear function of salinity The brine fraction r is given by o Sm The enthalpy of the combined ice brine system is given by E T r r Li CpoT 1 r CpiT 201 Substituting in for r and differentiating gives OE SimLi a i 202 OT T Tor ee Inside the snow we have ks Q BB 203 El The heat conduction in the upper part of the ice layer is 2k Qr A 204 YA These can be set equal to each other to solve for To T3 C T 7 205 e 14 Ck where aeh Substituting into 204 we get 2ki T T3 gt 206 Qs Qr hy 0 206 Note that in the absence of snow Ck becomes zero and we recover the formula for the no snow case in which T T2 At the bottom of the ice we have Qr a 207 i 42 Variable Value Definition b 3 0 factor E evaporation k 0 4 von Karman s constant v 1 8 x 10 m s kinematic viscosity of seawater P precipitation Pr 13 0 molecular Prandtl number Prt 0 85 turbulent Prandtl number So surface salinity Tio stress on the ocean from the ice Tao stress on the ocean from the wind T i
178. rite single precision restart fields OUT_DOUBLE Define to write double precision output fields INLINE_2DIO Define to read write 3D fields level by level 6 8 Important parameters The following is a list of the important parameters in the model These are in mod_ param F and many of them are read from the standard input file as described in 7 1 12 Lm Number of interior grid points in the direction Mm Number of interior grid points in the 7 direction N Number of grid points in the vertical NAT Number of active tracers usually 2 for temperature and salinity NBT Number of biological tracers This will depend on the ecosystem model used NST Number of sediment tracers NPT Number of passive tracers NT Total number of tracer fields NT NAT NBT NST NPT Nfloats Number of Lagrangian floats to track Nstation Number of stations for station output MTC Maximum number of tidal constituents NtileI Number of tiles in the direction NtileJ Number of tiles in the 7 direction 6 9 Domain decomposition ROMS supports serial OpenMP and MPI computations with the user choosing between them at compile time The serial code can also take advantage of multiple small tiles which can be sized to fit in cache All are accomplished through domain decomposition in the horizontal All of the horizontal operations are explicit with a relatively small footprint so the tiling is a logical choice Some goals in the parallel de
179. ror code set it will write out a restart record of the current model fields in case they are useful in diagnosing the trouble wclock_ off ends the built in timers and causes them to print out a report close_ io closes all open files so as to flush the buffers and put NetCDF files into a finished state 6 2 6 main3d This solves the full three dimensional equations described in It has siblings main2d for solving the depth integrated equations and main3d_ offline for reading files from a prior simulation and using them to advect the biological tracers or the Lagrangian floats The full version is shown in Fig Note that many subroutines are optional and only get called if the appropriate C preprocessor switches have been set The subroutines are described as follows set_ data time interpolates between the records read in by get_ data ini_zeta checks for wet dry cells if needed and initializes all the time levels of zeta ini_ fields initializes the 2 D velocities to match the vertical integral of the 3 D velocities making all the time levels match set__massflux computes horizontal mass fluxes Hon and Hav rho_ eos computes the nonlinear equation of state 49 set_ data ana_ vmix first step only rhs3d lmd_ vmix ini_zeta my25_ prestep bvf_mix ini_fields gls_ prestep hmixing set_ massflux lt loop gt omega rho_ eos step2d wvelocity diag step2d set_ zeta radiation_ stress step3d_uv set_ diags ca
180. rsion repository svn client software must be installed on your local machine If you are compiling subversion on your own be sure to build it with SSL support or you will not be able to download the ROMS source code Most Linux distributions come with subversion the command name is svn so shell commands may be used without installing additional software If your username on your local computer is not the same as your ROMS username you will need to pass the username lt username gt option to svn an example is given below The general form of subversion commands is svn action from to optional_qualifiers To check out the files from the ROMS repository trunk svn checkout https www myroms org svn src trunk MyDir where MyDir is the destination directory on your local computer Note the https rather than http If your username on your local computer is not the same as your ROMS username you will need to pass the username option to svn svn checkout username joe_roms https www myroms org svn src trunk MyDir You only check out once after that a hidden directory called svn exists to keep track of the source destination and a bunch of other information Your username and password will also be saved 149 13 Updates Now and again you might feel the urge to get up to speed with the latest changes that have been made to the ROMS repository When that happens simply go to the directory that was MyDir above and type svn u
181. s a user would change to the cpp switches 6 7 the ana_xx h files 6 5 and the ASCII input file 7 1 12 A query on the ROMS forum might be in order if you have something specific in mind If you do choose to add bits know that the makefile fragments called Module mk will attempt to compile anything with a F extension in the directories already populated with such code 7 2 Upwelling Downwelling Example The application which ROMS is configured to run as distributed is a wind driven upwelling and downwelling example described in Macks and Middleton 47 There is a shelf on each wall of a periodic channel and an along channel wind forcing which drives upwelling at one wall and downwelling at the other This problem depends on the Ekman layer so a surface stress is used with vertical viscosity The Ekman depth is estimated to be 9 m if A 0 01m s so the vertical grid spacing must resolve this The maximum depth is 150 m and our choice of the vertical grid parameters leads to a surface Az of 4 0 m 7 2 1 cppdefs h The C preprocessor variable UPWELLING is used for the upwelling configuration of the model The makefile will direct cppdefs h to include the file upwelling h define UV_ADV define UV_COR define UV_LDRAG define UV_VIS2 undef MIX_GEO_UV define MIX_S_UV define TS_U3HADVECTION define TS_C4VADVECTION undef TS_MPDATA define DJ_GRADPS define TS_DIF2 undef TS_DIF4 undef MIX_GEO_TS define MIX_S_TS
182. s approximation ALBEDO Define to use albedo equation for shortwave radiation DIURNAL_SRFLUX Define to impose the local diurnal cycle onto the shortwave radiation 58 General model configuration SOLVE3D Define to solve the 3 D primitive equations MASKING Define if there is land in the domain to be masked out BODYFORCE Define to apply the surface stress as a body force PROFILE Define for time profiling AVERAGES Define to write out time averaged model fields AVERAGES2 Define to write out secondary time averaged model fields AVERAGES_AKV Define to write out time averaged AKv AVERAGES_AKT Define to write out time averaged AKt AVERAGES_AKS Define to write out time averaged AKs AVERAGES_DETIDE Define to write out time averaged detided fields one method FILTERED Define to write out time averaged detided fields using a Lanczos filter AVERAGES_ FLUXES Define to write out time averaged surface fluxes AVERAGES_NEARSHORE Define to write out time averaged nearshore stresses AVERAGES_QUADRATIC Define to write out time averaged quadratic terms AVERAGES_BEDLOAD Define to write out time averaged bed load ICESHELF Define for ice shelf cavities SPHERICAL Define if lat lon coordinates rather than x y SPLINES Define for conservative parabolic spline reconstruction of vertical derivatives STATIONS Define to write out time series information at specific points in the model STATIONS_ CGRID Define if stations are on native C grid F
183. s of the primitive ocean variables and optionally the sediment variables mod_ parallel F This sets up some global variables such as Master which is true for the master thread or process It also initializes the internal ROMS profiling arrays mod_param F This contains the sizes of each grid used plus things like how many tidal con stituents are being used Many of these are read from the input files during initialization not known at compile time mod_scalars F This contains a large number of scalars i e values which don t have spatial dependence Some are fixed constants such as itemp referring to the temperature tracer Others could have a different value on each grid mod_sediment F If either SEDIMENT or BBL_ MODEL is defined this contains parame ters for the respective model mod_sources F If one of UV_PSOURCE TS_ PSOURCE or Q_ PSOURCE is defined this contains the variables used for point sources mod_ stepping F This contains the time stepping variables used to point to the relevant time level mod_storage F If PROPAGATOR is defined this module defines the work space for the Generalized Stability Theory GST Analysis package ARPACK mod_strings F This contains strings such as a title for the run the list of cpp options defined and the names of the sections of code being profiled mod_tides F If SSH_TIDES and or UV_TIDES is defined this will provide the storage for the tidal constituents 6 5 Functionals T
184. s parameter g acceleration of gravity h x y depth of sea floor below mean sea level H x y z vertical grid spacing V Vo molecular viscosity and diffusivity Km Kc vertical eddy viscosity and diffusivity P total pressure P pogz olz y z t dynamic pressure P po Po p x y 2 t total in situ density S x y z t salinity t time T x y z t potential temperature u v W the x y z components of vector velocity Y x y horizontal coordinates z vertical coordinate C x y t the surface elevation Table 1 The variables used in the description of the ocean model 3 2 Vertical boundary conditions The vertical boundary conditions can be prescribed as follows top z x y 1 Km u T mt Km ou T nit 0C_ Q Ko dz a 3 Ot and bottom z h z y Km de TF x y t Km ge me yt OC _ w vu Vh 0 Variable Description Qc surface concentration flux T T surface wind stress Ti T bottom stress Table 2 The variables used in the vertical boundary conditions for the ocean model The surface boundary condition variables are defined in Table 3 2 Since Qr is a strong function of the surface temperature we usually choose to compute Qr using the surface temperature and the atmospheric fields in an atmospheric bulk flux parameterization This bulk flux routine also computes the wind stress from the winds On the variable bottom z h x y the horiz
185. s similar to that in the momentum equations and is in fact simpler Values of T S and inside the land 14 masks such as point H in Fig 3 are set to zero after every time step This point would be used by the horizontal diffusion term for points B J and N This is corrected by setting the first derivative terms at points E I and L to zero to be consistent with a no flux boundary condition Note that the equation of state must be able to handle T S 0 since this is the value inside masked regions 4 2 3 Wetting and drying There is now an option to have wetting and drying in the model in which a cell can switch between being wet or being dry as the tides come in and go out for instance Cells which are masked out as in Fig 3 are never allowed to be wet however e In the case of wetting and drying a critical depth Derit is supplied by the user e The total water depth D h is compared to Deri If the water level is less than this depth no flux is allowed out of that cell Water can always flow in and resubmerge the cell e The wetting and drying only happens during the 2 D computations the 3 D computations see a depth of Derit in the dry areas e The ice component now checks for dry cells when computing the ice rheology 4 3 Time stepping overview While time stepping the model we have a stored history of the model fields at time n 1 an estimate of the fields at the current time n and we need t
186. s stable for longer time steps than many other schemes It is however unstable for advection terms A 5 Forward Backward Feedback RK2 FB One option for solving equation is a predictor corrector with predictor step cate e F u At yrth u Bg cu tr als 1 _ BGC At 238 and corrector step qu 14 L F u tA Flu At 239 urt u de a eG dl 1 ge G At Setting G 0 in the above it becomes a standard second order Runge Kutta scheme which is unstable for a non dissipative system Adding the 8 and e terms adds Forward Backward feedback to this algorithm and allows us to improve both its accuracy and stability The choice of 3 1 3 and e 2 3 leads to a stable third order scheme with Qmax 2 14093 75 121 A 6 LF TR and LF AM3 with FB Feedback Another possibility is a leapfrog predictor AR aa 2F u At n n n n n 240 YE u 42 1 28 G GC The 4 G C HAt and either a two time trapezoidal or a three time Adams Moulton corrector n 1 n 1 n 1 x 1 n n 1 T a Hut 5 t 27 F u F u At 241 wt mara 1 a 1 or 5 27 000 10067 At where the parameters 6 and e introduce FB feedback during both stages while y controls the type of corrector scheme Without FB feedback we have y 0 gt LF TR Cmax V2 max 4 B e 0 y 1 12 gt LF AM3 Qmax 1 5874 y 0 0804 gt max stability Qmax 1 5876 Keeping y 1
187. se it can follow these chained dependencies and do the right thing sfmakedepend assumes that all the files using and defining modules are in the same directory and are all in the list of files to be searched It seems that the industry has not settled on a practical way to deal with a separate modules directory anyway I sometimes include non existent files as a compile time consistency check ifndef PLOTS include must_define_PLOTS bogus include endif This program warns about include files it can t find but not if there is a bogus on the same line See the comments at the top of sfmakedepend for up to date information on the options I may someday get inspired to use a newer version of the getopt routine and rename the options to have names like SGI and Cray 148 I Subversion The ROMS source code is distributed using svn There are svn clients available for nearly every operating system and many resources available including an O Reilly book 9 PH just cover a few basics here taken in part from the ROMS wiki If you wish to maintain your own copy of ROMS in a source code repository you may want to investigate git which has the ability to download from an svn server It too has an O Reilly book 46 plus I enjoyed Swicegood s book 81 I 1 Overview Subversion is a tool for managing software development that keeps track of who modified what and allows the return to a previous version if changes d
188. sign of ROMS were e Minimize code changes e Don t hard code the number of processes 65 e MPI and OpenMP share the same basic structure e Don t break the serial optimizations e Same result as serial code for any number of processes e Portability able to run on any Unix system First some cpp options If we re compiling for MPI the option DMPI gets added to the argument list for cpp Then in globaldefs h we have if defined MPI define DISTRIBUTE Hendif The rest of the code uses DISTRIBUTE to identify distributed memory jobs The OpenMP case is more straightforward with D_OPENMP getting passed to cpp and _ OPENMP being the tag to check within ROMS The whole horizontal ROMS grid is shown in Fig The computations are done over the cells inside the darker line the cells are numbered 1 to Lm in the direction and 1 to Mm in the y direction Those looking ahead to running in parallel would be wise to include factors of two in their choice of Lm and Mm ROMS will run in parallel with any values of Lm and Mm but the computations might not be load balanced 6 9 1 ROMS internal numbers A domain with tiles is shown in Fig The overlap areas are known as ghost points or halo points Each tile is an MPI process an OpenMP thread or a discrete unit of computation in a serial run The tile contains enough information to time step all the interior points so the number of ghost points is dictated by the footprint of the
189. snapshot in a history file the averages files have been detided 0 o o 114 23 Ice concentration averaged over the month of April 1959 115 PENE 116 25 The o surfaces for the North Atlantic with a 9 0 0001 and b 0 b 90 8 and b 0 c 0 8 and b 1 d The actual values used in this domain were 0 5 A ea SENA AAA 124 List of Tables 1 The variables used in the description of the ocean model 9 2 The variables used in the vertical boundary conditions for the ocean model 9 3 The time stepping schemes used in the various ROMS versions a w t is the Courant number and w ck is the frequency for a wave component with wavenumber e e E eaten yee Su Me a ee te ee a 17 ee ee 39 Ae rok tence Meee a es See e ees Al 6 Ocean surface variables lt so resors sad idoti aka ka E a aA 43 T Frazil ice variables dae e a a e a E a G a a A 44 8 Variables used in computing the incoming radiation and latent and sensible heat 128 1 Introduction This user s manual for the Regional Ocean Modeling System ROMS describes the model equations and algorithms as well as additional user configurations necessary for specific applications This manual also describes the sea ice model that we are using Budgell 5 The principle attributes of the model are as follows General e e e J e e Horizontal e e e e Vertical Oo Oo Ice Primitive e
190. sources sinks Q_PSOURCE Define for mass point sources sinks Tracers The default horizontal and vertical advection is 4th order centered The 3rd order upstream split advection TS__ U3ADV_ SPLIT can be used to cor rect for the spurious diapycnal diffusion of the advection operator in terrain following coordinates If this is the case the advection operator is split in advective and diffusive components and several internal flags are activated in globaldefs h Notice that hori zontal and vertical advection of tracer is 4th order centered plus biharmonic diffusion to correct for spurious diapycnal mixing The total time dependent horizontal mixing coef ficient are computed in hmixing F It is also recommended to use the rotated mixing tensor along geopotentials MIX__GEO_ TS for the biharmonic operator TS_U3ADV_SPLIT Define for 3rd order upstream split tracer advection TS _A4HADVECTION Define for 4nd order Akima horizontal advection TS_C2HADVECTION Define for 2nd order centered horizontal advection TS_C4HADVECTION Define for 4rd order centered horizontal advection TS_MPDATA Define for recursive MPDATA 3D advection 50 TS_U38HADVECTION Define for 3nd order upstream horizontal advection TS_A4VADVECTION Define for 4nd order Akima vertical advection TS _C2VADVECTION Define for 2nd order centered vertical advection TS _C4VADVECTION Define for 4rd order centered vertical advection TS_SADVECTION Define for splines vertical advection
191. subs NtileX ng NtileE ng numthreads DO tile subs thread subs thread 1 1 1 CALL set_data ng TILE END DO END DO I OMP END PARALLEL DO What isn t obvious from this is that the argument TILE means different things depending on if we re using OpenMP or MPI ifdef DISTRIBUTE define TILE MyRank else define TILE tile endif Likewise NtileX also depends on whether or not we re using MPI ifdef DISTRIBUTE NtileX 1 Ngrids 1 ttelse NtileX 1 Ngrids NtileI 1 Ngrids endif In other words for MPI TILE becomes the process number and the loop is only executed once In looking at a typical routine that s called from main3d the routine is usually quite short calling a _ tile version of itself in which the actual work happens 69 Non periodic western tile central tile eastern tile IstrIstrU Iend IstnUIstr lend IstrUlstr lend X x gt X XK gt OK x IstrR IendR IstrR lendR IstrR IendR 0 1 2 3 4 5 6 7 8 Lm L Periodic western tile central tile eastern tile IstrnU Istr lend IstiUIstr lend IstrUlstr Tend XK XK XK XK XK XK IstrR lendR IstrR lendR IstrR lendR 1 2 3 4 5 6 7 8 Lm Figure 15 Some ROMS variables for tiles for both a periodic and non periodic case Shown are the variables in the i direction the j direction is similar SUBROUTINE set_data ng tile include tile h CALL set_data_tile ng tile amp amp LBi UBi LBj UBj amp amp IminS
192. t p t _ 1 5 FE F t Ad 233 This leapfrog trapezoidal time step is stable with respect to diffusion and it strongly damps the computational mode However the right hand side terms are computed twice per time step A 3 Third order Adams Bashforth AB3 The time step on SCRUM s full 3 D fields is done with a third order Adams Bashforth step It uses three time levels of the right hand side terms pt At t At where the coefficients a 3 and y are chosen to obtain a third order estimate of p t At We use 0F t BF t At yF t 2At 234 a Taylor series expansion p t ia At HE o t _ yl At n At At AA gl vee 235 120 where Ft lA 1 At I F t At 0 At age Me ae F t 2At g 2Atd 2Ar We find that the coefficients are 23 4 5 12 B gt 17 Q 3 This requires one time level for the physical fields and three time levels of the right hand side information and requires special treatment on startup A 4 Forward Backward In equation above we assume that multiple equations for any number of variables are time stepped synchronously For coupled equations we can actually do better by time stepping asyn chronously Consider these equations o a F u e 236 If we time step them alternately we can always be using the newest information CH FA arth u 4 OLA 237 This scheme is second order accurate and i
193. t currently supported Compilers This contains makefile fragments as described in Data Directories under here contain example forcing grid and initial condition NetCDF files There is also a directory containing the headers of these files in the format produced by ncdump CDL Lib The ARPACK and MCT libraries are needed by the data assimilation codes and by the coupled models respectively makefile This is the standard ROMS makefile as described in Master The ROMS main program is here in various forms for the forward model coupled models and others See ROMS These files are for the ocean model as opposed to other components of the coupled system 46 Bering Apps Bering_10k Atmosphere IL CGOA Adjoint Compilers m Circle c Bin Data L NEP Drivers Libs External makefile m Functionals Master m Include ROMS License _ ROMS txt User Modules Waves Nonlinear Biology Obsolete L Sediment Programs Representer L Sealce Tangent Utility Version Figure 9 ROMS directory structure Adjoint This is the adjoint of the forward model for data assimilation Bin Various shell and Perl scripts for use with the model Note that the sh files are actually csh scripts not sh scripts if it were up to me I d rename them all to csh Drivers The main program includes one of these files depending on ho
194. tails as they apply to you The directories shown here are Apps This directory contains a subdirectory for each of my personal applications The sub directory contains files used by that application the ROMS header file for setting cpp definitions the analytic formulations for fields computed in the model rather than read from files bottom heat flux of zero for instance and ASCII input files read by ROMS on startup to set things such as forcing file names and model time step Some of these applications are Bering This is a 4 km grid of the Bering Sea aligned with the Northeast Pacific grid but at three times the resolution Bering_10k This is a 10 km grid of the Bering Sea a subset of the Northeast Pacific domain with the same extent as the 4 km grid above CGOA This is a 3 km grid of the Gulf of Alaska It is a subset of the Northeast Pacific grid but at four times the resolution Circle This is a circular domain wave propagation problem with an analytic solution used as a test problem 39 NEP _ This is the Northeast Pacific domain covering the waters off the west coast of the US from California to the Bering Sea It is a rectangular domain at about 11 km resolution when viewed in a conformal conic projection with standard latitudes of 40 and 60 N Paul Budgell s applications are also here The application specific files included in the main trunk ROMS are elsewhere Atmosphere This directory is under development no
195. te out current date logical switch to read and plot coastlines and islands bottom and top map latitudes south values are negative left and right map longitudes west values are negative arid dat di arango scrum3 0 plot Palettes gs1 pal d2 kate plot d ice_rst nc scrum_rst_1 nc d2 kate arctic ui coasts coas x x Above FILE 1st 2nd 3rd 4th 5th 6th 7th xx x IREF Sec 1 0 eww Ne efault cnt gridpak arctic_grid_2 nc t dat NAMES line input variables ID file line input color palette file line input contour parameters line input primary NetCDF file line input secondary NetCDF file line input grid NetCDF file line input coastlines file ondary or reference field option Overlay horizontal grid no secondary or reference field to plot plot field overlay from primary file plot field overlay from secondary file primary secondary file field subtraction Day0 DayN field subtraction 118 IPROJ Map Projections option Some of the values for the projection Pole and rotation angle are suggested Check the NCAR manual for details 1 Cylindrical equidistant PLON 0 PLAT 0 ROTA 0 2 Mercator PLON PLAT 0 ROTA 0 3 Lambert conformal conic PLON PLAT ROTA 4 Stereographic azimuthal PLON PLAT 90 or 90 ROTA 0 IVINC JVINC vector grid sampling If either value is negative the velocity vectors at drawn at PSI points Otherwise if
196. te system which essentially flattens out the variable bottom at z h x y Such o coordinate systems have long been used with slight appropriate modification in both meteorology and oceanography e g Phillips and Freeman et al 19 To proceed we make the coordinate transformation and See Appendix B for the form of used here Also see Shchepetkin and McWilliams 2005 for a discussion about the nature of this form of and how it differs from that used in SCRUM In the stretched system the vertical coordinate o spans the range 1 lt lt 0 we are therefore left with level upper a 0 and lower 1 bounding surfaces The chain rules for this transformation are o 0 1 Oz o r Ox mn H Ox Oo o 0 1 Oz o Oy N0y NH X0y 00 481 0 10 Oz X0z 0r H 00 do As a trade off for this geometric simplification the dynamic equations become somewhat more complicated The resulting dynamic equations are after dropping the carats Ou _ 00 gp Oz Oc 1 Km Ou Ot JUTU MUS Ox 2 Gan HOF H 00 where an 8 10 Ov E Od gp Oz Oc 1 Kym Ov po gt Vv Dy eae ee Vu Dy 2 A T E Og H v Fo 9 C _1 0 Ke ac tv FE C tet Pr 10 p p T S P 11 Oo Po OH 0 Hzu O H v HQ MRE hag ae 13 where y u v Q P 0 v V ET T Oy T ra The vertical velocity in coordinates is Q x y o t a o E
197. tegrated time step It is called in a loop over all the short time steps first as a predictor step then as a corrector step step3d_uv completes the time step for the three dimensional velocities omega computes the Q vertical velocity my25__corstep performs the corrector step for turbulent kinetic energy and length scale prog nostic variables tke and gls 57 and 20 my25__corstep performs the corrector step for turbulent kinetic energy and length scale prog nostic variables tke and gls 85 biology computes the changes to the biological tracers due to biological activity using one of several options for the ecosystem model sediment computes changes to the sediment tracers 87 step3d_t completes the tracer time step ice_frazil computes the frazil ice growth if any Not in the trunk code step_ floats time steps the Lagrangian floats 51 6 3 Initialization checkdefs Reports on which C preprocessor variables have been 7 defined and checks their con sistency ana_ grid Computes the grid internally ana_mask Computes the land mask internally get_ grid Reads in the curvilinear coordinate arrays as well as f and h from a grid NetCDF file set_ scoord Sets and initializes relevant variables associated with the vertical transformation to nondimensional o coordinate described in Appendix B set_ weights Sets the barotropic time step average weighting function metrics Computes the metric term combinations whi
198. term The second term of the LMD scheme s surface boundary layer formulation is the non local transport term y which can play a significant role in mixing during surface cooling events This is a redistribution term included in the tracer equation separate from the diffusion term and is written as o Ky az T LMD base their formulation for non local scalar transport on a parameterization for pure free convection from Mailh t and Benoit 48 They extend this parameterization to cover any unstable surface forcing conditions to give 143 wo wTR sO 144 Yr o 144 for temperature and we 145 18 sloh for salinity other scalar quantities with surface fluxes can be treated similarly LMD argue that although there is evidence of non local transport of momentum as well the form the term would take is unclear so they simply specify ym 0 32 The interior scheme The interior scheme of Large McWilliams and Doney estimates the viscos ity coefficient by adding the effects of several generating mechanisms shear mixing double diffusive mixing and internal wave generated mixing vold v vf v 146 Shear generated mixing The shear mixing term is calculated using a gradient Richardson number formulation with viscosity estimated as Vo Rig lt 0 vi 4 voll Rig Rio 0 lt Rig lt Rio 147 0 Rig gt Rio where vo is 5 0 x 1073 Rig 0 7 Double diffusive processes The second component of the
199. the file and let Subversion know that you have resolved the conflict gt svn resolved User Functionals ana_hmixcoef h Resolved conflicted state of User Functionals ana_hmixcoef h 1 5 2 Copying a file onto your working file If you get a conflict and decide that you want to throw out your changes you can merely copy one of the temporary files created by Subversion over the file in your working copy Let s say you want to use the new revision from the repository gt cd User Functionals gt ls ana_hmixcoef hx ana_hmixcoef h ana_hmixcoef h mine ana_hmixcoef h r280 ana_hmixcoef h r291 gt cp ana_hmixcoef h r291 ana_hmixcoef h gt svn resolved ana_hmixcoef h Resolved conflicted state of ana_hmixcoef h 1 5 3 Punting Using svn revert If you get a conflict and upon examination decide that you want to throw out your changes and start your edits again just revert your changes gt cd User Functionals gt svn revert ana_hmixcoef h Reverted ana_hmixcoef h Note that when you revert a conflicted file you don t have to run svn resolved 152 References 1 10 11 12 14 15 16 17 J S Allen P A Newberger and J Federiuk Upwelling circulation on the oregon continental shelf part i Response to idealized forcing J Phys Oceanogr 25 1843 1866 1995 A Arakawa and V R Lamb Methods of computational physics volume 17 pages 174 265 Academic Press
200. the ratio of the surface buoyancy flux to that at the entrainment depth be a constant Thus the entrainment flux at the bottom of the boundary layer under such conditions should be independent of the stratification at that depth Without a turbulent shear term in the denominator of the bulk Richardson number calculation the estimated boundary layer depth is too shallow and the diffusivity at the entrainment depth is too low to obtain the necessary entrainment flux Thus by adding a turbulent shear term proportional to the stratification in the denominator the calculated boundary layer depth will be deeper and will lead to a high enough diffusivity to satisfy the empirical rule of convection 128 turbulent velocity scale To estimate wy where x is m momentum or s any scalar throughout the boundary layer surface layer similarity theory is utilized Following an argument by Troen and Mahrt 84 Large et al estimate the velocity scale as Kx O En where is the surface layer stability parameter defined as z L x is a non dimensional flux profile which varies based on the stability of the boundary layer forcing The stability parameter used in this equation is assumed to vary over the entire depth of the boundary layer in stable and neutral conditions In unstable conditions it is assumed only to vary through the surface layer which is defined as eh sp where e is set at 0 10 Beyond this depth is set equal to its value at ehgpy The f
201. these terms R and Rogon equations and become Du 4 o Duu i o Dut Dft Ot mn oE n On m mn loose gt wy p R Po a Ye oe 35 m ow n 0 mn mn Dv 0 Duv 2 O Duv i Dft t mn 0 n On m mn ew E a gD amp D 1 aod ony P Ram Dry Tp 36 n m m On mn mn When time stepping the model we compute the right hand sides for equations and as well as the right hand sides for equations and 36 The vertical integral of the 3 D right hand sides are obtained and then the 2 D right hand sides are subtracted The resulting fields are the slow forcings R and Roux This was found to be the easiest way to retain the baroclinic contributions of the non linear terms such as uu UU Uslow 19 The model is time stepped from time n to time n 1 by using short time steps on equations 35 36 and 30 Equation is time stepped first so that an estimate of the new D is available for the time rate of change terms in equations and 36 A third order predictor corrector time stepping is used In practice we actually time step all the way to time n dtfast x M while maintaining weighted averages of the values of u v and The averages are used to replace the values at time n 1 in both the baroclinic and barotropic modes and for recomputing the vertical grid spacing H Fig 6 shows one option for how these weights might look The primary weights am are used to compute
202. tnu2 02 visc2 Akt_bak 01 Akt_bak 02 Akv_bak rdrg rdrg2 Zob Vtransform Vstretching theta_s theta_b Tcline rho0 dstart time_ref Tnudg 01 Tnudg 02 Znudg M2nudg M3nudg obcfac TO So to standard output Switch to create a new output NetCDF file s Number of timesteps between the writing fields into history file Starting timestep for the accumulation of output time averaged data Number of timesteps between the writing of time averaged data into averages file Starting timestep for the accumulation of output time averaged diagnostics data Number of timesteps between the writing of time averaged data into diagnostics file Horizontal harmonic mixing coefficient m2 s for tracer 01 temp Horizontal harmonic mixing coefficient m2 s for tracer 02 salt Horizontal harmonic mixing coefficient m2 s for momentum Background vertical mixing coefficient m2 s for tracer 01 temp Background vertical mixing coefficient m2 s for tracer 02 salt Background vertical mixing coefficient m2 s for momentum Linear bottom drag coefficient m s Quadratic bottom drag coefficient Bottom roughness m S coordinate transformation equation S coordinate stretching function S coordinate surface control parameter S coordinate bottom control parameter S coordinate surface bottom layer width m used in vertical coordinate stretching Mean density kg m3 for Boussinesq approximation Time sta
203. ts Some MPI environments require that you submit your job with cd PBS_O_WORKDIR mpirun np 32 oceanM ocean_benchmark3 in while others need aprun np 32 oceanM ocean_benchmark3 in You just have to find out from the locals If all goes according to plan ROMS will create both a collection of NetCDF files and a verbose text file on standard out Chapter 8 describes one way to view the gridded NetCDF files Other tools that have been used include and If things don t go according to plan the text output file is your friend Examine it carefully If it fails on the UPWELLING problem you can compare your output to that in 2 6 Warnings and bugs ROMS is not a large program by some standards but it is still complex enough to require some effort to use effectively Some specific things to be wary of include It is recommended that you use 64 bits of precision rather than 32 bits The code must be run through the C preprocessor before it is compiled This can occasionally be dangerous especially with the newer ANSI C versions of cpp Potential problems are listed in Appendix F The gnu cpp with the traditional flag is known to work well The vertical v coordinate was chosen as being a sensible way to handle variations in the water depth as seen in the coastal oceans Changes to the code have allowed us to expand the well behaved range of depths and the range of values for THETA_S plus there are some new vertical coordinate
204. umed to be in the current directory The C preprocessor has an equivalent include statement include file h We are using the C preprocessor style of include because many of the ROMS include files are not pure Fortran and must be processed by cpp F 2 Macro substitution A macro definition has the form define name replacement text where name would be replaced with replacement text throughout the rest of the file This was used in SCRUM as a reasonably portable way to get 64 bit precision define BIGREAL real 8 It is customary to use uppercase for cpp macros the C preprocessor is case sensitive It is also possible to define macros with arguments as in define av2 al a2 5 Cal a2 although this is riskier than the equivalent statement function av2 al a2 5 al a2 The statement function is preferable because it allows the compiler to do type checking and because you don t have to be as careful about using enough parentheses The third form of macro has no replacement text at all define MASKING In this case MASKING will evaluate to true in the conditional tests described below 129 F 3 Conditional inclusion It is possible to control which parts of the code are seen by the Fortran compiler by the use of cpp s conditional inclusion For example the statements Hifdef MASKING include mask h Hendif MASKING ifdef MASKING c c Apply Land Sea mask slipperiness c do j 1
205. ur profile or login files to globally set them Here is an example for tcsh setenv NETCDF_INCDIR usr local netcdf4 include setenv NETCDF_LIBDIR usr local netcdf4 lib setenv HDF5_LIBDIR usr local hdf5 lib The ksh bash equivalent is export NETCDF_INCDIR usr local netcdf4 include export NETCDF_LIBDIR usr local netcdf4 lib export HDF5_LIBDIR usr local hdf5 lib 2 4 3 Build scripts If you have more than one application or more than one compiler you will get tired of editing the makefile One option is to have a makefile for each configuration then type make f makefile circle_pgi for instance Another option of keeping track of the user defined choices in a build script The advantage is that updates to the build scripts are less frequent than updates to the makefile There are now two of these scripts in the ROMS Bin directory build sh which is surprisingly a csh script and build bash The build scripts use environment variables to provide values for the list above overwriting those found in the ROMS makefile Just as in the multiple makefile option you will need as many copies of the build script as you have applications The scope of these variables is local to the build script allowing you to compile different applications at the same time from the same sources as long as each SCRATCH_ DIR is unique Both scripts have the same options j N Compile in parallel using N cpus omit argument for all available CPUs
206. ve the form J K 06 154 Oo gona we 33 where represents one of u v or C and K is the corresponding vertical viscous or diffusive coefficient This is timestepped using a semi implicit Crank Nicholson scheme with a weighting of 0 5 on the old timestep and 0 5 on the new timestep Specifically the equation of motion for can be written as 3 5 OH K 7 155 a g NO Te ae where Rg represents all of the forcing terms other than the vertical viscosity or diffusion Since we want the diffusion term to be evaluated partly at the current timestep n and partly at the next timestep n 1 we introduce the parameter and rewrite equation 155 as 0 H 0 K K d 1 1 pp mn 5 be Oe E r Gi The discrete form of equation 156 is Artigntl H gn a A K K Zk k zrk kin n k 1 n n Ai mnRg As Pky1 0 Haa ok am A Ky n 1 n 1 Kr 1 n 1 n 1 Forth opty eto asn where k is used as the vertical level index This can be reorganized so that all the terms involving tl are on the left and all the other terms are on the right The equation for oor will contain terms involving the neighbors above and below prti and a which leads to a set of coupled equations with boundary conditions for the top and bottom The general form of these equations is AR BLA Cp Dp 158 where the boundary conditions are written into the coefficients for the end points
207. w you are running the model The forward model is in nl_ocean h External ROMS reads an ASCII file on startup Here are examples for various applica tions also examples of the optional files for extra components such as a sediment model Functionals The file analytical F can include one or more code bits for the analytic specification of for instance the initial conditions Here are examples for the supported model test problems Include Each application has a header file with C preprocessor options for that applica tion For instance the UPWELLING case has the include file upwelling h containing C preprocessor options for its periodic channel domain The full list of available options is in cppdefs h License_ ROMS txt The open source license under which ROMS is copyrighted Modules The ROMS data structures are now in Fortran 90 module files located here Nonlinear The routines used by the nonlinear forward model are here implementing the physics described in Biology The files for the ecosystem parts of the forward model are here 47 mpi_ init initialize_ parallel lt loop gt ROMS _ initialize wclock__on get_ data ROMS_ run inp_ par main3d or main2d ROMS _finalize mod_ arrays mpi_ finalize initial c get_ data if trouble wrt_rst wclock off L close io Figure 10 ROMS main structure Sediment The files for the sediment parts of the forward model are he
208. wdir_ eval omega set_ filter ccsm_ flux my25__corstep set_ avg bulk_ flux gls__corstep set_avg2 ncep_ flux biology output bblm m sediment exit if last set_ vbc step3d_t step done set_ tides ice_ frazil seaice step_ floats Figure 11 Flow chart of the model main program diag computes some global sums prints them and checks them to see if they are sensible If not it stops the model run radiation_stress computes the radiation stresses due to wave current interactions 53 and 54 cawdir_ eval computes a 24 hour mean albedo at the marine surface Not in the trunk code ccsm_ flux computes the surface fluxes from the atmosphere based on a marine boundary layer This version comes from CCSM and is reputed to do better outside of the tropics Not in the trunk code bulk_ flux computes the surface fluxes from the atmosphere based on a marine boundary layer This version comes from COARE version 3 0 LJ and 58 ncep_ flux computes the surface fluxes from the NCEP atmospheric model Not in the trunk code bblm compute the bottom stresses from one of three bottom boundary layer models set_vbc computes the surface and bottom fluxes and stresses that aren t computed elsewhere set vertical boundary conditions set_ tides computes the tidal boundary conditions from the tidal constituents seaice runs the sea ice model described in It changes the surface boundary conditions for the ocean and therefore gets called before the
209. with traditional Fortran include statements I later wrote a variant of it for use with the C preprocessor called sfmakedepend The latest version of sfmakedepend does the job of both programs and also searches for the dependencies introduced by Fortran 90 modules It is used by the Makefiles described in It recursively searches for Fortran style includes for instance if file f has the statement include commons h the line file o commons h will be added to the bottom of the Makefile This tells make that file o depends on commons h as well as file f and to recompile file f whenever commons h is modified It likewise searches source files for C style includes such as include commons h and adds the corresponding dependencies to the Makefile It has several options including s required for Fortran compilers which will not invoke the C preprocessor for you In this case the above dependency line would become file o commons h file f commons h letting make know that the C preprocessor must be rerun on file F whenever commons h is updated When using the C preprocessor you can ask it to search directories other than the current directory Likewise sfmakedepend can be instructed to search other directories with I dir options Note that it is legal to have more than one I dir option as in sfmakedepend I usr local include I home me include F Fortran 90 introduces some interesting dependencies Two compilers I have ac
210. y ANA_ TOBC or they can be read from a boundary NetCDF file There is logic in globaldefs h by which ROMS decides whether or not it needs to read a boundary file 7 1 11 Model forcing a Winds and thermal fluxes There are two different ways to apply a wind forcing as a surface momentum flux in the vertical viscosity term or as a body force over the upper water column Usually we set the vertical coordinate parameters to retain some resolution near the surface and apply the fluxes as boundary conditions to the vertical viscosity diffusivity In either case the surface and bottom fluxes are either defined analytically read from a forcing file or computed inside ROMS using a bulk flux formula from the appropriate atmospheric fields air temperature and winds for instance You must either edit the appropriate ana__xxx h or create a NetCDF forcing file in the format expected by ROMS Note that it is quite common to put the surface variables in the forcing file while having an analytic bottom heat flux ROMS now has the capability of reading in a list of forcing fields it can be convenient to have one file per field rather than stuffing tides winds river inputs all into one file In the past our vertical resolution was relatively coarse and the vertical viscosity would have to have been unreasonably large for us to resolve the surface Ekman layer If that is your situation define BODYFORCE in cppdefs h and provide a value for levsfre in o
211. y solving Ss 1 342B 1841 42 G 42 1 6A B7 121 Sm 1 941 42Gp Ss Gr 184i 94142 G Ai 1 3C1 641 B7 122 2 2 Gr min Zz 0 028 123 Sq 0 415m 124 The constants are set to Aj A2 Bi Bo C1 E1 Ez 0 92 0 74 16 6 10 1 0 08 1 8 1 33 The quantities q and q l are both constrained to be no smaller than 1078 while 1 is set to be no larger than 0 53q N 4 11 2 The Large McWilliams and Doney parameterization The vertical mixing parameterization introduced by Large McWilliams and Doney is a versatile first order scheme which has been shown to perform well in open ocean settings Its design facilitates experimentation with additional or modified representations of specific turbulent processes Surface boundary layer The Large McWilliams and Doney scheme LMD matches separate parameterizations for vertical mixing of the surface boundary layer and the ocean interior A formulation based on boundary layer similarity theory is applied in the water column above a calculated boundary layer depth hsa This parameterization is then matched at the interior with schemes to account for local shear internal wave and double diffusive mixing effects Viscosity and diffusivities at model levels above a calculated surface boundary layer depth hsp are expressed as the product of the length scale hsp a turbulent velocity scale wy and a non dimensional shape function Vz hspiW2 7 Gz o
212. y starting the line with a G 1 2 Implicit Rules Make has some rules already built in For fortran you might be able to get away with OBJS main o init o plot o model OBJS FC LDFLAGS o model 0BJS 134 as your whole Makefile Make will automatically invoke its default Fortran compiler possibly f77 or g77 with whatever default compile options it happens to have FFLAGS One built in rule often looks like CC CFLAGS c lt which says to compile c files to o files using the compiler CC and options CFLAGS We can write our own suffix rules in this same style The only thing to watch for is that make by default has a limited set of file extentions that it knows about Let s write our Makefile using a suffix rule Cray version F90 90 F9OFLAGS 03 LDFLAGS F90 FOOFLAGS c lt model main o init o plot o F90 LDFLAGS o model main o init o plot o clean rm o core G 1 3 Dependencies There may be additional dependencies beyond the source gt object ones In our little example all our source files include a file called commons h If commons h gets modified to add a new variable everything must be recompiled Make won t know that unless you tell it include dependencies main o commons h init o commons h plot o commons h Fortran 90 introduces module dependencies as well See SH for how we automatically generate this dependency information The most c
213. yes ISOVAL iso surface value to process see below VLWD vector line width 1 0 for default VLSCL vector length scale 1 0 for default IVINC vector grid sampling in the X direction 1 for default JVINC vector grid sampling in the Y direction 1 for default IREF secondary or reference field option see below IDOVER overlay field identification for IREF 1 2 only LEVOVER level of the overlay field set to 0 if same as current FLDLEV RMIN overlay field minimum value to consider 0 0 for default RMAX overlay field maximum value to consider 0 0 for default 117 10 0 LGRID 4 IPROJ 60 0 PLON 90 0 PLAT 0 0 ROTA 0 LMSK 1 NPAGE T READGRD F PLTLOGO T WRTHDR T WRTBLAB T WRTRANG T WRTFNAM T WRTDATE T CST 50 0 50 0 110 0 80 0 d2 kate plot v Desired longitude latitude grid spacing degrees map projection see below projection Pole longitude west values are negative projection Pole latitude south values are negative projection rotation angle clockwise degrees flag to color mask land 0 no 1 yes number of plots per page currently 1 2 or 4 logical switch to read in positions from grid NetCDF file logical switch draw Logo logical switch to write out the plot header titles logical switch to write out the plot bottom title logical switch to write out data range values and CI logical switch to write out input primary filename logical switch to wri

Technical Manual for a Coupled Sea-Ice/Ocean Circulation

Contents

Download Pdf Manuals

Related Search

Related Contents