Nested 1/150° configuration crashes with neg dp

[Alfredo Terrazas]

Hello everyone,

I am working on a very high-resolution HYCOM configuration in the southern Gulf of Mexico with a horizontal resolution of 1/150°. This is the first time I have worked at such a high resolution, and I am encountering some difficulties with nesting.

My base experiment (expt_01.0), which runs without nesting, completes successfully (1-year run). However, when I run the nested configuration using boundary conditions from a parent simulation at 1/50° resolution, the model crashes after a few days with a "neg dp" error.

I have already experimented with reducing the time step, but the problem persists.

The model was compiled with:

setenv OCN_MISC "-DMASSLESS_1MM -DRDNEST_MASK -DLATBDT_NPLINE3"

I would like to ask specifically about the nesting mask. I am attaching the mask used for the nested domain and would appreciate any feedback on whether it looks reasonable or if there are obvious issues that could lead to negative layer thicknesses near the boundaries.

Has anyone experienced similar problems when nesting from 1/50° into 1/150°? Are there any particular checks you would recommend for the nesting mask or boundary treatment at this resolution?

I am also attaching the blkdat.input from the working non-nested experiment (expt_01.0) for reference and also from the nested experiment (expt_01.2).

Thank you for any suggestions.

[Alan Wallcraft]

The two blkdat.input's have differences other than the open boundary.

Run exact twins and see if the problem persists.

I suggest starting with momtyp=2, which is more stable that momtyp=4.

Alan.

[Alfredo Terrazas]

Hi Alan,

Thank you for your suggestion.

Following your recommendation, I performed a twin experiment (expt_02.0) so that the only intended difference from expt_01.0 is the use of the open boundary conditions.

The base experiment (expt_01.0) runs successfully with baroclinic/barotropic time steps of 100/2 s. Both expt_01.0 and expt_02.0 are already using momtyp=2.

For the nested experiment (expt_02.0), I found that these time steps were too large. To complete the first month, I had to reduce them to 40/1 s. Even then, during the second month the simulation crashed after only a few days with a neg dp error.

I also tested both the original (100/2 s) and the reduced (40/1 s) time steps for the second month, but the simulation failed in both cases.

I am attaching:

• the blkdat.input files for expt_01.0 and expt_02.0,

• the log files for months 1 and 2 of expt_02.0.

At this point, the only intended difference between the two experiments is the use of the open boundary conditions:

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

diff expt_02.0/blkdat.input expt_01.0/blkdat.input
6c6
<  020    'iexpt ' = experiment number x10
---
>  010    'iexpt ' = experiment number x10
121,122c121,122
<    0.125    'bnstfq' = number of days between baro nesting archive input
<    0.125    'nestfq' = number of days between 3-d  nesting archive input
---
>    0.0    'bnstfq' = number of days between baro nesting archive input
>    0.0    'nestfq' = number of days between 3-d  nesting archive input
124,125c124,125
<  40.0    'baclin' = baroclinic time step (seconds), int. divisor of 86400
<   1.0   'batrop' = barotropic time step (seconds), int. div. of baclin/2
---
>  100.0    'baclin' = baroclinic time step (seconds), int. divisor of 86400
>   2.0   'batrop' = barotropic time step (seconds), int. div. of baclin/2
229c229
<    2      'lbflag' = lateral barotropic bndy flag (0=none, 1=port, 2=input)
---
>    0      'lbflag' = lateral barotropic bndy flag (0=none, 1=port, 2=input)

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

 I would appreciate any suggestions on what else I should check.

Thank you.

[Alan Wallcraft]

Never reduce batrop, unless it is inconsistent with baclin, because a smaller value sis less stable.  So batrop=2 should work, as should baclin=90.0 (say).

Do you:

a) have the bathymetry in the relax zone interpolated from the coarse grid

b) retain the sea points "landward" of the open boundary in the nest archives

If so, this should work.

One thing you could try is retaining the 3-D nesting only: nestfq=0.125

Alan.

[Alfredo Terrazas]

Hi Alan,

Thank you for your suggestions.

I checked the bathymetry preparation again. The parent bathymetry was first interpolated from the coarse grid to the fine grid using isuba_topog. After that, I merged it with my higher-resolution bathymetry using topo_merge, preserving the parent bathymetry within the northern relaxation zone, as follows:

  IF     =     1,     1,
  IL     =  1020,  1020,
  JF     =   372,   356,
  JL     =   385,   371,
! ---      scale       - multiplier for 2nd bathymetry in box
! ---                     =-9.0; use boxscl linearly varying multiplier
! ---                     < 0.0; use 1st bathymetry, but merge land
! ---                     = 0.0; use 1st bathymetry only (default)
! ---                     = 1.0; use 2nd bathymetry only
! ---                     = 0.0-1.0; use fraction of each bathymetry
! ---                     = 2.0; use 2nd bathymetry where 1st is land
! ---                     = 3.0; use 2nd bathymetry where 1st is near land
  SCALE  =   1.0, -9.0,
! ---      boxscl(1,:) - scale factor for (if,jf) between 0.0 and 1.0
! ---      boxscl(2,:) - scale factor for (il,jf) between 0.0 and 1.0
! ---      boxscl(3,:) - scale factor for (il,jl) between 0.0 and 1.0
! ---      boxscl(4,:) - scale factor for (if,jl) between 0.0 and 1.0
  BOXSCL =   1.0,   1.0,   1.0,   1.0,
             0.0,   0.0,   1.0,   1.0,  
 &END
'E-o-D'

Also the relaxation mask was generated using the following regions:

   1      'if    ' = first i point of sub-region (<=0 to end)
1020      'il    ' = last  i point of sub-region
 372      'jf    ' = first j point of sub-region
 386      'jl    ' = last  j point of sub-region
   1.0    'efoldA' = bottom left  e-folding time in days
   1.0    'efoldB' = bottom right e-folding time in days
   0.1    'efoldC' = top    right e-folding time in days
   0.1    'efoldD' = top    left  e-folding time in days
   1      'if    ' = first i point of sub-region (<=0 to end)
1020      'il    ' = last  i point of sub-region
 356      'jf    ' = first j point of sub-region
 371      'jl    ' = last  j point of sub-region
   9.0    'efoldA' = bottom left  e-folding time in days
   9.0    'efoldB' = bottom right e-folding time in days
   1.0    'efoldC' = top    right e-folding time in days
   1.0    'efoldD' = top    left  e-folding time in days
  -1      'if    ' = first i point of sub-region (<=0 to end)
'E-o-D'
./rmu_linear

I also compared the merged bathymetry against the bathymetry obtained directly from isuba_topog. The attached figure shows the maximum bathymetry difference for each model row. As expected, the differences decrease through the transition zone (J=356–371) and become exactly zero within the relaxation zone adjacent to the open boundary (J=372–385). This suggests that the bathymetry being used is consistent with the regions defined in topo_merge and the relaxation mask

Finally, I tested the configuration using only 3-D nesting (nestfq=0.125), but unfortunately the model still crashes with the same neg dp error.

Do these topo_merge and relaxation mask definitions look consistent with what you had in mind? If so, do you have any suggestions on what else I should check?

Thank you very much for your help.

[Alan Wallcraft]

The relaxation zone is J=356–385 which is appropriate for 3x nesting, but your bathymetry is only from the outer model for J=372–385.

Extend the region identical to the outer bathymetry.  The transition zone can remain as 15 grid point, but it will be further into the domain.

Alan.

[Alfredo Terrazas]

Hi Alan,

Thank you for your suggestion.

I modified the bathymetry merge so that the bathymetry from the parent grid now extends over the entire relaxation zone (J=356–385), and I moved the 15-point transition farther into the interior (J=341–355), as you suggested.

I kept the relaxation mask unchanged, since it already covers J=356–385. After creating the new bathymetry, I regenerated the nesting archive files. 

Unfortunately, the model still crashes with the same neg dp error, after a few days of simulation. I also tried only with nestfq=0.125

Do you have any suggestions on what else I should investigate? Could the problem be related to the nesting archive files or the boundary forcing itself?