[Dart-dev] DART/trunk Revision: 10818

[dart at ucar.edu]

thoar at ucar.edu
2016-12-19 14:39:50 -0700 (Mon, 19 Dec 2016)
512
Incorporating solutions from Dec 15 in response to crashing on yellowstone:

romain escudier <re217 at envsci.rutgers.edu>
AttachmentsDec 15 (4 days ago)

to me, dart, Enrique, Raphael 
Hi Tim,

We finally managed to make it work on yellowstone! I corrected a few things 
on model_mod.f90 (new version attached to the email) to address the problems 
you mentioned. It is also working with debug mode and bound-check options. 
When you have time, can you look at it to include it in the next revision of DART code.




Modified: DART/trunk/models/ROMS/model_mod.f90
===================================================================
--- DART/trunk/models/ROMS/model_mod.f90	2016-12-19 18:52:00 UTC (rev 10817)
+++ DART/trunk/models/ROMS/model_mod.f90	2016-12-19 21:39:50 UTC (rev 10818)
@@ -4071,6 +4071,9 @@
 !> (2)after we find the box which is closest to the given point,
 !>    we can calculate distances to the cells within the box.
 !> bjchoi 2014/08/07
+!>
+!> @todo This WHOLE ROUTINE can be replaced with something similar 
+!>       (and much faster) from POP.
 !-------------------------------------------------------------
 
 subroutine get_reg_box_indices(lon, lat, LON_al, LAT_al, x_ind, y_ind, istatus)
@@ -4110,10 +4113,13 @@
 ! boxLON and boxLAT
 allocate( boxLON(nxbox,nybox) )
 allocate( boxLAT(nxbox,nybox) )
-allocate( rtemp(Nx*Ny) )
-allocate( itemp(Nx*Ny) )
-allocate( jtemp(Nx*Ny) )
 
+! Temporary variables that are only used on the space of the boxes (nxbox*nybox)
+! or the space of 4 boxes ((2*num_cell+1)*(2*num_cell+1))
+allocate( rtemp(max(nxbox*nybox,(2*num_cell+1)*(2*num_cell+1))) )
+allocate( itemp(max(nxbox*nybox,(2*num_cell+1)*(2*num_cell+1))) )
+allocate( jtemp(max(nxbox*nybox,(2*num_cell+1)*(2*num_cell+1))) )
+
 rtemp = huge(rtemp)   ! ultimately want the minimum value
 itemp = -1            ! should not matter
 jtemp = -1            ! should not matter
@@ -4144,18 +4150,9 @@
 jj = nybox
 j = (num_cell/2) + (jj-1)*num_cell
 if (j .GT. Ny) j = Ny
+boxLON(ii,jj) = LON_al(i,j)
+boxLAT(ii,jj) = LAT_al(i,j)
 
-!> @todo boxLON(nxbox,nybox) and boxLAT(nxbox,nybox) are not getting set
-!> TJH          37          17   302.859687426721        52.1531690036829     
-!> TJH          37          18   301.879518471865        53.1934976695757     
-!> TJH          37          19   1.797693134862316E+308  1.797693134862316E+308
-
-!> @todo hardcode a fix for now - I know this is not correct with tilted domains
-boxLon(nxbox,nybox) = boxLon(nxbox,nybox-1) - &
-                     (boxLon(nxbox,nybox-1) - boxLon(nxbox,nybox-2))
-boxLat(nxbox,nybox) = boxLat(nxbox,nybox-1) + &
-                     (boxLat(nxbox,nybox-1) - boxLat(nxbox,nybox-2))
-
 !TJH: these have no need to be inside the loops.
 
 latrad = lat*DEG2RAD
@@ -4167,7 +4164,6 @@
 t = 0
 do ii=1,nxbox
 do jj=1,nybox
-!   write(*,*)'TJH', nxbox, nybox, ii, jj, boxLON(ii,jj), boxLAT(ii,jj)
     lon_dif = (lon - boxLON(ii,jj)) * part1 - part2
     lat_dif = (lat - boxLat(ii,jj))*111.13295_r8
     t = t+1
@@ -4186,10 +4182,10 @@
 ! the neighboring 4 model grid points are within this range.
 ii = itemp(loc_box)
 jj = jtemp(loc_box)
-i_start = max((ii-1)*num_cell - num_cell/2, 1) !
-i_end   = min(    ii*num_cell + num_cell/2,Nx) ! span 2 cells
-j_start = max((jj-1)*num_cell - num_cell/2, 1)
-j_end   = min(    jj*num_cell + num_cell/2,Ny)
+i_start = max(ii*num_cell - num_cell, 1) !
+i_end   = min(ii*num_cell + num_cell,Nx) ! span 2 cells
+j_start = max(jj*num_cell - num_cell, 1)
+j_end   = min(jj*num_cell + num_cell,Ny)
 
 ! now, calculate the distance between the obs and close candiate points
 
@@ -4212,68 +4208,63 @@
 !find minimum distances
 min_location = minloc(rtemp(1:t))
 loc1 = min_location(1)
+ii = itemp(loc1)
+jj = jtemp(loc1)
 
-!> @todo FIXME dies here with bounds checking turned on and you attempt
-!>       to interpolate something outside the domain. This WHOLE ROUTINE
-!>       can be replaced with something similar from POP.
-
-if ( (jtemp(loc1)+1) > size(LON_al,2) ) then
-   write(string1,*)'ERROR in finding matching grid point for'
-   write(string2,*)'longitude ',lon,' latitude ',lat
-   write(string3,*)'jtemp(loc1)+1 ',jtemp(loc1)+1,' is > size(LON_al,2) ',size(LON_al,2)
-   call error_handler(E_ERR,'get_reg_box_indices:', string1, &
-              source, revision, revdate, text2=string2, text3=string3)
+if ( ii == 1 .or. ii == Nx .or. jj == 1 .or. jj == Ny ) then
+   istatus = 1
+   return
 endif
 
-!check point location
-tp1 = checkpoint(LON_al(itemp(loc1)  ,jtemp(loc1)  ), &
-                 LAT_al(itemp(loc1)  ,jtemp(loc1)  ), &
-                 LON_al(itemp(loc1)+1,jtemp(loc1)+1), &