Seems like the m file gcircle.m for grid refinement is not provided in the Rutgers matlab respository. Is this available or have i missed it?
thanks,
john
missing routine
Re: missing routine
also seems like G.uniform is undefined in bottom of grid_metrics.m
- arango
- Site Admin
- Posts: 1367
- Joined: Wed Feb 26, 2003 4:41 pm
- Location: DMCS, Rutgers University
- Contact:
Re: missing routine
The gcircle.m script is in the matlab/utility directory. The G.uniform should be set-up in get_roms_grid.m. Do you have the latest version of this script? or perhaps it is shadowed by an old version in another disrectory in your path. If not, it missed something about this. I always get the uniform structure field in idealized and realistic grid structure.
Re: missing routine
I've had issues with missing G.uniform flags being set too... if I remember correctly, it was in one of the nested grid setup utilities (coarse2fine.m, maybe?), though I don't have the code in front of me right now to check. My domain is non-idealized and non-uniform, for the record.
Re: missing routine
So i will take one, but that leaves one for you.
- I found gcircle in the Utility dir. did not have all the dirs in my path. user error.
- But the .uniform is missing in the structure. Here is what I did:
F=coarse2fine('USeast_grd19.nc','Car_nest5.nc',5,411,549,210,259);
This is to create the 5x child grid called Car_nest5.nc.
The call to coarse2fine goes well until it reaches line496:
[F.pm, F.pn, F.dndx, F.dmde]=grid_metrics(F, GreatCircle);
here we are going to grid_metrics to compute pm, pn etc. But in grid_metrics, we need at the bottom:
if (~G.uniform),
dndx(2:L,2:M) = 0.5.*(1.0./pn(3:Lp,2:M ) - 1.0./pn(1:Lm,2:M ));
dmde(2:L,2:M) = 0.5.*(1.0./pm(2:L ,3:Mp) - 1.0./pm(2:L ,1:Mm));
end
But G.uniform (really it is F.uniform) does not exist, because the F.uniform needed to be created in coarse2fine. But coarse2fine does not determine G.uniform. It is only determined in get_grid, but the creation of the Fine grid structure "F" does not call get_grid, it is created. So somewhere we need to compute F.uniform before the call to grid_metrics. I am not sure the best way to do this, but what i did was add to grid_metric:
%jcw add
if (~exist('G.uniform','var'))
G.uniform = false;
if (length(unique(pm(:))) == 1 && ...
length(unique(pn(:))) == 1),
G.uniform = true;
end
end
if (~G.uniform),
dndx(2:L,2:M) = 0.5.*(1.0./pn(3:Lp,2:M ) - 1.0./pn(1:Lm,2:M ));
dmde(2:L,2:M) = 0.5.*(1.0./pm(2:L ,3:Mp) - 1.0./pm(2:L ,1:Mm));
end
This worked, but i know you will do something else to fix this. I will move on with my testing for now.
-john
- I found gcircle in the Utility dir. did not have all the dirs in my path. user error.
- But the .uniform is missing in the structure. Here is what I did:
F=coarse2fine('USeast_grd19.nc','Car_nest5.nc',5,411,549,210,259);
This is to create the 5x child grid called Car_nest5.nc.
The call to coarse2fine goes well until it reaches line496:
[F.pm, F.pn, F.dndx, F.dmde]=grid_metrics(F, GreatCircle);
here we are going to grid_metrics to compute pm, pn etc. But in grid_metrics, we need at the bottom:
if (~G.uniform),
dndx(2:L,2:M) = 0.5.*(1.0./pn(3:Lp,2:M ) - 1.0./pn(1:Lm,2:M ));
dmde(2:L,2:M) = 0.5.*(1.0./pm(2:L ,3:Mp) - 1.0./pm(2:L ,1:Mm));
end
But G.uniform (really it is F.uniform) does not exist, because the F.uniform needed to be created in coarse2fine. But coarse2fine does not determine G.uniform. It is only determined in get_grid, but the creation of the Fine grid structure "F" does not call get_grid, it is created. So somewhere we need to compute F.uniform before the call to grid_metrics. I am not sure the best way to do this, but what i did was add to grid_metric:
%jcw add
if (~exist('G.uniform','var'))
G.uniform = false;
if (length(unique(pm(:))) == 1 && ...
length(unique(pn(:))) == 1),
G.uniform = true;
end
end
if (~G.uniform),
dndx(2:L,2:M) = 0.5.*(1.0./pn(3:Lp,2:M ) - 1.0./pn(1:Lm,2:M ));
dmde(2:L,2:M) = 0.5.*(1.0./pm(2:L ,3:Mp) - 1.0./pm(2:L ,1:Mm));
end
This worked, but i know you will do something else to fix this. I will move on with my testing for now.
-john
Re: missing routine
There is an error in the nesting. If i try a refinement ratio of 5, then there is an out of bounds error in nesting.F, :
"Fortran runtime error: Index '-4' of dimension 1 of array 'bry_contact' below lower bound of -3"
In correct_tracer_tile, this error is trying to access
TFF=TFF+BRY_CONTACT(iwest,dgcr)%Tflux(Jbf,k,itrc)
the variable Jbf is -4, but the lower allocated limit on it is -3.
The issue is that the tracer correction flux is being computed at one grid cell outside the parent boundary, and we dont need to be computing this correction flux here. We need to use the Jbc_min to limit where we start.
Here are the gory details:
We have:
TFF=0.0_r8
Jedge=Jo+(Jbc-Jbc_min)*RefineScale(ngf)
DO jsum=-half,half
Jbf=Jedge+jsum
DO k=1,N(ngf)
write(*,*) 'west ',jbf,k,itrc,jedge,jbc,jbc_min,jo
TFF=TFF+BRY_CONTACT(iwest,dgcr)%Tflux(Jbf,k,itrc)
Jbc actually starts out one less that Jbc_min, so we get
Jedge = 3 +(-1)*5 = -2
jsum = -2,2
Jbf=-2+-2 = -4 and this is too low.
This actually works for refinement of 3, as
Jedge = 2 +(-1)*3 = -1
jsum = -1,1
Jbf=-1+-1 = -2 and this is ok, but not really. we dont want to be comuting here.
We need to limit this some how such as replace
Jbc=BRY_CONTACT(iwest,rgcr)%Jb(j)
with
Jbc=max(BRY_CONTACT(iwest,rgcr)%Jb(j),Jbc_min)
or something like that.
I can do some more testing.
-john
"Fortran runtime error: Index '-4' of dimension 1 of array 'bry_contact' below lower bound of -3"
In correct_tracer_tile, this error is trying to access
TFF=TFF+BRY_CONTACT(iwest,dgcr)%Tflux(Jbf,k,itrc)
the variable Jbf is -4, but the lower allocated limit on it is -3.
The issue is that the tracer correction flux is being computed at one grid cell outside the parent boundary, and we dont need to be computing this correction flux here. We need to use the Jbc_min to limit where we start.
Here are the gory details:
We have:
TFF=0.0_r8
Jedge=Jo+(Jbc-Jbc_min)*RefineScale(ngf)
DO jsum=-half,half
Jbf=Jedge+jsum
DO k=1,N(ngf)
write(*,*) 'west ',jbf,k,itrc,jedge,jbc,jbc_min,jo
TFF=TFF+BRY_CONTACT(iwest,dgcr)%Tflux(Jbf,k,itrc)
Jbc actually starts out one less that Jbc_min, so we get
Jedge = 3 +(-1)*5 = -2
jsum = -2,2
Jbf=-2+-2 = -4 and this is too low.
This actually works for refinement of 3, as
Jedge = 2 +(-1)*3 = -1
jsum = -1,1
Jbf=-1+-1 = -2 and this is ok, but not really. we dont want to be comuting here.
We need to limit this some how such as replace
Jbc=BRY_CONTACT(iwest,rgcr)%Jb(j)
with
Jbc=max(BRY_CONTACT(iwest,rgcr)%Jb(j),Jbc_min)
or something like that.
I can do some more testing.
-john
- arango
- Site Admin
- Posts: 1367
- Joined: Wed Feb 26, 2003 4:41 pm
- Location: DMCS, Rutgers University
- Contact:
Re: missing routine
Yes, the uniform field in the grid structure was added recently. The fix in coarse2fine.m and fine2coarse.m is very simple. ITS value should be available and inherited from the donor grid.
Thank you for bringing this to my attention. Please update.
Thank you for bringing this to my attention. Please update.
- arango
- Site Admin
- Posts: 1367
- Joined: Wed Feb 26, 2003 4:41 pm
- Location: DMCS, Rutgers University
- Contact:
Re: missing routine
Okey, all my refinement test cases have a ratio of 1:3. I will need to have a test case with a 1:5 ratio to reproduce this problem and check in the debugger for a fix to the out-of-bounds error. It shouldn't that difficult to fix. This is only used in coupled routines.
I think that something is wrong in the formula to compute Ibc and Jbc. We always have the same number of contact points, regardless of what is the refinement ratio.
I will check on this after I create a 1:5 refinement application.
I think that something is wrong in the formula to compute Ibc and Jbc. We always have the same number of contact points, regardless of what is the refinement ratio.
I will check on this after I create a 1:5 refinement application.