Application setup. Please help!

Discussion of how to use ROMS on different regional and basin scale applications.

Moderators: arango, robertson

Post Reply
Message
Author
Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

Application setup. Please help!

#1 Unread post by Diego »

Preamble:
* I am new to ROMS - new to linux.
* I got ROMS up & running and ran several test cases
* I read most articles in WikiROMS, most tutorials / Webinars, most posts in this forum and the 2000 SCRUM manual.
* I started with (and successfully ran) a simple 5 by 5 grid application with periodic boundaries.
* I started scaling up to a coastal grid when... hell broke loose! ...I keep causing ROMS to blow-up within the first time-step. I need help!
* Payment: In return for your help, I plan to pay with a nice series of WikiROMS tutorials... so that the next generation of ROMS users don't have to ask too many questions in this forum.

-------------------

Step one
As I said above, I started with a 5 by 5 grid with East-West and North-South periodic boundaries that I got from Katia Fennel. She uses this grid as a simple 1-D model to do quick tests. It seemed a good idea to start simple and work my way up. As I said, this first try ran fine, here are cpp definitions that I used for compilation:

Code: Select all

#define UV_ADV                   /* use to turn ON or OFF advection terms  */
#define UV_COR                   /* use to turn ON or OFF Coriolis term    */
#define UV_QDRAG                 /* use to turn ON or OFF quadratic bottom friction */
#define DJ_GRADPS                /* use if splines density Jacobian (Shchepetkin, 2000) */
#define UV_VIS2                  /* use to turn ON or OFF harmonic horizontal mixing */
#define MIX_S_UV                 /* momentum mixing on s-surfaces */
#define TS_DIF2                  /* use to turn ON or OFF harmonic horizontal mixing */
#define MIX_GEO_TS               /* tracer mixing on constant z surfaces */
#define TS_U3HADVECTION          /* use if 3rd-order upstream horiz. advection */
#define TS_C4VADVECTION          /* use if 4th-order centered vertical advection */
#define TS_MPDATA                /* use if recursive MPDATA 3D advection  */
#define SOLAR_SOURCE             /* use if solar radiation source term */
#define NONLIN_EOS               /* use if using nonlinear equation of state */
#define SALINITY                 /* use if having salinity */
#define SPLINES                  /* use to activate parabolic splines reconstruction */
#define AVERAGES                 /* use if writing out time-averaged data */
#define AVERAGES_FLUXES          /* use if writing out time-averaged fluxes */
#define AVERAGES_AKV             /* use if writing out time-averaged AKv */
#define AVERAGES_AKT             /* use if writing out time-averaged AKt */
#define SOLVE3D                  /* use if solving 3D primitive equations */
#define EW_PERIODIC              /* use if East-West periodic boundaries */
#define NS_PERIODIC              /* use if North-South periodic boundaries */
#define MY25_MIXING              /* use if Mellor/Yamada Level-2.5 closure */

#ifdef MY25_MIXING                 
# define N2S2_HORAVG             /* use if Large et al. (1994) interior closure */
# define KANTHA_CLAYSON          /* use if Kantha and Clayson stability function */
#endif                

#undef ECOSIM
#if defined ECOSIM || defined BIO_FASHAM
# undef ANA_BIOLOGY              /* use if analytical biology initial conditions */
# define ANA_SPFLUX              /* use if analytical surface passive tracers fluxes */
# define ANA_BPFLUX              /* use if analytical bottom passive tracers fluxes */
# define ANA_CLOUD               /* use if analytical cloud fraction */
#endif

#define BULK_FLUXES              /* use if bulk fluxes computation */

#ifdef BULK_FLUXES
# define EMINUSP                 /* use if computing E-P */
# define LONGWAVE                /* use if computing net longwave radiation */
# undef ANA_RAIN                 /* use if analytical rain fall rate */
#else
# define ANA_SMFLUX              /* use if analytical surface momentum stress */
# define ANA_STFLUX              /* use if analytical surface temperature flux */
#endif

#if defined BULK_FLUXES || defined ECOSIM
# undef ANA_CLOUD
# undef PAPA_CLM                 /* ? */
#endif
#define ANA_SSFLUX                /* use if analytical surface salinity flux */
#define ANA_BSFLUX                /* use if analytical bottom salinity flux */
#define ANA_BTFLUX                /* use if analytical bottom temperature flux */
Initialization and forcing: I initialized the model with the same value of: salt, temp, u, ubar, v, vbar, ocean_time and zeta... for each of the cells of the 5 by 5 grid. I forced the model with time-series of Pair, Qair, Tair, Uwind, Vwind, cloud, rain, swrad and data_time. I got the initialization and forcing NetCDF files from Katja.

Step 2: Varying the grid size
I modified the Matlab code that Katja used to create the 5by5 grid to do a 25by40 grid... Then I wrote Matlab code to do the corresponding initialization NetCDF file (i.e. that match the new grid size). The larger grid also ran fine. Then, I played with different initialization values and with different grid sizes, but with the same forcing file. All tests ran successfully.

Step 3: Realistic Grid and masking
I tried SeaGRID but my Matlab just keep crashing and not linking it. So I tried GRIDGEN but I got stuck with Python... so I decided to write my own Matlab code (following Katja's 5by5-grid as an example): I added the capability to rotate the grid, but I only allow for rectangular domains and squared grid cells. I also included bathymetry. Below is an image of the grid. The geographical location is a fjord in Nova Scotia called Ship Harbour. This first grid is probably too coarse, but it is just for testing.

Image


The next figure is just a zoom of the figure above. I included, for reference, the rho, U, V and PSI grids. All graphs were generated with John Wilkin's version of pcolor (I was going to put a link to the download, but I can't find where I got it).
Image



Then I created the masks. The rho_mask was created with the bathymetry (depths > 0 are land). Then I created the U, V and PSI masks from the RHO mask as follows:

Code: Select all

%This is a Matlab code snipet...

%u mask
for j = 1:length(mask_rho(:,1));
    for i = 1:length(mask_rho(1,:))-1;
        if mask_rho(j,i)==0 | mask_rho(j,i+1)==0;
            mask_u(j,i) = 0;
        else
            mask_u(j,i) = 1;
        end
    end
end

%v mask
for i = 1:length(mask_rho(1,:));
    for j = 1:length(mask_rho(:,1))-1;
        if mask_rho(j,i)==0 | mask_rho(j+1,i)==0;
            mask_v(j,i) = 0;
        else
            mask_v(j,i) = 1;
        end
    end
end

%psi mask
for j = 1:length(mask_rho(:,1))-1;
    for i = 1:length(mask_rho(1,:))-1;
        if mask_u(j,i)==0 | mask_u(j+1,i)==0 | mask_v(j,i)==0 | mask_v(j,i+1)==0;
            mask_psi(j,i) = 0;
        else
            mask_psi(j,i) = 1;
        end
    end
end





Below is the RHO mask for the same (zoomed) domain of the figure above (blue=land=0 and red=ocean=1):
Image




I recompiled ROMS with all the cpp's from above but including also masking and turning off the periodic boundaries.

Code: Select all

#undef EW_PERIODIC              /* use if East-West periodic boundaries */
#undef NS_PERIODIC              /* use if North-South periodic boundaries */
#define MASKING                       /* use if analytical masking is enabled */

ROMS run
I used the same forcing NetCDF file used in "step one" above... but I created a new initialization NetCDF file, where parameters where defined at the rho points. ROMS blew up in the first time step

I went back to the working 25 by 40 grid... and I did some basic experimentation:

* Turn off periodic boundaries.
RESULT: ROMS blows up

* Turn periodic boundaries back on and add a 5by5 island in the middle on the grid.
RESULT: ROMS blows up

* Same as above but with: #define MASKING
RESULT: ROMS blows up

* Modify as follows:

Code: Select all

#define WESTERN_WALL 	/*use if Western edge, closed wall condition*/
#define NORTHERN_WALL 	/*use if Northern edge, closed wall condition*/
#define SOUTHERN_WALL 	/*use if Southern edge, closed wall condition*/
#define EASTERN_WALL 	/*use if Eastern edge, closed wall condition*/
RESULT: ROMS blows up



At this point is pretty obvious that I can keep going with this trial-an-error system for a few decades (or at least until November 2008 when the ROMS manual comes out). It is time to ask for help!

As I mentioned above, I plan to make tutorials for WikiROMS where I will include the "cleaned up" outcome of this process of setting up a realistic application. This will be very useful for next generations of users.

Thanks in Advance!


Diego

jcwarner
Posts: 1172
Joined: Wed Dec 31, 2003 6:16 pm
Location: USGS, USA

#2 Unread post by jcwarner »

just a quick thought: are any of your depth (h) values == 0 ?

Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

#3 Unread post by Diego »

I just double checked... and No. None of my depths are 0.

User avatar
arango
Site Admin
Posts: 1347
Joined: Wed Feb 26, 2003 4:41 pm
Location: DMCS, Rutgers University
Contact:

#4 Unread post by arango »

Usually, if ROMS blow-up in the first time step is either because a set-up problem in the grid, initial conditions, and forcing or a time-step CFL violation. This implies that the model is blowing-up in the barotropic engine (step2d). You probable are no resolving the gravity wave physics with your current time step. Then, you need to reduce your baroclinic time-step (DT) or increase the number of barotropic steps within a baroclinic step (NDTFAST).

The CFL condition for gravity wave (fast) physics is governed by dt < SQRT( (dx^2+dy^2)/(g*h(x,y) ).

Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

#5 Unread post by Diego »

Using the formula above... I calculated dt = 7.4 seconds

I adjusted (to be in the safe side) as follows:

Code: Select all

DT == 1
NDTFAST == 100
Unfortunately, ROMS is still blowing up. That means that likely my problem is in the grid, initial conditions or forcing. I will see if I can find what's the problem...

more soon.

Diego

User avatar
kate
Posts: 4088
Joined: Wed Jul 02, 2003 5:29 pm
Location: CFOS/UAF, USA

#6 Unread post by kate »

Since it is blowing up in the first step, this is when I would go to a debugger. Perhaps first look for compiler options to dump core when you get Inf or NaN instead of silently running with garbage until the next check in diag. There were some recent articles about it here: http://www.arsc.edu/support/news/HPCnews.shtml
Check out issues 376 and 377.

RubenDiez-Lazaro

#7 Unread post by RubenDiez-Lazaro »

Diego:

I have a very similar problem, and i follow similar steps from scratch to today...

My application run right at open ocean, but when i put land, the model blow up.....
I think that the error is on the grid stuff; I ran some applications with analytical versions for all data except for the grid, and they blew up....

I try SeaGRID too, and, like in your case, just crashing , so I wrote a C code for grid generation, based on "gridgen" from Pavel Sakov.

I suspect that the cause of the problem are the so named "rvalues": What are the "Maximum grid stiffness ratios" (rx0 and rx1) at the roms output in your case?

I try to smooth the rvales using the "smoothgrid.m" script from the roms-tools (roms-agrif) project, but it seems not work for my case...

Now i am on and "trial-an-error" phase....

What about collaborate for find where is our mistake???

You can get one of my grids at http://web.usc.es/~rdlazaro/roms_stuff/ ... d_probe.nc

Do you see any strange thing on this grid??

Regards

Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

#8 Unread post by Diego »

Hi Ruben,

Joining efforts makes sense... Although I may be too new to ROMS to be helpful.

I looked at your grid and I plotted mask_rho (figure below). I can see some "one-pixel-lakes" (arrow A) that I suspect could cause ROMS to blow up when trying to numerically advect things in or out of those isolated pixels. Arrow B shows some "Rias" that are not connected to the ocean body... which also may cause troubles. Finally, you have a few "one-pixel-islands" (arrow C), I don't know if they are problem or not. I would try to manually mask (or de-mask) all this irregularities and try running it again.

Image





I spent some time checking my grid and initial-conditions files, and I couldn't find problems (although, I am not sure what to look for). I will try smoothing the bathymetry of my grid... I'll let you know soon if it worked or not.


I thought about trying a test case that uses a grid file and a initial-conditions file... After a successful run of such test case... I was planning to replicate the grid & init files with my own code to see it still runs. The problem is that, out of the test cases that explored within ROMS/Include (e.g. nj_bight.h), I couldn't find the associated grid, initial and forcing files.

Maybe a more senior ROMS-user could direct us where to find a (simple) ROMS application, with its associated grid/initial/forcing files.


Kate, Thanks for the link to the ARSC Newsletters! I have been reading the suggested articles... they are bit thick for my background... so it is taking me a while to get through them... but now worries, I'll get there :D

Thanks you all for your replies!

More soon...

Diego

User avatar
kate
Posts: 4088
Joined: Wed Jul 02, 2003 5:29 pm
Location: CFOS/UAF, USA

#9 Unread post by kate »

Ruben, you should definitely try running editmask in Matlab. I don't think the lakes should cause trouble, but one grid waterways are hard to visualize so that if your troubles are there, it might be hard to diagnose. That's why there's a recommendation out there to simply mask the one-grid bays and the disconnected lakes. Note that the island pointed to by the upper C arrow has a one-grid bay in it. The tiny islands do not cause trouble in my experience.

Diego - poke Dave Robertson for wiki rights!

RubenDiez-Lazaro

#10 Unread post by RubenDiez-Lazaro »

Thanks you very much for the indications.

The grid i have sent is the "virgin" one, without any post manipulation....

Some of the things i have tried was remove one-cell bays, islands and internal lakes.
This is same grid after manual modification using "editmask":
http://web.usc.es/~rdlazaro/roms_stuff/ ... itmask1.nc

The application also blew up using this "editmaskared" grid...

Another thing i have tried is use the smoothing routines from the "roms-tools" (the script named "smoothgrid.m") instead the mine one, but it also blow up...

This is the "editmaskared" grid smoothed with the "smoothgrid.m" script:
http://web.usc.es/~rdlazaro/roms_stuff/ ... othgrid.nc
that also blew up.

And, finally, this is my application configuration: http://web.usc.es/~rdlazaro/roms_stuff/mog.h

Do you see any other mistake or strange thing in this files?

Another time, thank very much by the help

Regards

User avatar
wilkin
Posts: 875
Joined: Mon Apr 28, 2003 5:44 pm
Location: Rutgers University
Contact:

#11 Unread post by wilkin »

A clarification: are any of your depths h>=0 (not just ==0)?

Your nested "for" loops of mask generating code (tsk tsk - no need for the loops in Matlab) is a bit hard to decipher, and it is essential that the u,v masks be consistent with rho. Editmask will enforce this for you by calling uvp_masks (no for loops) (be careful if you want to run uvp_mask outside of editmask - Hernan's conventions for reading netcdf variables is transposed from the way snctools reads the data).

Check your log file: At initialization ROMS reports min/max values of all the initial conditions. Are any values unreasonable? Even if they fall on land they will cause trouble because ROMS computes on all points and applies the mask afterward.
John Wilkin: DMCS Rutgers University
71 Dudley Rd, New Brunswick, NJ 08901-8521, USA. ph: 609-630-0559 jwilkin@rutgers.edu

Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

#12 Unread post by Diego »

John,

Thanks for the clarification... YES! all my "above-sea-level" cells are >0. Mountains go up 50 m while my deepest point only goes to -20 m. As soon as I click "submit" I will start my experimentation with Editmask.


Also, when I got your reply I was analyzing grid files I found at DATA/ROMS/Grid/

The variable h in the file inlet_test_grid.nc has positive depths... while mine are negative (except for above-water land). I am not sure if that can be a problem or not.

Thanks!!!!

more soon

User avatar
wilkin
Posts: 875
Joined: Mon Apr 28, 2003 5:44 pm
Location: Rutgers University
Contact:

#13 Unread post by wilkin »

Sorry, I should have asked "any depths h<=0' which means land above z=0.

Unless you have WET_DRY defined (which you do not want to be experimenting with as a novice) you should probably set min(h) to something like a little more than the maximum tidal variation you expect.
John Wilkin: DMCS Rutgers University
71 Dudley Rd, New Brunswick, NJ 08901-8521, USA. ph: 609-630-0559 jwilkin@rutgers.edu

jcwarner
Posts: 1172
Joined: Wed Dec 31, 2003 6:16 pm
Location: USGS, USA

#14 Unread post by jcwarner »


Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

#15 Unread post by Diego »

Thanks jcwarner!

(1) I change my depths... now they are positive (1 to 23 m).
(2) I use editmask to erase 1-pixel bays... also this should have setup mask_u, mask_v and mask_psi correctly.

Still blows up :(

CLICK Here for the log file, which ends with... " MAIN: Abnormal termination: BLOWUP."

CLICK Here for the ocean_application.in file that I used.

CLICK Here for the header.h file that I used to compile ROMS.

CLICK Here to download the grid.nc file that I used.

CLICK Here to download the Initial_conditions.nc file that I used.

CLICK Here to download the Forcing.nc file that I used.

PS. Thanks to Dave for giving me posting-rights for WikiROMS... a nice synthesis of all this process will go there as soon I get this up and running.

bberx
Posts: 9
Joined: Thu Oct 11, 2007 1:33 pm
Location: Marine Scotland - Science

#16 Unread post by bberx »

Hi Diego,

I've been very keenly following this topic because I have similar intentions. Anyway, my computer is now all set, so today I've also been attempting my own set-up. Did you activate wetting-drying? Am doing my run while writing this, but will let you know if mine doesn't blow up. Luckily it is almost weekend time here :D

Bee

User avatar
wilkin
Posts: 875
Joined: Mon Apr 28, 2003 5:44 pm
Location: Rutgers University
Contact:

#17 Unread post by wilkin »

Your initial u,v (=0.05) is inconsistent with your initial ubar,vbar (=0 everywhere). Probably not a good idea.
John Wilkin: DMCS Rutgers University
71 Dudley Rd, New Brunswick, NJ 08901-8521, USA. ph: 609-630-0559 jwilkin@rutgers.edu

Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

#18 Unread post by Diego »

John,

I made a new initial_conditions.nc file where u, v, ubar and vbar are 0 everywhere. ROMS still BLOWS UP.

(huff... I'm running out of coffee :()

User avatar
wilkin
Posts: 875
Joined: Mon Apr 28, 2003 5:44 pm
Location: Rutgers University
Contact:

#19 Unread post by wilkin »

Here it is:
Step 3: Realistic Grid and masking
I tried SeaGRID but my Matlab just keep crashing and not linking it. So I tried GRIDGEN but I got stuck with Python... so I decided to write my own Matlab code ...
The stdout log shows

Code: Select all

 Minimum X-grid spacing, DXmin =  3.83055437E-05 km
 Maximum X-grid spacing, DXmax =  3.83055437E-05 km
 Minimum Y-grid spacing, DYmin =  5.39956803E-05 km
 Maximum Y-grid spacing, DYmax =  5.39956803E-05 km
but your X grid spacing is not really 3.83 centimeters!

You have incorrectly calculated the pm and pn grid metric factors which are the reciprocal of the cell dimensions centered on each rho point.

Working with the lon/lats in your grid file I calculate a grid size for rho cell 0,0 more like 100 m, which implies pm=0.01. Your pm(0,0)=26.
John Wilkin: DMCS Rutgers University
71 Dudley Rd, New Brunswick, NJ 08901-8521, USA. ph: 609-630-0559 jwilkin@rutgers.edu

Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

#20 Unread post by Diego »

John,

Thanks for spending the time reviewing my case! I really appretiate it!

In my code, I made the grid with cells that are exactly 100 by 100 m. Then I rotate the grid (while in meters) and finally I calculate the lat lon of every cell point with respect of the bottom-left corner of the grid (using formulas to convert meter to degrees Lat and meter to degrees Lon as a function of lat). Somehow my code is screwing things up... I think I will give it one last try to fix it (before going back to try to get GRIDGEN working)... just because I can see the use of a (very non-fancy) matlab-based grid-generation code that runs in most Matlabs (i.e. few dependencies) for the ROMS community... a quick and dirty of sorts.

Thanks again John!

more soon...

nganju
Posts: 82
Joined: Mon Aug 16, 2004 8:47 pm
Location: U.S. Geological Survey, Woods Hole
Contact:

#21 Unread post by nganju »

From the output file: your barotropic courant number (>5) seems high, i think <1 is more appropriate...perhaps your grid resolution is too fine (dx = 3.8 e -5 km = 0.038 m ???????). something is wrong with your grid spacing, i assume?

nganju
Posts: 82
Joined: Mon Aug 16, 2004 8:47 pm
Location: U.S. Geological Survey, Woods Hole
Contact:

#22 Unread post by nganju »

sorry, i was late to the dance

Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

#23 Unread post by Diego »

ROMS is Running!!

After a long battle with my matlab code that generates the grid, ROMS is running without blowing up (many many Thanks to John Wilkin!! for his hints, tricks and tips). The joy is big! Unfortunately, I cannot open the output .nc files using matlab's ncload... hence I haven't been able to see what came out of the simulation.

Matlab error message:

Code: Select all

??? varget1:  cannot have empty set in input position 4.


Error in ==> mexcdf53 at 9
	[varargout{:}] = feval('mexnc', varargin{:});

Error in ==> ncmex at 139
	[varargout{:}] = feval(fcn, varargin{:});

Error in ==> ncvar.subsref at 247
			[result, status] = ncmex('varget1', ncid(self), varid(self), ...

Error in ==> netcdf.subsref at 52
         result = subsref(result, s);

Error in ==> ncload at 30
   assignin('caller', varargin{i}, f{varargin{i}}(:))
I can open the output .nc files using "Panoply", so the variables are there. But ncload is been fuzzy.

Anyway, I haven't done any digging in this matter. I just thought share the joy of the running ROMS before continuing :D

more soon...

User avatar
wilkin
Posts: 875
Joined: Mon Apr 28, 2003 5:44 pm
Location: Rutgers University
Contact:

#24 Unread post by wilkin »

I recommend snctools in matlab to read/write netcdf files.
http://sourceforge.net/project/showfile ... _id=130914

But first make sure you have the correct mexnc for your architecture
http://sourceforge.net/project/showfile ... _id=119464
John Wilkin: DMCS Rutgers University
71 Dudley Rd, New Brunswick, NJ 08901-8521, USA. ph: 609-630-0559 jwilkin@rutgers.edu

RubenDiez-Lazaro

#25 Unread post by RubenDiez-Lazaro »

Diego:

Can you put a link to your (working) grid, please??

I would like try it and compare it with mine....

Thanks in advance

Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

#26 Unread post by Diego »

Ruben,

Yes I will... even better, I am planning to also post my grid-generating code (with fixes from many folks in this forum... thanks a lot!). It is nothing fancy... but it may be useful for some to have a "quick-&-dirty" alternative to Seagrid and Gridgen for quick tests.

The problem is that my NetCDF utilities were working a bit funny... I was not able to read some .nc files... However, in the process of trying to fix it, I screw it up more... now I can't do anything with .nc files (hence can't do final tests on my code). I am trying to fix it as we speak... I have the feeling the problem is related to the fact that I have Matlab R14 service pack 1... and available binaries are only R14...

More in ~1 hr.

Diego

Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

#27 Unread post by Diego »

Ruben, sorry for taking so long... I was able to get my Matlab MEXCDF back to "not-too-bad"... but I still have troubles opening some files (but others are fine).

Anyway... here is the test grid file that works (i.e. doesn't make ROMS blow-up). Also you may want to take a look to roms_make_simple_grid.m, which is my matlab code to "quick-and-dirty" generate the grid and initial-conditions .nc files to input into ROMS. You can use it to make you own grid (cross fingers that works for you) or perhaps even more useful... you can look at the code to get ideas on how to fix your own C code. I tried to keep very detailed comments throughout my code so that it is easy to follow. I also tried to code it as self-sustained as possible... the only extra thing you should need, is the proper mexnc for your Matlab version / platform.

Inputs are easy... grid width, length, resolution, rotation-angle, lat/lon of corner, and a .mat file with depths (with their lat/lon).

Output are both, grid and initial-conditions .nc files... plus a screen display of what to copy-paste into the ocean.in file (Lm, Mn N) including DT calculated using Arango's formula posted earlier in this thread :)

The main purpose of this code is to help newbies (like me) pass the hump of setting up a "slightly realistic" application so that we can continue learning ROMS. However, if you plan to use this code for a more "realistic" application (which I might)... I strongly suggest you to (1) read the ASSUMPTIONS section in the code's help and (2) wait for people in this forum to try it out and criticize it. Here are a test_bathymetry file and a test_coastline file for you to play around.

On another note, I was able to see the results of my simulation (using the fixed grid) with help of another computer... because mine can't open those particular .nc output files. Since the forcing and parameterization were just a (non-realistic) test... the output is, non-surprisingly, completely non-realistic... yet, beautiful! I hope I can post some results soon.

Also, I should start redacting a WikiROMS post with what was learned up to now... this to give chance for things to be added (or taken out) before posting them on the Wiki.

RubenDiez-Lazaro

#28 Unread post by RubenDiez-Lazaro »

Diego:

Did you made any changes in the variable names definitions?
Seems i need the "varinfo_mgV3.dat" file in order to use same variable definitions that you....
Could do you upload it? Thanks

Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

#29 Unread post by Diego »

Ruben,

Here is the varinfo_mgV3.dat file.

Let me know if it works.

Diego

ocecept
Posts: 42
Joined: Tue Jan 08, 2008 3:57 pm
Location: Universidade Federal do Ceará
Contact:

#30 Unread post by ocecept »

Hi Diego;

First I would like to thank you for what you have been doing here.

I tried to use your roms_make_simple_grid.m to create a grid for South Australia, without success. My grid is larger (930 Km x 846 Km) than your test case and located in the south hemisphere.

The problem appears to be related with the Normal radius of curvature that results in a smaller value to the length of one degree of longitude.

If I make:

latdist = 111120;
equarad = 6378137;
polarad = 6356752.3142;
normradius = (equarad^2) ./ (sqrt(((equarad*cos(lat))^2) + ((polarad*sin(lat))^2)));
length_deg_long=(pi/180) * cos(lat) * normradius)

and I set the same value of latitude used in your test case (lat=44.812) the length of one degree of longitude is equal 75288.51, but try to set the lat equal 45 and you will get 58620.82. Compare with the values in http://en.wikipedia.org/wiki/Longitude and you will see that something is odd.

I was not able to find the formula to the Normal radius of curvature, so I changed the length of one degree of longitude [ (pi/180) * cos(lat) * normradius ] in your program by a fixed value that I get from the m_lldist (m_map package) and it works fine.

I'm not sure if I doing the right thing and I don't know how it will affect your grid.


Many thanks

Carlos

Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

#31 Unread post by Diego »

Carlos,

I'm glad this thread (and associated code) has been helpful...

In roms_make_simple_grid, I create the grid in a plane coordinate (in meters), then rotate the grid... and finally I estimated the lat/lon's for each grid point with respect to the bottom-left corner (user-specified parameter). I struggled to do this last step... but since Ruben was waiting for my reply, I decided do a quick fix and add a note in the .m file help stating the limitations:

"ASSUMPTIONS: This script should NOT be used to generate grids over large domains that span over many degrees of Latitude. The assumption of constant Coriolis over the grid domain and constant length of a degree of longitude (as a function of Latitude at the bottom-left corner of the grid), are violated if the grid spans over many degrees of Latitude."

Now that you found a way to make it work for big grids... I hope you don't mind if I include your fix in the code. However, I am not sure where you found m_lldist, I download the set of matlab tools from the "ROMS Processing Packages" page (matlab.tar.gz @ http://www.myroms.org/index.php?page=Processing )... I found the m_map directory but there is no m_lldist. Could you tell me where you found it? :?

Thanks a lot Carlos!

Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

#32 Unread post by Diego »

Aha! Found it!

M_Map is a mapping package for Matlab developed by R. Pawlowicz... you can find a instruction manual and download link at: http://www.eos.ubc.ca/~rich/map.html

It is bit awkward to find the different Matlab utilities that exist for ROMS (they are scattered around). You can visit WikiROMS Tools for a list of some of the utilities.... in there there is a link to ROMS_AGRIF / ROMSTOOLS which has another set of tools (including M_map).

I shall start code fixing...

ocecept
Posts: 42
Joined: Tue Jan 08, 2008 3:57 pm
Location: Universidade Federal do Ceará
Contact:

#33 Unread post by ocecept »

Hi Diego;

That's the m_map package. Sorry for the incomplete information.

If you REALLY want your program working for a large area with my solution (and m_map), change:

((pi/180) * cos(lat) * normradius) );

by

(1000*m_lldist([lon lon+1],[lat lat])));

at lines 145;163; 183 and 202

I hope the program will work. Maybe you will have a problem here since your dxdy WILL NOT be constant in the role grid for large grids.

To solve the problem with the Coriolis parameter:

just do:

OMEGA=7.29E-5;

for ii=1:Mp;
for jj=1:Lp;
f(ii,jj)= 2*OMEGA*sin(lat_rho(ii,jj).*pi/180.);
end
end


It will work!

Sorry about the extra work, but I think your program will be another solution to generate grids.

Many thanks

Carlos :D :D

Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

#34 Unread post by Diego »

I have been uploading tutorials to WikiROMS (e.g. https://www.myroms.org/wiki/index.php/MEXNC).

Things were going soothly with tutorial-writing until I ran into an interesting problem... I just realize that the Matlab files that I have been using to create/read/write/manipulate/plot ROMS files are NOT quite the same as the ones that are publicly available at http://www.myroms.org/index.php?page=Processing

The Matlab files that I use (which I got from Katja, and she got them from someone else) have the same name and ALMOST the same content as the publicly available files. The small differences are bug fixes that make the files that I have to work, while the publicly available simply keep crashing (my main problem being with editmask.m)

Maybe the "ROMS Processing Packages Directory Listing" is an old repository and we should not be getting files from there...

If that is the case, where should I send people (in the tutorials) to get their ROMS Matlab Tools?

:?

RubenDiez-Lazaro

#35 Unread post by RubenDiez-Lazaro »

Diego, Carlos:

I try the Diego's grid. It outputs some strange things: the time marks in the output of roms is changed by asterisks ('*') characters...

Very thanks for submit your matlab script. At first time your point of view is similar at mine.

I decide submit my code (mygrid_nc) for make grids in the roms way in the hope someone play with it and found what is broken: The grid generated seems ok, but it not work....

The code is based in the csa interpolation library and the gridgen grid generation library from Pavel Sakov
http://www.marine.csiro.au/~sak007/


In order to compile it, you need:

1- A working installation of the netcdf library (with development files).

2- Compile the csa code from Pavel and put the "libcsa.a" file where the compiler could find it (for example, /usr/local/lib/), and the "csa.h" file where the compiler could find it (for example, /usr/local/include/)

3- Compile the gridgen code from Pavel and put the "libgridgen.a" file where the compiler could find it (for example, /usr/local/lib/), and the "gridgen.h" file where the compiler could find it (for example, /usr/local/include/)

4- Get etopo2 and/or etopo5 batimetry data. etopo2 (in netcdf format) can be downloaded from the roms-tools page:
ftp://ftp.ifremer.fr/ifremer/ird/romsto ... 004.tar.gz
and the etopo5 from:
http://www.ngdc.noaa.gov/mgg/global/rel ... PO/ETOPO5/
(you need the etopo5.dat for big endian machines and the "etopo5.dos" for little endian machines).

5- Uncompress the file downloaded from:
http://web.usc.es/~rdlazaro/roms_stuff/ ... nc_c03.tgz
go to the "mygrid_nc_c03" directory and edit the "defs.h" in order to fix configuration stuff. Play attention to the paths to the etop2 and etopo5 files: put the values in your case.
See too the LAND_UMBRAL, MIN_DEPTH and LAND_UMBRAL values...
As last steep, type "make"...

USAGE:
./mygrid_nc_c [-p #.#] [-n #] [-o STR] [-s #] [-N] [-r #.#] [-t STR] in_file

OPTIONS:
[-p #.#] --> precision for gridgen
[-n #] --> newton for gridgen
[-o STR] --> out file name
[-s #] --> number of pass by 2D noise filter
[-N] --> do not pass by rfilter (for rvalue ensure)
[-r #.#] --> max rvalue allowed
[-t STR] --> batimetry font. Posible values:
"etopo2" --> NetCDF etopo2 batimetry
"etopo5" --> Raw etopo5 batimetry

NOTES:
- rfilter (fix for bad rvalues) not yet implemented.
- default batimetry data: etopo2
- The Examples.in dir content some examples of input data. The format for input data is a simple ascii text, equal to the format used by Pavel's gridgen...

EXAMPLE:
In the mygrid_nc_c03 directory type:

./mygrid_nc_c -o roms_grid_probe.nc -s 2 Examples.in/grid-l2.in

It uses the "Examples.in/grid-l2.in" input file for generate the grid. Output filename will be "roms_grid_probe.nc" ( -o option), and will make two passes by a gaussian filter ( -s option). The batimetry used will be the default one (etopo2).

-----------------
Please consider a beta version suitable only for beta-testers.
Some few comments in the code are still in spanish: sorry....

Diego
Posts: 36
Joined: Mon Dec 03, 2007 3:30 pm
Location: Dalhousie University - Dept. of Oceanography

#36 Unread post by Diego »

Hi All,

As promised, I summarized what I learned throughout this thread in 3 WikiROMS tutorials. Even though I did my best to test / check all I posted... I would appreciate if more senior ROMS users check this tutorials. I just don't want to spread wrong information. I also thought it may be useful to set a unique thread to discuss/comment each tutorial. A link in the forum thread will direct traffic to the tutorial (WikiROMS) and a link in the tutorial will redirect discussion back to the forum. Lets experiment and see if this works. Anyway... here are the tutorials:

(1) How to download/install MEXNC, SNCTOOLS and ROMS tool-kit for MATLAB
Post comments & questions for this tutorial: HERE

(2) Tutorial to grid-generation using EASYGRID
After I fixed my matlab script to create grids... I though it will be nice to share it with others. So I named it easygrid and I put it up in Google-code for free and open download.
Post comments & questions for this tutorial: HERE

(3) How to setup Tidal forcing in a Fjord
Post comments & questions for this tutorial: I am not sure if I should start a new thread or if I should direct discussion to an existing thread like THIS one.


Cheers!

Diego

hodamahyar
Posts: 16
Joined: Wed Jan 11, 2017 3:25 pm
Location: k.n.toosi University of Technology

Re: Application setup. Please help!

#37 Unread post by hodamahyar »

Hello.
first:i generate ROMS grid by its example in the file (Fjord_bathy.mat and Fjord_coast.m). but it show me the grid in 2D. how can i have 3D plot? anything should change in the romsgrid.mat code?
second:what is ocean.in file?where can i find it?
third:would you please send me mexnc,snctools,ROMS Matlab tool-kit,ocean.in files for me in a zip file?

Post Reply