WET_DRY and PERFECT_RESTART in Ice branch

Discussion about modeling ice with ROMS

Moderators: arango, robertson

Post Reply
Message
Author
Pysh
Posts: 30
Joined: Tue Nov 29, 2011 3:51 pm
Location: Hydrometcenter of Russia

WET_DRY and PERFECT_RESTART in Ice branch

#1 Unread post by Pysh »

I run model with flags
#define PERFECT_RESTART
#define WET_DRY
#define SOLVE3D
......................................
#define ICE_MODEL
#ifdef SOLVE3D
# ifdef ICE_MODEL
# define ICE_THERMO
# define ICE_MK
# undef ICE_ALB_EC92
# define ICE_SMOOTH
# define ICE_MOMENTUM
# define ICE_MOM_BULK
# define ICE_EVP
# define ICE_ADVECT
# define ICE_SMOLAR
# define ICE_UPWIND
# define ICE_BULK_FLUXES
# define ICE_SHOREFAST
# endif
#endif
After restart in the subroutine ice_advect ice concentrations are equal to 0

Code: Select all

!
! step number 1 in mpdata:
!
      DO j=Jmin,Jmax
        DO i=Imin,Imax
!
          ar(i,j)=1.0_r8/omn(i,j)
          aif(i,j)=(scr(i,j,liold)-dtice(ng)*(aflxu(i+1,j)-aflxu(i,j)   &
     &        +aflxv(i,j+1)-aflxv(i,j))*ar(i,j))
          aif(i,j) = aif(i,j)*rmask(i,j)
          aif(i,j) = aif(i,j)*rmask_wet(i,j)
        END DO
      END DO 
.........................................................
      DO j=Jstr,Jend
        DO i=Istr,Iend
            scr(i,j,linew) = aif(i,j)
        END DO
      END DO
because of rmask_wet(i,j) = 0 for all cells
In order to avoid this situation I've made changes
in file checkvars.F

Code: Select all

#ifdef PERFECT_RESTART
!
!  Determine perfect restart fields to read from input NetCDF file.
!
        IF (((model.eq.0).or.(model.eq.iNLM)).and.(nrrec(ng).ne.0)) THEN
# ifdef SOLVE3D
          get_var(idRu3d)=.TRUE.
          get_var(idRv3d)=.TRUE.

          get_var(idRwet)=.TRUE.    !!!!!!!!!!!  included
          get_var(idUwet)=.TRUE.    !!!!!!!!!!!  included
          get_var(idVwet)=.TRUE.    !!!!!!!!!!!  included

# endif
          get_var(idRzet)=.TRUE.
          get_var(idRu2d)=.TRUE.
          get_var(idRv2d)=.TRUE.
# if defined GLS_MIXING || defined MY25_MIXING
          get_var(idMtke)=.TRUE.
          get_var(idMtls)=.TRUE.
          get_var(idVmLS)=.TRUE.
          get_var(idVmKK)=.TRUE.
#  ifdef GLS_MIXING
          get_var(idVmKP)=.TRUE.
#  endif
# endif 
in file get_state.F

Code: Select all

............................................................................

!
!  Read in nonlinear free-surface (m).
!
        IF (get_var(idFsur)) THEN
          foundit=find_string(var_name, n_var, TRIM(Vname(1,idFsur)),   &
     &                        varid)
          IF (Perfect2D) THEN
            gtype=var_flag(varid)*r3dvar
          ELSE
            gtype=var_flag(varid)*r2dvar
          END IF
          IF (Perfect2D) THEN
            status=nf_fread3d(ng, IDmod, ncname, ncINPid,               &
     &                        Vname(1,idFsur), varid,                   &
     &                        InpRec, gtype, Vsize,                     &
     &                        LBi, UBi, LBj, UBj, 1, 3,                 &
     &                        Fscl, Fmin, Fmax,                         &
# ifdef MASKING
     &                        GRID(ng) % rmask,                         &
# endif
     &                        OCEAN(ng) % zeta)
          ELSE
            status=nf_fread2d(ng, IDmod, ncname, ncINPid,               &
     &                        Vname(1,idFsur), varid,                   &
     &                        InpRec, gtype, Vsize,                     &
     &                        LBi, UBi, LBj, UBj,                       &
     &                        Fscl, Fmin, Fmax,                         &
# ifdef MASKING
     &                        GRID(ng) % rmask,                         &
# endif
     &                        OCEAN(ng) % zeta(:,:,Tindex))
          END IF
          IF (status.ne.nf90_noerr) THEN
            IF (Master) THEN
              WRITE (stdout,30) string, TRIM(Vname(1,idFsur)), InpRec,  &
     &                          TRIM(ncname)
            END IF
            exit_flag=2
            ioerror=status
            RETURN
          ELSE
            IF (Master) THEN
              WRITE (stdout,70) TRIM(Vname(2,idFsur)), Fmin, Fmax
            END IF
          END IF
        END IF
# if defined MASKING   && \ 
  defined PERFECT_RESTART
  
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!  start of inclusion
!
!  Read in wet/dry mask in RHO-points.
!
        IF (get_var(idRwet)) THEN
          foundit=find_string(var_name, n_var, TRIM(Vname(1,idRwet)),   &
     &                        varid)
          gtype=var_flag(varid)*r2dvar
          status=nf_fread2d(ng, IDmod, ncname, ncINPid,               &
     &                        Vname(1,idRwet), varid,                   &
     &                        InpRec, gtype, Vsize,                     &
     &                        LBi, UBi, LBj, UBj,                       &
     &                        Fscl, Fmin, Fmax,                         &
     &                        GRID(ng) % rmask,                         &
     &                        GRID(ng) % rmask_wet)
          IF (status.ne.nf90_noerr) THEN
            IF (Master) THEN
              WRITE (stdout,30) string, TRIM(Vname(1,idRwet)), InpRec,  &
     &                          TRIM(ncname)
            END IF
            exit_flag=2
            ioerror=status
            RETURN
          ELSE
            IF (Master) THEN
              WRITE (stdout,70) TRIM(Vname(2,idRwet)), Fmin, Fmax
            END IF
          END IF
        END IF
!
!  Read in wet/dry mask in U-points.
!
        IF (get_var(idUwet)) THEN
          foundit=find_string(var_name, n_var, TRIM(Vname(1,idUwet)),   &
     &                        varid)
          gtype=var_flag(varid)*u2dvar
          status=nf_fread2d(ng, IDmod, ncname, ncINPid,               &
     &                        Vname(1,idUwet), varid,                   &
     &                        InpRec, gtype, Vsize,                     &
     &                        LBi, UBi, LBj, UBj,                       &
     &                        Fscl, Fmin, Fmax,                         &
     &                        GRID(ng) % umask,                         &
     &                        GRID(ng) % umask_wet)
          IF (status.ne.nf90_noerr) THEN
            IF (Master) THEN
              WRITE (stdout,30) string, TRIM(Vname(1,idUwet)), InpRec,  &
     &                          TRIM(ncname)
            END IF
            exit_flag=2
            ioerror=status
            RETURN
          ELSE
            IF (Master) THEN
              WRITE (stdout,70) TRIM(Vname(2,idUwet)), Fmin, Fmax
            END IF
          END IF
        END IF
!
!  Read in wet/dry mask in V-points.
!
        IF (get_var(idVwet)) THEN
          foundit=find_string(var_name, n_var, TRIM(Vname(1,idVwet)),   &
     &                        varid)
          gtype=var_flag(varid)*v2dvar
          status=nf_fread2d(ng, IDmod, ncname, ncINPid,               &
     &                        Vname(1,idVwet), varid,                   &
     &                        InpRec, gtype, Vsize,                     &
     &                        LBi, UBi, LBj, UBj,                       &
     &                        Fscl, Fmin, Fmax,                         &
     &                        GRID(ng) % vmask,                         &
     &                        GRID(ng) % vmask_wet)
          IF (status.ne.nf90_noerr) THEN
            IF (Master) THEN
              WRITE (stdout,30) string, TRIM(Vname(1,idVwet)), InpRec,  &
     &                          TRIM(ncname)
            END IF
            exit_flag=2
            ioerror=status
            RETURN
          ELSE
            IF (Master) THEN
              WRITE (stdout,70) TRIM(Vname(2,idVwet)), Fmin, Fmax
            END IF
          END IF
        END IF
 !!!!!!!!!!!!!!!!!!!!!!!!!!   end of inclusion
# endif
!
!  Read in nonlinear RHS of free-surface.
!
Then ice concentration (and others ice values) after restart save its values

In the history file some values are multiplied on rmask_io, umask_io, vmask_io. But theses masks are calculated in subroutine wetdry_tile (from step2d) after first output of history values (from subroutine output). To avoid this situation I've changed file set_masks.F

Code: Select all

...................................................................................
      CALL set_masks_tile (ng, tile, model,                             &
     &                     LBi, UBi, LBj, UBj,                          &
     &                     IminS, ImaxS, JminS, JmaxS,                  &
# ifdef UV_PSOURCE
     &                     Msrc(ng), Nsrc(ng),                          &
     &                     SOURCES(ng) % Isrc,                          &
     &                     SOURCES(ng) % Jsrc,                          &
     &                     SOURCES(ng) % Dsrc,                          &
# endif
     &                     GRID(ng) % pmask,                            &
     &                     GRID(ng) % rmask,                            &
     &                     GRID(ng) % umask,                            &
     &                     GRID(ng) % vmask,                            &

# if defined PERFECT_RESTART                                               !!!!  included
     &                     GRID(ng) % rmask_wet,                        &  !!!!  included
     &                     GRID(ng) % umask_wet,                        &  !!!!  included
     &                     GRID(ng) % vmask_wet,                        &  !!!!  included
# endif                                                                    !!!!  included

# ifdef OUTFLOW_MASK
     &                     GRID(ng) % mask_outflow,                     &
# endif
............................................................................................
     SUBROUTINE set_masks_tile (ng, tile, model,                       &
     &                           LBi, UBi, LBj, UBj,                    &
     &                           IminS, ImaxS, JminS, JmaxS,            &
# ifdef UV_PSOURCE
     &                           Msrc, Nsrc,                            &
     &                           Isrc, Jsrc, Dsrc,                      &
# endif
     &                           pmask, rmask,                          &
     &                           umask, vmask,                          &
# if defined PERFECT_RESTART                                                !!!!  included  
     &                           rmask_wet,                             &  !!!!  included
     &                           umask_wet, vmask_wet,                  &   !!!!  included
# endif                                                                    !!!!  included 

# ifdef OUTFLOW_MASK
     &                           mask_outflow,                          &
# endif
#
............................................................................................................
      real(r8), intent(in) :: pmask(LBi:,LBj:)
      real(r8), intent(in) :: rmask(LBi:,LBj:)
      real(r8), intent(in) :: umask(LBi:,LBj:)
      real(r8), intent(in) :: vmask(LBi:,LBj:)

# if defined PERFECT_RESTART                                   !!!!  included  
      real(r8), intent(in) :: rmask_wet(LBi:,LBj:)                !!!!  included  
      real(r8), intent(in) :: umask_wet(LBi:,LBj:)               !!!!  included  
      real(r8), intent(in) :: vmask_wet(LBi:,LBj:)               !!!!  included  
# endif                                                         !!!! included   
....................................................................

!!!!!!!!!!!!!!!  start of inclusion
# if defined PERFECT_RESTART
      rmask_io(LBi:UBi,LBj:UBj) = 0
      umask_io(LBi:UBi,LBj:UBj) = 0
      vmask_io(LBi:UBi,LBj:UBj) = 0

      DO j=Jstr,JendR
        DO i=Istr,IendR
          pmask_io(i,j)=pmask(i,j)
        END DO
      END DO
      DO j=JstrR,JendR
        DO i=IstrR,IendR
          rmask_io(i,j)=rmask(i,j) * rmask_wet(i,j)
        END DO
      END DO
      DO j=JstrR,JendR
        DO i=Istr,IendR
          umask_io(i,j)=umask(i,j) * umask_wet(i,j)
        END DO
      END DO
      DO j=Jstr,JendR
        DO i=IstrR,IendR
          vmask_io(i,j)=vmask(i,j) * vmask_wet(i,j)
        END DO
      END DO
# else
!!!!!!!!!!!!!!!  end of inclusion


      DO j=Jstr,JendR
        DO i=Istr,IendR
          pmask_io(i,j)=pmask(i,j)
        END DO
      END DO
      DO j=JstrR,JendR
        DO i=IstrR,IendR
          rmask_io(i,j)=rmask(i,j)
        END DO
      END DO
      DO j=JstrR,JendR
        DO i=Istr,IendR
          umask_io(i,j)=umask(i,j)
        END DO
      END DO
      DO j=Jstr,JendR
        DO i=IstrR,IendR
          vmask_io(i,j)=vmask(i,j)
        END DO
      END DO
# endif                                                 !!!! included

# ifdef OUTFLOW_MASK
And

Values Pair in input forcing file have units Pascal. In the model its values are multiplied by 0.01 ( for millibar units). But in history and average output files should be to return to Pascal units
I've made some corrections in files wrt_his.F and wrt_avg.F

Code: Select all

.....................................
# if defined BULK_FLUXES || defined ECOSIM || defined ATM_PRESS
!
!  Write out surface air pressure.
!
      IF (Hout(idPair,ng)) THEN
        scale=100.0_r8                    !!!!!! changed
        gtype=gfactor*r2dvar   
......................................    
.......................................
#  if defined BULK_FLUXES || defined ECOSIM || defined ATM_PRESS
!
!  Write out surface air pressure.
!
      IF (Aout(idPair,ng)) THEN
        scale=100.0_r8                   !!!!!! changed
        gtype=gfactor*r2dvar
..........................................
And. I'm not sure, but changes in file ice_mk.h perhaps is true

Code: Select all

...........................................
#ifdef CASPIAN_XXX
          hh = h(i,j)+Zt_avg1(i,j)
          IF (hh.LT.1.0_r8) THEN
            fac_shflx = hh
!            fac_shflx = 0.0_r8
          ELSE
            fac_shflx = 1.0_r8
          END IF
#else
          fac_shflx = 1.0_r8
#endif
#ifdef ICE_SHOREFAST             !!!! included
          fac_sf = 1.0_r8        !!!! included
#endif                           !!!! included

#ifdef ICESHELF
          IF (zice(i,j).eq.0.0_r8) THEN
#endif
            IF(ai(i,j,linew).LE.min_a(ng)) THEN
               stflx(i,j,itemp) = qao_n(i,j)*fac_shflx
            ELSE
#ifdef ICE_SHOREFAST  
............................................
Without this correction value fac_sf is not determined


And
I have not change, but in file set_avg.F compilator make warning about fac in line 2563. It is not determined. May be should be rfac(i,j) instead of fac?

..............................................

Code: Select all

       IF (Aout(idTauiw,ng)) THEN
          DO j=JstrR,JendR
            DO i=IstrR,IendR
              AVERAGE(ng)%avgutau_iw(i,j)=                              &
     &                      fac*AVERAGE(ng)%avgutau_iw(i,j)                      !!!!!!!!!!!!!!!!!!!!!!!!
            END DO
          END DO
          IF (EWperiodic(ng).or.NSperiodic(ng)) THEN
            CALL exchange_r2d_tile (ng, tile,                           &
     &                              LBi, UBi, LBj, UBj,                 &
     &                              AVERAGE(ng)%avgutau_iw)
#  ifdef DISTRIBUTE
            CALL mp_exchange2d (ng, tile, iNLM, 1,                      &
     &                          LBi, UBi, LBj, UBj,                     &
     &                          NghostPoints,                           &
     &                          EWperiodic(ng), NSperiodic(ng),         &
     &                          AVERAGE(ng)%avgutau_iw)
#  endif
          END IF
        END IF
................................................................
If my corections are not true excuse me, If true, will be usefull to know it :)
Now I run model with its

Thanks

Pysh
Posts: 30
Joined: Tue Nov 29, 2011 3:51 pm
Location: Hydrometcenter of Russia

Re: WET_DRY and PERFECT_RESTART in Ice branch

#2 Unread post by Pysh »

And

in file set_vbc.F

Code: Select all

!***********************************************************************
      SUBROUTINE set_vbc (ng, tile)
!***********************************************************************
!
      USE mod_param
      USE mod_grid
      USE mod_forces
      USE mod_ocean
      USE mod_stepping
...............................
#  ifndef BBL_MODEL
     &                   FORCES(ng) % bustr,                            &
     &                   FORCES(ng) % bvstr,                            &
#  endif
     &                   FORCES(ng) % stflx,                            &
     &                   FORCES(ng) % qao_n,                            &  !!!! included
     &                   FORCES(ng) % btflx)
#  ifdef PROFILE
      CALL wclock_off (ng, iNLM, 6)
#  endif

      RETURN
      END SUBROUTINE set_vbc
!
!***********************************************************************
      SUBROUTINE set_vbc_tile (ng, tile,                                &
     &                         LBi, UBi, LBj, UBj,                      &
     &                         IminS, ImaxS, JminS, JmaxS,              &
     &                         nrhs,                                    &
     &                         Hz,                                      &
...............................................
#  ifndef BBL_MODEL
     &                         bustr, bvstr,                            &
#  endif
     &                         stflx, qao_n, btflx)                  !!!! changed
............................................
#   ifndef BBL_MODEL
      real(r8), intent(inout) :: bustr(LBi:UBi,LBj:UBj)
      real(r8), intent(inout) :: bvstr(LBi:UBi,LBj:UBj)
#   endif
      real(r8), intent(inout) :: stflx(LBi:UBi,LBj:UBj,NT(ng))
      real(r8), intent(inout) :: btflx(LBi:UBi,LBj:UBj,NT(ng))
#  endif
      real(r8), intent(inout) :: qao_n(LBi:UBi,LBj:UBj)             !!!! included    
........................................................
#  ifdef QCORRECTION
!
!-----------------------------------------------------------------------
!  Add in flux correction to surface net heat flux (degC m/s).
!-----------------------------------------------------------------------
!
! Add in net heat flux correction.
!
      DO j=JstrR,JendR
        DO i=IstrR,IendR
          stflx(i,j,itemp)=stflx(i,j,itemp)+                            &
     &                     dqdt(i,j)*(t(i,j,N(ng),nrhs,itemp)-sst(i,j))
         qao_n(i,j) = -stflx(i,j,itemp) * rho0 * Cp                        !!!! included
        END DO
      END DO
#  endif

    
Otherwise after qcorrection in subroutine set_vbc correction will be droped in subroutine ice_thermo (from seaice)

Thanks

User avatar
kate
Posts: 4091
Joined: Wed Jul 02, 2003 5:29 pm
Location: CFOS/UAF, USA

Re: WET_DRY and PERFECT_RESTART in Ice branch

#3 Unread post by kate »

Thank you! :D I added PERFECT_RESTART to the ice code, but never checked it with WET_DRY. I will see about incorporating these changes.

Pysh
Posts: 30
Joined: Tue Nov 29, 2011 3:51 pm
Location: Hydrometcenter of Russia

Re: WET_DRY and PERFECT_RESTART in Ice branch

#4 Unread post by Pysh »

Hi

After last code correction some inaccuracies are appeared, I believe

in set_state.F

string 78
# if defined PERFECT_RESTART && defined WET_DRY
string 128
# if defined PERFECT_RESTART && defined WET_DRY
string 182
# if defined PERFECT_RESTART && defined WET_DRY
string 226
# if defined PERFECT_RESTART && defined WET_DRY
string 278
# if defined PERFECT_RESTART && defined WET_DRY

in set_avg.F
string 4838
AVERAGE(ng)%avgutau_iw(i,j)= &
& rfac(i,j)*AVERAGE(ng)%avgutau_iw(i,j)


Best regards

User avatar
kate
Posts: 4091
Joined: Wed Jul 02, 2003 5:29 pm
Location: CFOS/UAF, USA

Re: WET_DRY and PERFECT_RESTART in Ice branch

#5 Unread post by kate »

Thanks, you're right about set_avg.F. I don't have a set_state.F though, nor do your fixes match get_state.F. Could you please elaborate?

Pysh
Posts: 30
Joined: Tue Nov 29, 2011 3:51 pm
Location: Hydrometcenter of Russia

Re: WET_DRY and PERFECT_RESTART in Ice branch

#6 Unread post by Pysh »

Hi

In my remarks about inaccuracies there is an inaccuracy :D

Neither in the file set_state.F, nor in the file get_state.F. In the file set_masks.F

Sorry

Best regards

User avatar
kate
Posts: 4091
Joined: Wed Jul 02, 2003 5:29 pm
Location: CFOS/UAF, USA

Re: WET_DRY and PERFECT_RESTART in Ice branch

#7 Unread post by kate »

Indeed, thanks! My delay in responding is due to a week at Tai Chi camp in Canada.

Post Reply