New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Bathymetry problems at 0.1 deg #99

Open
aekiss opened this Issue Jun 11, 2018 · 46 comments

Comments

Projects
None yet
5 participants
@aekiss
Contributor

aekiss commented Jun 11, 2018

I'm getting crashes like this from time to time:

FATAL from PE  3457:  Error: temperature out of range with value     3.081776472849E+02 at (i,j,k) = (2153,2093,  2),  (lon,lat,dpt) = (  -64.7500,   64.3310,    1.7294 m)

always at the same i,j=2153,2093 but with k=2 or 3. This is at the coast on the southern side of the mouth of Cumberland Sound on the east coast of Baffin Island, Canada. This is marked by the red circle in the figures below. It looks like /g/data3/hh5/tmp/cosima/bathymetry/topog_latest.nc has a large shallow estuary with a 1-cell mouth on Hall Peninsula, whereas in /g/data3/hh5/tmp/cosima/bathymetry/gebco.nc this is a series of islands and fjords.

I think we could just close off the "estuary". @russfiedler what do you think? Would you have a chance to look at this sometime soon?

screen shot 2018-06-12 at tue 12-6 9 08am

screen shot 2018-06-12 at tue 12-6 8 57am

screen shot 2018-06-12 at tue 12-6 8 58am 1

@PaulSpence

This comment has been minimized.

Show comment
Hide comment
@PaulSpence

PaulSpence Jun 12, 2018

PaulSpence commented Jun 12, 2018

@russfiedler

This comment has been minimized.

Show comment
Hide comment
@russfiedler

russfiedler Jun 12, 2018

@PaulSpence Yes, it should be removed. In fact, the point at i=2153, j=2093 is non advective and therefore illegal so nothing was entering nor leaving that estuary anyway except by horizontal diffusion. I'll need to check why it wasn't flagged.

russfiedler commented Jun 12, 2018

@PaulSpence Yes, it should be removed. In fact, the point at i=2153, j=2093 is non advective and therefore illegal so nothing was entering nor leaving that estuary anyway except by horizontal diffusion. I'll need to check why it wasn't flagged.

@russfiedler

This comment has been minimized.

Show comment
Hide comment
@russfiedler

russfiedler Jun 12, 2018

Ok, I've found a logic error with how non advective cells are detected in my code. There are special cases like the one above that are missed. There may be other cases too. Hold everything...

russfiedler commented Jun 12, 2018

Ok, I've found a logic error with how non advective cells are detected in my code. There are special cases like the one above that are missed. There may be other cases too. Hold everything...

@russfiedler

This comment has been minimized.

Show comment
Hide comment
@russfiedler

russfiedler Jun 12, 2018

Ouch! 10K potholes at depth. Never let me be in charge of road maintenance! New bathymetry at /g/data3/hh5/tmp/cosima/bathymetry/topog_12_06_2018.baffin.nc. I haven't updated the topog_latest.nc symlink yet in case it breaks things somewhere.

russfiedler commented Jun 12, 2018

Ouch! 10K potholes at depth. Never let me be in charge of road maintenance! New bathymetry at /g/data3/hh5/tmp/cosima/bathymetry/topog_12_06_2018.baffin.nc. I haven't updated the topog_latest.nc symlink yet in case it breaks things somewhere.

@PaulSpence

This comment has been minimized.

Show comment
Hide comment
@PaulSpence

PaulSpence Jun 12, 2018

PaulSpence commented Jun 12, 2018

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 12, 2018

Contributor

Thanks @russfiedler, that looks much better! I'll give it a whirl tomorrow.

Do you think I'll need to do anything to fix up the restarts (i.e. did you make any extra ocean)? Or should I just give it a go and see?

screen shot 2018-06-12 at tue 12-6 7 32pm

Contributor

aekiss commented Jun 12, 2018

Thanks @russfiedler, that looks much better! I'll give it a whirl tomorrow.

Do you think I'll need to do anything to fix up the restarts (i.e. did you make any extra ocean)? Or should I just give it a go and see?

screen shot 2018-06-12 at tue 12-6 7 32pm

@russfiedler

This comment has been minimized.

Show comment
Hide comment
@russfiedler

russfiedler Jun 13, 2018

No, there's no extra ocean anywhere. It's all filled. You may have to make eta_t and eta_t_bar zero over the new land points in ocean_barotropic.res.nc. The thickness restart may have to be updated but I think you can get away without doing that. I'll update the topog_latest.nc link to point to the new bathymetry.

russfiedler commented Jun 13, 2018

No, there's no extra ocean anywhere. It's all filled. You may have to make eta_t and eta_t_bar zero over the new land points in ocean_barotropic.res.nc. The thickness restart may have to be updated but I think you can get away without doing that. I'll update the topog_latest.nc link to point to the new bathymetry.

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 13, 2018

Contributor

The run is in the queue. I didn't tweak any restarts as the SSH in the spurious estuary is significantly positive (1-2m, presumably because it can't drain any accumulated runoff from any of its three advectively disconnected basins) so hopefully it will run - we'll see. There are two other new land points (one one the south side of Hall Peninsula and one on Victoria Island) that might also cause grief - see https://github.com/aekiss/notebooks/blob/master/bathymetry-closeup-baffin.ipynb
screen shot 2018-06-13 at wed 13-6 3 53pm

Contributor

aekiss commented Jun 13, 2018

The run is in the queue. I didn't tweak any restarts as the SSH in the spurious estuary is significantly positive (1-2m, presumably because it can't drain any accumulated runoff from any of its three advectively disconnected basins) so hopefully it will run - we'll see. There are two other new land points (one one the south side of Hall Peninsula and one on Victoria Island) that might also cause grief - see https://github.com/aekiss/notebooks/blob/master/bathymetry-closeup-baffin.ipynb
screen shot 2018-06-13 at wed 13-6 3 53pm

@russfiedler

This comment has been minimized.

Show comment
Hide comment
@russfiedler

russfiedler Jun 13, 2018

I knew about the point on the south side of Hall peninsula but didn't spot the othe.

I don't like that piece of coast near Victoria Island. Before or after! I wanted to avoid 2 cell inlets like that. If it leads (or has lead) to problems I'd be happy to sort it out.

russfiedler commented Jun 13, 2018

I knew about the point on the south side of Hall peninsula but didn't spot the othe.

I don't like that piece of coast near Victoria Island. Before or after! I wanted to avoid 2 cell inlets like that. If it leads (or has lead) to problems I'd be happy to sort it out.

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 13, 2018

Contributor

Is 2140,2318 another non-advective point still in /g/data3/hh5/tmp/cosima/bathymetry/topog_12_06_2018.baffin.nc?
It looks like it and eta_t seems anomalous also.
screen shot 2018-06-13 at wed 13-6 11 08pm
screen shot 2018-06-13 at wed 13-6 11 10pm 2

Contributor

aekiss commented Jun 13, 2018

Is 2140,2318 another non-advective point still in /g/data3/hh5/tmp/cosima/bathymetry/topog_12_06_2018.baffin.nc?
It looks like it and eta_t seems anomalous also.
screen shot 2018-06-13 at wed 13-6 11 08pm
screen shot 2018-06-13 at wed 13-6 11 10pm 2

@russfiedler

This comment has been minimized.

Show comment
Hide comment
@russfiedler

russfiedler Jun 14, 2018

No, it isn't non advective cell, however it is an isolated bay since there is no flow across the inlet. The code to remove these bays was only naively diffusing from the open ocean T cells and not checking if the pairs of points to the north and south-east/west (and other combinations) would allow flow through. I'll sort this fairly quickly.

russfiedler commented Jun 14, 2018

No, it isn't non advective cell, however it is an isolated bay since there is no flow across the inlet. The code to remove these bays was only naively diffusing from the open ocean T cells and not checking if the pairs of points to the north and south-east/west (and other combinations) would allow flow through. I'll sort this fairly quickly.

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 14, 2018

Contributor

OK thanks @russfiedler. In gebco it looks like it's basically closed so I think we should just fill it in. Here's a closeup:
screen shot 2018-06-14 at thu 14-6 12 11pm

Contributor

aekiss commented Jun 14, 2018

OK thanks @russfiedler. In gebco it looks like it's basically closed so I think we should just fill it in. Here's a closeup:
screen shot 2018-06-14 at thu 14-6 12 11pm

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 14, 2018

Contributor

I've found a couple more like that, will post soon.

Contributor

aekiss commented Jun 14, 2018

I've found a couple more like that, will post soon.

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 14, 2018

Contributor

This bay mouth at 2049,640 doesn't look advectively closed, but eta_t looks odd. From gebco it looks like the mouth should be open. Maybe just leave it as-is?
screen shot 2018-06-14 at thu 14-6 12 24pm
screen shot 2018-06-14 at thu 14-6 12 32pm 1
screen shot 2018-06-14 at thu 14-6 12 24pm 1

Contributor

aekiss commented Jun 14, 2018

This bay mouth at 2049,640 doesn't look advectively closed, but eta_t looks odd. From gebco it looks like the mouth should be open. Maybe just leave it as-is?
screen shot 2018-06-14 at thu 14-6 12 24pm
screen shot 2018-06-14 at thu 14-6 12 32pm 1
screen shot 2018-06-14 at thu 14-6 12 24pm 1

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 14, 2018

Contributor

This spot around 1573,1801 behind Vancouver Island also looks a bit odd in eta_t. Very constricted but not advectively closed, right?
screen shot 2018-06-14 at thu 14-6 12 48pm
screen shot 2018-06-14 at thu 14-6 12 46pm 1
screen shot 2018-06-14 at thu 14-6 12 48pm 1

Contributor

aekiss commented Jun 14, 2018

This spot around 1573,1801 behind Vancouver Island also looks a bit odd in eta_t. Very constricted but not advectively closed, right?
screen shot 2018-06-14 at thu 14-6 12 48pm
screen shot 2018-06-14 at thu 14-6 12 46pm 1
screen shot 2018-06-14 at thu 14-6 12 48pm 1

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 14, 2018

Contributor

@russfiedler I've run a check for other cells with non-advective open edges like the one at 2140,2318 above.

Turns out there are 102 of them (or perhaps more - I haven't done the i=0 edge) - see https://github.com/aekiss/notebooks/blob/master/non-advective.ipynb

The only ones that dam the only exit of an embayment are:

  • 2140,2318 (already noted above)
  • 3080,2599
  • 2123,2609
  • 2133,2624

So these should be fixed as the dammed bays might dry out or hold back runoff. They're all small so I think they should just be filled in as it makes restarts easier.

screen shot 2018-06-14 at thu 14-6 5 43pm
screen shot 2018-06-14 at thu 14-6 5 45pm
screen shot 2018-06-14 at thu 14-6 5 46pm

Contributor

aekiss commented Jun 14, 2018

@russfiedler I've run a check for other cells with non-advective open edges like the one at 2140,2318 above.

Turns out there are 102 of them (or perhaps more - I haven't done the i=0 edge) - see https://github.com/aekiss/notebooks/blob/master/non-advective.ipynb

The only ones that dam the only exit of an embayment are:

  • 2140,2318 (already noted above)
  • 3080,2599
  • 2123,2609
  • 2133,2624

So these should be fixed as the dammed bays might dry out or hold back runoff. They're all small so I think they should just be filled in as it makes restarts easier.

screen shot 2018-06-14 at thu 14-6 5 43pm
screen shot 2018-06-14 at thu 14-6 5 45pm
screen shot 2018-06-14 at thu 14-6 5 46pm

@russfiedler

This comment has been minimized.

Show comment
Hide comment
@russfiedler

russfiedler Jun 14, 2018

russfiedler commented Jun 14, 2018

@russfiedler

This comment has been minimized.

Show comment
Hide comment
@russfiedler

russfiedler Jun 15, 2018

Ok. bathymetry updated. Many/most of those 100 choked faces are deliberate. They aren't completely non advective but they do restrict flows in straits and around islands without the need to put land in.

russfiedler commented Jun 15, 2018

Ok. bathymetry updated. Many/most of those 100 choked faces are deliberate. They aren't completely non advective but they do restrict flows in straits and around islands without the need to put land in.

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 15, 2018

Contributor

Awesome, many thanks @russfiedler. I'll give it a go.

Contributor

aekiss commented Jun 15, 2018

Awesome, many thanks @russfiedler. I'll give it a go.

@aekiss aekiss changed the title from Bathymetry problem in Cumberland Sound (Baffin Island) at 0.1 deg to Bathymetry problems at 0.1 deg Jun 15, 2018

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 16, 2018

Contributor

I'm having trouble getting the new bathymetry running.

I set /short/v45/aek156/access-om2/input/mom_01deg/topog.nc to a copy of the new corrected topography /g/data3/hh5/tmp/cosima/bathymetry/topog_13_06_2018.baffin.nc and also used https://github.com/aekiss/notebooks/blob/master/fix-restarts.ipynb
to set eta_t and eta_t_bar to zero inside the new land points in /short/v45/aek156/access-om2/control/01deg_jra55_ryf/archive/restart387/ocean/ocean_barotropic.res.nc.00??

However this run crashed on initialization, with 7 cores giving Caught signal 8 (Floating point exception) with this backtrace

[r2693:446  :0] Caught signal 8 (Floating point exception)
==== backtrace ====
 2 0x000000000005a64c mxm_handle_error()  /var/tmp/OFED_topdir/BUILD/mxm-3.6.3104/src/mxm/util/debug/debug.c:641
 3 0x000000000005a7bc mxm_error_signal_handler()  /var/tmp/OFED_topdir/BUILD/mxm-3.6.3104/src/mxm/util/debug/debug.c:616
 4 0x0000000000032510 killpg()  ??:0
 5 0x000000000044f432 atmo_boundary_layer()  /short/v45/aek156/CHUCKABLE/aekiss/access-om2/src/cice5/build_auscom_3600x2700_2000p/ice_atmo.f90:314
 6 0x0000000000424d5d get_i2o_fluxes()  /short/v45/aek156/CHUCKABLE/aekiss/access-om2/src/cice5/build_auscom_3600x2700_2000p/cpl_forcing_handler.f90:558
 7 0x000000000040f6bf ice_step()  /short/v45/aek156/CHUCKABLE/aekiss/access-om2/src/cice5/build_auscom_3600x2700_2000p/CICE_RunMod.f90:367
 8 0x000000000040f6bf cice_run()  /short/v45/aek156/CHUCKABLE/aekiss/access-om2/src/cice5/build_auscom_3600x2700_2000p/CICE_RunMod.f90:180
 9 0x000000000040c37d icemodel()  /short/v45/aek156/CHUCKABLE/aekiss/access-om2/src/cice5/build_auscom_3600x2700_2000p/CICE.f90:57
10 0x000000000040c31e main()  ??:0
11 0x000000000001ed1d __libc_start_main()  ??:0
12 0x000000000040c229 _start()  ??:0
===================

(more details are in /short/v45/aek156/access-om2/control/01deg_jra55_ryf/work-7652169)

/short/v45/aek156/CHUCKABLE/aekiss/access-om2/src/cice5/build_auscom_3600x2700_2000p/ice_atmo.f90 lines 306-320 are

      do ij = 1, icells
         i = indxi(ij)
         j = indxj(ij)

      !------------------------------------------------------------
      ! define some more needed variables
      !------------------------------------------------------------

         TsfK       = Tsf(i,j) + Tffresh     ! surface temp (K)
         qsat       = qqq * exp(-TTT/TsfK)   ! saturation humidity (kg/m^3)
         ssq        = qsat / rhoa(i,j)       ! sat surf hum (kg/kg)

         thva(ij)   = potT(i,j) * (c1 + zvir * Qa(i,j)) ! virtual pot temp (K)
         delt(i,j)  = potT(i,j) - TsfK       ! pot temp diff (K)
         delq(i,j)  = Qa(i,j) - ssq          ! spec hum dif (kg/kg)

and the floating point exception happens at TsfK = Tsf(i,j) + Tffresh ! surface temp (K).
Hard to see how that line would be a problem unless Tsf(i,j) is a signalling nan?

I guess this means I need to tweak some of the cice restart files? eg do I need to update iceumask(nj, ni) in /short/v45/aek156/access-om2/control/01deg_jra55_ryf/archive/restart387/ice/iced.0036-01-01-00000.nc to exclude the newly-masked points? (is this needed so that indxi, indxj work correctly?)

Grateful for any advice, I'm really just guessing here.

Contributor

aekiss commented Jun 16, 2018

I'm having trouble getting the new bathymetry running.

I set /short/v45/aek156/access-om2/input/mom_01deg/topog.nc to a copy of the new corrected topography /g/data3/hh5/tmp/cosima/bathymetry/topog_13_06_2018.baffin.nc and also used https://github.com/aekiss/notebooks/blob/master/fix-restarts.ipynb
to set eta_t and eta_t_bar to zero inside the new land points in /short/v45/aek156/access-om2/control/01deg_jra55_ryf/archive/restart387/ocean/ocean_barotropic.res.nc.00??

However this run crashed on initialization, with 7 cores giving Caught signal 8 (Floating point exception) with this backtrace

[r2693:446  :0] Caught signal 8 (Floating point exception)
==== backtrace ====
 2 0x000000000005a64c mxm_handle_error()  /var/tmp/OFED_topdir/BUILD/mxm-3.6.3104/src/mxm/util/debug/debug.c:641
 3 0x000000000005a7bc mxm_error_signal_handler()  /var/tmp/OFED_topdir/BUILD/mxm-3.6.3104/src/mxm/util/debug/debug.c:616
 4 0x0000000000032510 killpg()  ??:0
 5 0x000000000044f432 atmo_boundary_layer()  /short/v45/aek156/CHUCKABLE/aekiss/access-om2/src/cice5/build_auscom_3600x2700_2000p/ice_atmo.f90:314
 6 0x0000000000424d5d get_i2o_fluxes()  /short/v45/aek156/CHUCKABLE/aekiss/access-om2/src/cice5/build_auscom_3600x2700_2000p/cpl_forcing_handler.f90:558
 7 0x000000000040f6bf ice_step()  /short/v45/aek156/CHUCKABLE/aekiss/access-om2/src/cice5/build_auscom_3600x2700_2000p/CICE_RunMod.f90:367
 8 0x000000000040f6bf cice_run()  /short/v45/aek156/CHUCKABLE/aekiss/access-om2/src/cice5/build_auscom_3600x2700_2000p/CICE_RunMod.f90:180
 9 0x000000000040c37d icemodel()  /short/v45/aek156/CHUCKABLE/aekiss/access-om2/src/cice5/build_auscom_3600x2700_2000p/CICE.f90:57
10 0x000000000040c31e main()  ??:0
11 0x000000000001ed1d __libc_start_main()  ??:0
12 0x000000000040c229 _start()  ??:0
===================

(more details are in /short/v45/aek156/access-om2/control/01deg_jra55_ryf/work-7652169)

/short/v45/aek156/CHUCKABLE/aekiss/access-om2/src/cice5/build_auscom_3600x2700_2000p/ice_atmo.f90 lines 306-320 are

      do ij = 1, icells
         i = indxi(ij)
         j = indxj(ij)

      !------------------------------------------------------------
      ! define some more needed variables
      !------------------------------------------------------------

         TsfK       = Tsf(i,j) + Tffresh     ! surface temp (K)
         qsat       = qqq * exp(-TTT/TsfK)   ! saturation humidity (kg/m^3)
         ssq        = qsat / rhoa(i,j)       ! sat surf hum (kg/kg)

         thva(ij)   = potT(i,j) * (c1 + zvir * Qa(i,j)) ! virtual pot temp (K)
         delt(i,j)  = potT(i,j) - TsfK       ! pot temp diff (K)
         delq(i,j)  = Qa(i,j) - ssq          ! spec hum dif (kg/kg)

and the floating point exception happens at TsfK = Tsf(i,j) + Tffresh ! surface temp (K).
Hard to see how that line would be a problem unless Tsf(i,j) is a signalling nan?

I guess this means I need to tweak some of the cice restart files? eg do I need to update iceumask(nj, ni) in /short/v45/aek156/access-om2/control/01deg_jra55_ryf/archive/restart387/ice/iced.0036-01-01-00000.nc to exclude the newly-masked points? (is this needed so that indxi, indxj work correctly?)

Grateful for any advice, I'm really just guessing here.

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 16, 2018

Contributor

Actually it's crashing in ice_step, ie in timestepping, not initialisation. I don't know whether it crashed on the first timestep but it didn't get to the end of day 1. Walltime Used: 00:06:30

Contributor

aekiss commented Jun 16, 2018

Actually it's crashing in ice_step, ie in timestepping, not initialisation. I don't know whether it crashed on the first timestep but it didn't get to the end of day 1. Walltime Used: 00:06:30

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 17, 2018

Contributor

I tried updating iceumask(nj, ni) in /short/v45/aek156/access-om2/control/01deg_jra55_ryf/archive/restart387/ice/iced.0036-01-01-00000.nc to exclude the newly-masked points but it crashed with exactly the same backtrace.

Contributor

aekiss commented Jun 17, 2018

I tried updating iceumask(nj, ni) in /short/v45/aek156/access-om2/control/01deg_jra55_ryf/archive/restart387/ice/iced.0036-01-01-00000.nc to exclude the newly-masked points but it crashed with exactly the same backtrace.

@russfiedler

This comment has been minimized.

Show comment
Hide comment
@russfiedler

russfiedler Jun 18, 2018

I can't see these. the 01deg_jra55_ryf directory is under the x77 group.

drwxr-s--- 459 aek156 x77 28672 Jun 15 00:32 01deg_jra55_ryf

You'll definitely need to update the Tmask as well. That's where fluxes are calculated.

russfiedler commented Jun 18, 2018

I can't see these. the 01deg_jra55_ryf directory is under the x77 group.

drwxr-s--- 459 aek156 x77 28672 Jun 15 00:32 01deg_jra55_ryf

You'll definitely need to update the Tmask as well. That's where fluxes are calculated.

@PaulSpence

This comment has been minimized.

Show comment
Hide comment
@PaulSpence

PaulSpence Jun 18, 2018

PaulSpence commented Jun 18, 2018

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 18, 2018

Contributor

Sorry @russfiedler, you should have access to /short/v45/aek156/access-om2/control/01deg_jra55_ryf/archive/restart387 now.
As far as I can see, the CICE restart /short/v45/aek156/access-om2/control/01deg_jra55_ryf/archive/restart387/ice/iced.0036-01-01-00000.nc only has a u mask, which surprised me. Or did you mean the Tmask somewhere else?

Contributor

aekiss commented Jun 18, 2018

Sorry @russfiedler, you should have access to /short/v45/aek156/access-om2/control/01deg_jra55_ryf/archive/restart387 now.
As far as I can see, the CICE restart /short/v45/aek156/access-om2/control/01deg_jra55_ryf/archive/restart387/ice/iced.0036-01-01-00000.nc only has a u mask, which surprised me. Or did you mean the Tmask somewhere else?

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 18, 2018

Contributor

Thanks @PaulSpence that sounds like a foolproof method. I'll do that if I can't get a simpler shortcut to work.

Contributor

aekiss commented Jun 18, 2018

Thanks @PaulSpence that sounds like a foolproof method. I'll do that if I can't get a simpler shortcut to work.

@russfiedler

This comment has been minimized.

Show comment
Hide comment
@russfiedler

russfiedler Jun 18, 2018

@aekiss There's a kmt.nc file that has a wet dry mask. I presume this is how cice knows where open water is. It still has the old land mask.

russfiedler commented Jun 18, 2018

@aekiss There's a kmt.nc file that has a wet dry mask. I presume this is how cice knows where open water is. It still has the old land mask.

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 18, 2018

Contributor

Ah OK thanks @russfiedler I'll fix that up.

Contributor

aekiss commented Jun 18, 2018

Ah OK thanks @russfiedler I'll fix that up.

@ofa001

This comment has been minimized.

Show comment
Hide comment
@ofa001

ofa001 Jun 18, 2018

ofa001 commented Jun 18, 2018

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 18, 2018

Contributor

That seems to have done the trick - it's timestepping now.

Contributor

aekiss commented Jun 18, 2018

That seems to have done the trick - it's timestepping now.

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Jun 30, 2018

Contributor

This has improved stability - it can now do every month with dt=540s :-)

For future reference, here's a summary of what was done to fix the topography and change it partway through a run:

  1. topog.nc fixed by @russfiedler
  • old bathymetry: /g/data3/hh5/tmp/cosima/bathymetry/topog_12_10_17_yenesei.nc
  • fixed version: /g/data3/hh5/tmp/cosima/bathymetry/topog_13_06_2018.baffin.nc
  • choked cell faces detected by non-advective.ipynb
  • non-advective cells, 10,000 potholes, and 4 choked cell faces removed
  • no new ocean cells added
  • further info on changes: bathymetry-closeup-baffin.ipynb
  1. used fix-restarts.ipynb to make these consistent with the new topog.nc so that we could restart from a previous run that used the old bathymetry:
  • mask in /short/v45/aek156/access-om2/input/mom_01deg/ocean_mask.nc (though it seemed to run fine without this)
  • iceumask in CICE restart file /g/data3/hh5/tmp/cosima/access-om2-01/01deg_jra55v13_ryf8485_spinup6/restart387/ice/iced.0036-01-01-00000.nc (not sure if this was really necessary)
  • kmt in CICE restart file /g/data3/hh5/tmp/cosima/access-om2-01/01deg_jra55v13_ryf8485_spinup6/restart387/ice/kmt.nc (fixes floating point exception)
  • eta_t and eta_t_bar in MOM restart files /short/v45/aek156/access-om2/control/01deg_jra55_ryf/archive/restart387/ocean/ocean_barotropic.res.nc.00{17,22,24} (to avoid "free surface penetrating rock" error)
Contributor

aekiss commented Jun 30, 2018

This has improved stability - it can now do every month with dt=540s :-)

For future reference, here's a summary of what was done to fix the topography and change it partway through a run:

  1. topog.nc fixed by @russfiedler
  • old bathymetry: /g/data3/hh5/tmp/cosima/bathymetry/topog_12_10_17_yenesei.nc
  • fixed version: /g/data3/hh5/tmp/cosima/bathymetry/topog_13_06_2018.baffin.nc
  • choked cell faces detected by non-advective.ipynb
  • non-advective cells, 10,000 potholes, and 4 choked cell faces removed
  • no new ocean cells added
  • further info on changes: bathymetry-closeup-baffin.ipynb
  1. used fix-restarts.ipynb to make these consistent with the new topog.nc so that we could restart from a previous run that used the old bathymetry:
  • mask in /short/v45/aek156/access-om2/input/mom_01deg/ocean_mask.nc (though it seemed to run fine without this)
  • iceumask in CICE restart file /g/data3/hh5/tmp/cosima/access-om2-01/01deg_jra55v13_ryf8485_spinup6/restart387/ice/iced.0036-01-01-00000.nc (not sure if this was really necessary)
  • kmt in CICE restart file /g/data3/hh5/tmp/cosima/access-om2-01/01deg_jra55v13_ryf8485_spinup6/restart387/ice/kmt.nc (fixes floating point exception)
  • eta_t and eta_t_bar in MOM restart files /short/v45/aek156/access-om2/control/01deg_jra55_ryf/archive/restart387/ocean/ocean_barotropic.res.nc.00{17,22,24} (to avoid "free surface penetrating rock" error)

@aekiss aekiss closed this Jun 30, 2018

aekiss added a commit to OceansAus/01deg_jra55_ryf that referenced this issue Jul 3, 2018

Fix topog.nc and use dt=450s for Jan year 36 onwards
/short/v45/aek156/access-om2/input/mom_01deg/topog.nc
was the same as
/g/data3/hh5/tmp/cosima/bathymetry/topog_12_10_17_yenesei.nc
until the end of run 387 (Dec year 35).

However this had potholes and problems on Baffin Island and elsewhere
OceansAus/access-om2#99

So from run 388 (Jan year 36) onwards
/short/v45/aek156/access-om2/input/mom_01deg/topog.nc
was changed to a copy of this corrected topography
/g/data3/hh5/tmp/cosima/bathymetry/topog_13_06_2018.baffin.nc

Restarts also needed tweaking - see README-restart387.txt
@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Aug 22, 2018

Contributor

@russfiedler - I've noticed that the topog.nc depth histograms are "gappy", ie continuous apart from particular depth ranges (perhaps between each cell bottom and cell bottom plus minimum partial cell?) where it seems depths have been adjusted:
depth-hist

If you look on a log scale this ends up giving quantised depths in shallow water, ie depths can only take a few different values: [plot updated]
depth-hist-log2

...so that topog.nc has terraces on shelves. This is particularly noticeable on the very broad Bering shelf:
bad_vorticity_topog

I may have misunderstood but I thought topog.nc was supposed to vary continuously in height, so that mom can represent the real topography as closely as possible with partial cells. If topog.nc is terraced like this, it seems to remove the benefit of partial cells in shallow water. Was topog.nc quantized onto the vertical grid in order to remove potholes...?

This may also have a bearing on issue #107.

Contributor

aekiss commented Aug 22, 2018

@russfiedler - I've noticed that the topog.nc depth histograms are "gappy", ie continuous apart from particular depth ranges (perhaps between each cell bottom and cell bottom plus minimum partial cell?) where it seems depths have been adjusted:
depth-hist

If you look on a log scale this ends up giving quantised depths in shallow water, ie depths can only take a few different values: [plot updated]
depth-hist-log2

...so that topog.nc has terraces on shelves. This is particularly noticeable on the very broad Bering shelf:
bad_vorticity_topog

I may have misunderstood but I thought topog.nc was supposed to vary continuously in height, so that mom can represent the real topography as closely as possible with partial cells. If topog.nc is terraced like this, it seems to remove the benefit of partial cells in shallow water. Was topog.nc quantized onto the vertical grid in order to remove potholes...?

This may also have a bearing on issue #107.

@aekiss aekiss reopened this Aug 22, 2018

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Aug 22, 2018

Contributor

oops, I needed more bins for the log histogram. This one shows the quantization in shallow water clearly:
depth-hist-log2

Contributor

aekiss commented Aug 22, 2018

oops, I needed more bins for the log histogram. This one shows the quantization in shallow water clearly:
depth-hist-log2

@russfiedler

This comment has been minimized.

Show comment
Hide comment
@russfiedler

russfiedler Aug 22, 2018

Good points.

  1. Partial cells can be no thinner than the top cell or a fraction of the cell. 25% from memory. This may be a partial explanation.
  2. Yes Pothole elimination will quantize the heights significantly as they are set to the no deeper than the deepest of the 4 surrounding kmu levels.
    I actually think it's the roughness in the number of levels that you're seeing rin the vorticity rather than than a quanization issue.

russfiedler commented Aug 22, 2018

Good points.

  1. Partial cells can be no thinner than the top cell or a fraction of the cell. 25% from memory. This may be a partial explanation.
  2. Yes Pothole elimination will quantize the heights significantly as they are set to the no deeper than the deepest of the 4 surrounding kmu levels.
    I actually think it's the roughness in the number of levels that you're seeing rin the vorticity rather than than a quanization issue.
@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Aug 22, 2018

Contributor

There's quantization in the underlying GEBCO data since it is at 1m resolution, but quantization in topog.nc is coarser and seems to be based on the vertical grid. Even with the restrictions some of the GEBCO depth levels should be showing up in topog.nc below ~40m.
depth-hist-topog-gebco

Contributor

aekiss commented Aug 22, 2018

There's quantization in the underlying GEBCO data since it is at 1m resolution, but quantization in topog.nc is coarser and seems to be based on the vertical grid. Even with the restrictions some of the GEBCO depth levels should be showing up in topog.nc below ~40m.
depth-hist-topog-gebco

@russfiedler

This comment has been minimized.

Show comment
Hide comment
@russfiedler

russfiedler Aug 22, 2018

Maybe I should have just looked at the metadata in the topography file than trying to remember it!

float depth(yy, xx) ;
	depth:missing_value = -1.e+30f ;
	depth:long_name = "depth" ;
	depth:units = "m" ;
	depth:lakes_removed = "yes" ;
	depth:minimum_depth = 10.43281f ;
	depth:minimum_levels = 7 ;
	depth:min_thick = 10.f ;
	depth:min_frac = 0.2f ;

min_thick is the important one here. Partial cells don't kick in till till the vertical thicknesses become greater than 10m which is at about 100m.

russfiedler commented Aug 22, 2018

Maybe I should have just looked at the metadata in the topography file than trying to remember it!

float depth(yy, xx) ;
	depth:missing_value = -1.e+30f ;
	depth:long_name = "depth" ;
	depth:units = "m" ;
	depth:lakes_removed = "yes" ;
	depth:minimum_depth = 10.43281f ;
	depth:minimum_levels = 7 ;
	depth:min_thick = 10.f ;
	depth:min_frac = 0.2f ;

min_thick is the important one here. Partial cells don't kick in till till the vertical thicknesses become greater than 10m which is at about 100m.

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Aug 23, 2018

Contributor

Ah OK thanks. I guess that min_thick value was inherited from the previous vertical grid, which had 10m levels at the surface?

Contributor

aekiss commented Aug 23, 2018

Ah OK thanks. I guess that min_thick value was inherited from the previous vertical grid, which had 10m levels at the surface?

@aidanheerdegen

This comment has been minimized.

Show comment
Hide comment
@aidanheerdegen

aidanheerdegen Aug 23, 2018

Contributor

I thought we had decided min_thick should be the thickness of the top layer?

So in this case would be about 1m.

Contributor

aidanheerdegen commented Aug 23, 2018

I thought we had decided min_thick should be the thickness of the top layer?

So in this case would be about 1m.

@russfiedler

This comment has been minimized.

Show comment
Hide comment
@russfiedler

russfiedler Aug 23, 2018

I'm surprised that I chose the min_thick to be as large as 10m. I don't think chosing the thickness to be 1 would make too much differenc (it's essentially redundant below 30m). I believe the main problem here is that the topography is rough in the levels sense.

I'm plotting where you have jumps of >1 in X and Y in the number of levels. This has nothing to do with partial cells and we can see that the rough patches pretty much align with what you're seeing in the vorticity.

x_roughness
y_roughness

russfiedler commented Aug 23, 2018

I'm surprised that I chose the min_thick to be as large as 10m. I don't think chosing the thickness to be 1 would make too much differenc (it's essentially redundant below 30m). I believe the main problem here is that the topography is rough in the levels sense.

I'm plotting where you have jumps of >1 in X and Y in the number of levels. This has nothing to do with partial cells and we can see that the rough patches pretty much align with what you're seeing in the vorticity.

x_roughness
y_roughness

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Aug 30, 2018

Contributor

Just noting here a correction to #99 (comment), from @russfiedler today - the GEBCO bathy is averaged over each MOM cell so this average is very much smoother than the 1m resolution of GEBCO. This means we will get partial cells right up to the 2nd layer if we set min_thick=1.1m to match the top layer.

Contributor

aekiss commented Aug 30, 2018

Just noting here a correction to #99 (comment), from @russfiedler today - the GEBCO bathy is averaged over each MOM cell so this average is very much smoother than the 1m resolution of GEBCO. This means we will get partial cells right up to the 2nd layer if we set min_thick=1.1m to match the top layer.

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Sep 5, 2018

Contributor

/g/data3/hh5/tmp/cosima/bathymetry/topog_05_09_2018_1m_partial.nc looks way better, thanks @russfiedler

screen shot 2018-09-05 at wed 5-9 5 22pm

screen shot 2018-09-05 at wed 5-9 5 21pm

Contributor

aekiss commented Sep 5, 2018

/g/data3/hh5/tmp/cosima/bathymetry/topog_05_09_2018_1m_partial.nc looks way better, thanks @russfiedler

screen shot 2018-09-05 at wed 5-9 5 22pm

screen shot 2018-09-05 at wed 5-9 5 21pm

@aekiss

This comment has been minimized.

Show comment
Hide comment
@aekiss

aekiss Sep 7, 2018

Contributor

@russfiedler would you be able to determine roughly how many new wet cells would be created with this bathymetry? These would require filling in with some sort of extrapolation if we want to continue an existing run with this bathy.

It looks like lots of locations are made deeper, and a smaller number get shallower:
screen shot 2018-09-07 at fri 7-9 1 41pm

screen shot 2018-09-07 at fri 7-9 1 23pm

Contributor

aekiss commented Sep 7, 2018

@russfiedler would you be able to determine roughly how many new wet cells would be created with this bathymetry? These would require filling in with some sort of extrapolation if we want to continue an existing run with this bathy.

It looks like lots of locations are made deeper, and a smaller number get shallower:
screen shot 2018-09-07 at fri 7-9 1 41pm

screen shot 2018-09-07 at fri 7-9 1 23pm

@russfiedler

This comment has been minimized.

Show comment
Hide comment
@russfiedler

russfiedler Sep 10, 2018

@aekiss I would estimate 232838 wet T cells were added and 231862 wet u cells were added (roughly...). No cells are filled in.

See topog_05_09_2018_1m_partial_levs.nc and topog_13_06_2018.baffin_levs.nc i /g/data3/hh5/tmp/cosima/bathymetry

russfiedler commented Sep 10, 2018

@aekiss I would estimate 232838 wet T cells were added and 231862 wet u cells were added (roughly...). No cells are filled in.

See topog_05_09_2018_1m_partial_levs.nc and topog_13_06_2018.baffin_levs.nc i /g/data3/hh5/tmp/cosima/bathymetry

@aidanheerdegen

This comment has been minimized.

Show comment
Hide comment
@aidanheerdegen

aidanheerdegen Sep 17, 2018

Contributor

Apologies if these have been covered before, but when testing the new fast collation on tenth output I noticed the following:

At 1928,2586 (-87.15, 85.2) it is only 2 cells wide, which might be contributing the to large salinity difference (9 in the bay, 25 outside (green colour)

screen shot 2018-09-17 at 4 25 52 pm
screen shot 2018-09-17 at 4 27 59 pm

Contributor

aidanheerdegen commented Sep 17, 2018

Apologies if these have been covered before, but when testing the new fast collation on tenth output I noticed the following:

At 1928,2586 (-87.15, 85.2) it is only 2 cells wide, which might be contributing the to large salinity difference (9 in the bay, 25 outside (green colour)

screen shot 2018-09-17 at 4 25 52 pm
screen shot 2018-09-17 at 4 27 59 pm

@aidanheerdegen

This comment has been minimized.

Show comment
Hide comment
@aidanheerdegen

aidanheerdegen Sep 17, 2018

Contributor

Maybe these are examples reproduce the real world well because there is little exchange through these straits? Don't know, but thought it was worth pointing out.

Another one, narrow strait at 2011,2394 (-78.85,77.1) not such a large difference in salinity (red is 32.5, green is 26.5, same scale as above)

screen shot 2018-09-17 at 4 28 53 pm

Contributor

aidanheerdegen commented Sep 17, 2018

Maybe these are examples reproduce the real world well because there is little exchange through these straits? Don't know, but thought it was worth pointing out.

Another one, narrow strait at 2011,2394 (-78.85,77.1) not such a large difference in salinity (red is 32.5, green is 26.5, same scale as above)

screen shot 2018-09-17 at 4 28 53 pm

@aidanheerdegen

This comment has been minimized.

Show comment
Hide comment
@aidanheerdegen

aidanheerdegen Sep 17, 2018

Contributor

The thin inlet has 2 or 3 choke points that are only 2 cells wide by the looks of it: 1703,2403 (-109.85,77.47). Could be contributing to the 5psu difference between the far end and the bay it comes off.

screen shot 2018-09-17 at 4 33 26 pm

Contributor

aidanheerdegen commented Sep 17, 2018

The thin inlet has 2 or 3 choke points that are only 2 cells wide by the looks of it: 1703,2403 (-109.85,77.47). Could be contributing to the 5psu difference between the far end and the bay it comes off.

screen shot 2018-09-17 at 4 33 26 pm

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment