I am attempting to run a one-year simulation of a modified version of the estuary test case:
1051200 ntimes
30.000 dt
20 ndtfast
I have modified bathymetry (ana_grid.h) as well as the tidal period from 12.0 to 12.42 (ana_fsobc.h). I have also attempted to add code (albeit inefficiently) to reflect SLR (ana_fsobc.h):
#elif defined YORK_ESTUARY
IF (LBC(ieast,isFsur,ng)%acquire.and. &
& DOMAIN(ng)%Eastern_Edge(tile)) THEN
cff=1.0_r8*SIN(2.0_r8*pi*time(ng)/(12.42_r8*3600.0_r8))
DO j=JstrR,JendR ! M2 Tide Period
BOUNDARY(ng)%zeta_east(j)=cff &
& +(((time(ng))/31536000_r8)*((4.44_r8)/1000_r8))
END DO ! SLR s/yr*mm/yr-m
END IF
When running at the original 10-day ESTUARY_TEST timescale (ntimes = 28800) zeta plots (MATLAB surfc) as expected with variation about MSL. However, with the 365-day run, zeta appears unstable as it becomes increasingly negative over the model domain. I have also run the longer simulation without my SLR code and the result is the same. Do I need to adjust other parameters to achieve stability or have I erred elsewhere?
APPARENT SIMULATION INSTABILITY
Re: APPARENT SIMULATION INSTABILITY
We developed the estuary test as a simple example, i think originally to compare along channel variations of mixing of salt and fresh water, back in the days of adding the GLS vertical mixing to the code. That test case has 2 boundary conditions. The east is just a river, with a uni-directional flow, and the west was a tide + the river.
If you want to ensure long-term mass flux through the system, then you need to be careful on the boundary conditions that you choose. I would suggest that you use the UV_Psource option on the upstream end, where there is just uni-directional flow. The uv_psource computes a mass flux on every barotropic time step, and ensures that the mass entering the system will be exactly what you prescribe. On the downstream end, you can use any of the m2 reduced or flather, which ever one works better. Mass will be conserved in the system. Can you give that a try?
If you want to ensure long-term mass flux through the system, then you need to be careful on the boundary conditions that you choose. I would suggest that you use the UV_Psource option on the upstream end, where there is just uni-directional flow. The uv_psource computes a mass flux on every barotropic time step, and ensures that the mass entering the system will be exactly what you prescribe. On the downstream end, you can use any of the m2 reduced or flather, which ever one works better. Mass will be conserved in the system. Can you give that a try?
Re: APPARENT SIMULATION INSTABILITY
Thank You John, will attempt to adjust bc's as recommended. Will re-post if I hit any snags along the way.
Re: APPARENT SIMULATION INSTABILITY
John,
I have tried the Flather and Reduced BCs with Psource at the western boundary (my river source); but, have a few additional questions.
In my Include file, I define both UV_COR and ANA_PSOURCE. However, my output file does not show that either of these C-preprocessing Options have been activated.
I specify analytically in ana_grid as follows:
#elif defined YORK_ESTUARY
Xsize=200000.0_r8 ! Short version Xsize=160km
Esize=3000.0_r8
depth=200.0_r8
f0=1.0E-04_r8 ! Mid-latitude f0=1.0E-04_r8
beta=0.0_r8
and quantify my river flow in ana_psource.h as follows:
#if defined YORK_ESTUARY
IF (Master.and.DOMAIN(ng)%SouthWest_Test(tile)) THEN
Nsrc=1
Dsrc(Nsrc)=0.0_r8
Isrc(Nsrc)=2
Jsrc(Nsrc)=1
END IF
and
# if defined YORK_ESTUARY
IF ((tdays(ng)-dstart).lt.0.5_r8) THEN
fac=1.0_r8+TANH((time(ng)-43200.0_r8)/43200.0_r8)
ELSE
fac=1.0_r8
END IF
DO is=1,Src
Qbar(is)=fac*667.0_r8 ! 667 m2/s
END DO
and
# if defined YORK_ESTUARY
IF (DOMAIN(ng)%NorthWest_Test(tile)) THEN ! JRA E to W
DO k=1,N(ng)
DO is=1,Nsrc
Tsrc(is,k,itemp)=T0(ng)
Tsrc(is,k,isalt)=0.0_r8
Tidal forcing at the eastern BC is specified analytically from ana_m2obc.h:
#if defined YORK_ESTUARY
IF (LBC(ieast,isUbar,ng)%acquire.and. &
& LBC(ieast,isVbar,ng)%acquire.and. &
& DOMAIN(ng)%Eastern_Edge(tile)) THEN
cff1=0.007_r8 ! JRA east end 0.4
! cff2=0.08_r8
! riv_flow=cff2*300.0_r8*0.05_r8 ! JRA depth was 5.0
tid_flow=cff1*3000.0_r8*200.0_r8 ! JRA depth was 10
my_area=0.0_r8
my_flux=0.0_r8
DO j=Jstr,Jend
cff=0.5_r8*(zeta(Istr ,j,knew)+h(Istr ,j)+ &
& zeta(Istr-1,j,knew)+h(Istr-1,j))/pn(Istr,j)
my_area=my_area+cff
END DO
my_flux=-tid_flow*SIN(2.0_r8*pi*time(ng)/ &
& (12.42_r8*3600.0_r8)) ! JRA M2 Tide Period
DO j=Jstr,Jend
BOUNDARY(ng)%ubar_east(j)=my_flux/my_area
BOUNDARY(ng)%vbar_east(j)=0.0_r8
END DO
END IF
1. Do I need to initialize all grid cells with f0 specified or is there another reason UV_COR is not reflected in my output file? Note, extreme f0 values don't appear to make a difference in my 400 x 6 grid tracer output after a 1 yr simulation at dt=30s.
2. Point source specified doesn't appear to impact tracer distribution (salinity or bed thickness) even with extreme values of source. Do I need to distribute Qbar across the western boundary? Assume Qbar is the depth integrated volume and would have to be reduced to reflect desired mass transport. Assume I define TS_PSOURCE only if I want to introduce suspended sediment flux from my river?
3. Finally, in my analytical tidal forcing is cff1 used elsewhere? I am somewhat confused why I shouldn't specific the tidal flow directly, e.g., tid_flow=4200.0_r8?
Current run is MPI with 1 x 16 tile if that is relevant. Appreciate your patience as always, believe it or not, I am getting better at this.
R/John
I have tried the Flather and Reduced BCs with Psource at the western boundary (my river source); but, have a few additional questions.
In my Include file, I define both UV_COR and ANA_PSOURCE. However, my output file does not show that either of these C-preprocessing Options have been activated.
I specify analytically in ana_grid as follows:
#elif defined YORK_ESTUARY
Xsize=200000.0_r8 ! Short version Xsize=160km
Esize=3000.0_r8
depth=200.0_r8
f0=1.0E-04_r8 ! Mid-latitude f0=1.0E-04_r8
beta=0.0_r8
and quantify my river flow in ana_psource.h as follows:
#if defined YORK_ESTUARY
IF (Master.and.DOMAIN(ng)%SouthWest_Test(tile)) THEN
Nsrc=1
Dsrc(Nsrc)=0.0_r8
Isrc(Nsrc)=2
Jsrc(Nsrc)=1
END IF
and
# if defined YORK_ESTUARY
IF ((tdays(ng)-dstart).lt.0.5_r8) THEN
fac=1.0_r8+TANH((time(ng)-43200.0_r8)/43200.0_r8)
ELSE
fac=1.0_r8
END IF
DO is=1,Src
Qbar(is)=fac*667.0_r8 ! 667 m2/s
END DO
and
# if defined YORK_ESTUARY
IF (DOMAIN(ng)%NorthWest_Test(tile)) THEN ! JRA E to W
DO k=1,N(ng)
DO is=1,Nsrc
Tsrc(is,k,itemp)=T0(ng)
Tsrc(is,k,isalt)=0.0_r8
Tidal forcing at the eastern BC is specified analytically from ana_m2obc.h:
#if defined YORK_ESTUARY
IF (LBC(ieast,isUbar,ng)%acquire.and. &
& LBC(ieast,isVbar,ng)%acquire.and. &
& DOMAIN(ng)%Eastern_Edge(tile)) THEN
cff1=0.007_r8 ! JRA east end 0.4
! cff2=0.08_r8
! riv_flow=cff2*300.0_r8*0.05_r8 ! JRA depth was 5.0
tid_flow=cff1*3000.0_r8*200.0_r8 ! JRA depth was 10
my_area=0.0_r8
my_flux=0.0_r8
DO j=Jstr,Jend
cff=0.5_r8*(zeta(Istr ,j,knew)+h(Istr ,j)+ &
& zeta(Istr-1,j,knew)+h(Istr-1,j))/pn(Istr,j)
my_area=my_area+cff
END DO
my_flux=-tid_flow*SIN(2.0_r8*pi*time(ng)/ &
& (12.42_r8*3600.0_r8)) ! JRA M2 Tide Period
DO j=Jstr,Jend
BOUNDARY(ng)%ubar_east(j)=my_flux/my_area
BOUNDARY(ng)%vbar_east(j)=0.0_r8
END DO
END IF
1. Do I need to initialize all grid cells with f0 specified or is there another reason UV_COR is not reflected in my output file? Note, extreme f0 values don't appear to make a difference in my 400 x 6 grid tracer output after a 1 yr simulation at dt=30s.
2. Point source specified doesn't appear to impact tracer distribution (salinity or bed thickness) even with extreme values of source. Do I need to distribute Qbar across the western boundary? Assume Qbar is the depth integrated volume and would have to be reduced to reflect desired mass transport. Assume I define TS_PSOURCE only if I want to introduce suspended sediment flux from my river?
3. Finally, in my analytical tidal forcing is cff1 used elsewhere? I am somewhat confused why I shouldn't specific the tidal flow directly, e.g., tid_flow=4200.0_r8?
Current run is MPI with 1 x 16 tile if that is relevant. Appreciate your patience as always, believe it or not, I am getting better at this.
R/John
Re: APPARENT SIMULATION INSTABILITY
well, the gov is back up and we are getting everything back up to speed here.
if you are making changes, and the output is not reflecting it, then you need to make sure you rebuilt it correctly and/or you are making changes to the correct files. make sure you are changing the ana* files in the right location.
if you are making changes, and the output is not reflecting it, then you need to make sure you rebuilt it correctly and/or you are making changes to the correct files. make sure you are changing the ana* files in the right location.