Changeset 8805


Ignore:
Timestamp:
2017-11-24T12:28:45+01:00 (3 years ago)
Author:
jchanut
Message:

Enhance4-freesurface. step 1: recover tracer conservation with split explicit free surface - #1959

Location:
branches/2017/dev_r8624_ENHANCE4_FREESURFACE/NEMOGCM/NEMO
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • branches/2017/dev_r8624_ENHANCE4_FREESURFACE/NEMOGCM/NEMO/NST_SRC/agrif_opa_update.F90

    r7646 r8805  
    500500         DO jj=j1,j2 
    501501            DO ji=i1,i2 
     502               ! Update corrective fluxes: 
     503               un_bf(ji,jj)  = un_bf(ji,jj) &  
     504                &           + tabres(ji,jj) / e2u(ji,jj)  - ub2_b(ji,jj) 
    502505               ub2_b(ji,jj) = tabres(ji,jj) / e2u(ji,jj) 
    503506            END DO 
     
    531534         DO jj=j1,j2 
    532535            DO ji=i1,i2 
     536               ! Update corrective fluxes: 
     537               vn_bf(ji,jj)  = vn_bf(ji,jj) &  
     538                &           + tabres(ji,jj) / e1v(ji,jj)  - vb2_b(ji,jj) 
    533539               vb2_b(ji,jj) = tabres(ji,jj) / e1v(ji,jj) 
    534540            END DO 
  • branches/2017/dev_r8624_ENHANCE4_FREESURFACE/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90

    r7753 r8805  
    209209            ! (used as a now filtered scale factor until the swap) 
    210210            ! ---------------------------------------------------- 
    211             IF( ln_dynspg_ts .AND. ln_bt_fw ) THEN    ! No asselin filtering on thicknesses if forward time splitting 
    212                e3t_b(:,:,1:jpkm1) = e3t_n(:,:,1:jpkm1) 
    213             ELSE 
    214                DO jk = 1, jpkm1 
    215                   e3t_b(:,:,jk) = e3t_n(:,:,jk) + atfp * ( e3t_b(:,:,jk) - 2._wp * e3t_n(:,:,jk) + e3t_a(:,:,jk) ) 
     211            DO jk = 1, jpkm1 
     212               e3t_b(:,:,jk) = e3t_n(:,:,jk) + atfp * ( e3t_b(:,:,jk) - 2._wp * e3t_n(:,:,jk) + e3t_a(:,:,jk) ) 
     213            END DO 
     214            ! Add volume filter correction: compatibility with tracer advection scheme 
     215            ! => time filter + conservation correction (only at the first level) 
     216            zcoef = atfp * rdt * r1_rau0 
     217            IF ( .NOT. ln_isf ) THEN   ! if no ice shelf melting 
     218               e3t_b(:,:,1) = e3t_b(:,:,1) - zcoef * ( emp_b(:,:) - emp(:,:) & 
     219                              &                      - rnf_b(:,:) + rnf(:,:) ) * tmask(:,:,1) 
     220            ELSE                     ! if ice shelf melting 
     221               DO jj = 1, jpj 
     222                  DO ji = 1, jpi 
     223                     ikt = mikt(ji,jj) 
     224                     e3t_b(ji,jj,ikt) = e3t_b(ji,jj,ikt) - zcoef * (  emp_b   (ji,jj) - emp   (ji,jj)  & 
     225                        &                                           - rnf_b   (ji,jj) + rnf   (ji,jj)  & 
     226                        &                                           + fwfisf_b(ji,jj) - fwfisf(ji,jj)  ) * tmask(ji,jj,ikt) 
     227                  END DO 
    216228               END DO 
    217                ! Add volume filter correction: compatibility with tracer advection scheme 
    218                ! => time filter + conservation correction (only at the first level) 
    219                zcoef = atfp * rdt * r1_rau0 
    220                IF ( .NOT. ln_isf ) THEN   ! if no ice shelf melting 
    221                   e3t_b(:,:,1) = e3t_b(:,:,1) - zcoef * ( emp_b(:,:) - emp(:,:) & 
    222                                  &                      - rnf_b(:,:) + rnf(:,:) ) * tmask(:,:,1) 
    223                ELSE                     ! if ice shelf melting 
    224                   DO jj = 1, jpj 
    225                      DO ji = 1, jpi 
    226                         ikt = mikt(ji,jj) 
    227                         e3t_b(ji,jj,ikt) = e3t_b(ji,jj,ikt) - zcoef * (  emp_b   (ji,jj) - emp   (ji,jj)  & 
    228                            &                                           - rnf_b   (ji,jj) + rnf   (ji,jj)  & 
    229                            &                                           + fwfisf_b(ji,jj) - fwfisf(ji,jj)  ) * tmask(ji,jj,ikt) 
    230                      END DO 
    231                   END DO 
    232                END IF 
    233             ENDIF 
     229            END IF 
    234230            ! 
    235231            IF( ln_dynadv_vec ) THEN      ! Asselin filter applied on velocity 
  • branches/2017/dev_r8624_ENHANCE4_FREESURFACE/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg_ts.F90

    r7831 r8805  
    10241024         un_adv(:,:) = zwx(:,:) * r1_hu_n(:,:) 
    10251025         vn_adv(:,:) = zwy(:,:) * r1_hv_n(:,:) 
     1026         ! 
     1027         un_bf(:,:) = 0._wp 
     1028         vn_bf(:,:) = 0._wp  
    10261029      ELSE 
    1027          un_adv(:,:) = z1_2 * ( ub2_b(:,:) + zwx(:,:) ) * r1_hu_n(:,:) 
    1028          vn_adv(:,:) = z1_2 * ( vb2_b(:,:) + zwy(:,:) ) * r1_hv_n(:,:) 
     1030         un_adv(:,:) = z1_2 * ( ub2_b(:,:) + zwx(:,:) - atfp * un_bf(:,:) ) * r1_hu_n(:,:) 
     1031         vn_adv(:,:) = z1_2 * ( vb2_b(:,:) + zwy(:,:) - atfp * vn_bf(:,:) ) * r1_hv_n(:,:) 
     1032         ! Update corrective fluxes for next time step: 
     1033         un_bf(:,:)  = atfp * un_bf(:,:) + (zwx(:,:) - ub2_b(:,:)) 
     1034         vn_bf(:,:)  = atfp * vn_bf(:,:) + (zwy(:,:) - vb2_b(:,:)) 
    10291035      END IF 
    10301036 
     
    12001206         CALL iom_get( numror, jpdom_autoglo, 'ub2_b'  , ub2_b  (:,:) )    
    12011207         CALL iom_get( numror, jpdom_autoglo, 'vb2_b'  , vb2_b  (:,:) )  
     1208         CALL iom_get( numror, jpdom_autoglo, 'un_bf'  , un_bf  (:,:) )    
     1209         CALL iom_get( numror, jpdom_autoglo, 'vn_bf'  , vn_bf  (:,:) )  
    12021210         IF( .NOT.ln_bt_av ) THEN 
    12031211            CALL iom_get( numror, jpdom_autoglo, 'sshbb_e'  , sshbb_e(:,:) )    
     
    12191227         CALL iom_rstput( kt, nitrst, numrow, 'ub2_b'   , ub2_b  (:,:) ) 
    12201228         CALL iom_rstput( kt, nitrst, numrow, 'vb2_b'   , vb2_b  (:,:) ) 
     1229         CALL iom_rstput( kt, nitrst, numrow, 'un_bf'   , un_bf  (:,:) ) 
     1230         CALL iom_rstput( kt, nitrst, numrow, 'vn_bf'   , vn_bf  (:,:) ) 
    12211231         ! 
    12221232         IF (.NOT.ln_bt_av) THEN 
  • branches/2017/dev_r8624_ENHANCE4_FREESURFACE/NEMOGCM/NEMO/OPA_SRC/DYN/sshwzv.F90

    r7753 r8805  
    252252      ENDIF 
    253253      !              !==  Euler time-stepping: no filter, just swap  ==! 
    254       IF(  ( neuler == 0 .AND. kt == nit000 ) .OR.    & 
    255          & ( ln_bt_fw    .AND. ln_dynspg_ts )      ) THEN  
     254      IF( neuler == 0 .AND. kt == nit000 ) THEN 
    256255         sshb(:,:) = sshn(:,:)                              ! before <-- now 
    257256         sshn(:,:) = ssha(:,:)                              ! now    <-- after  (before already = now) 
  • branches/2017/dev_r8624_ENHANCE4_FREESURFACE/NEMOGCM/NEMO/OPA_SRC/oce.F90

    r7646 r8805  
    4545   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::                              hvr_e  !: inverse of v-depth 
    4646   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ub2_b  , vb2_b           !: Half step fluxes (ln_bt_fw=T) 
     47   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   un_bf  , vn_bf           !: Asselin filter corrective fluxes (ln_bt_fw=T) 
    4748#if defined key_agrif 
    4849   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) ::   ub2_i_b, vb2_i_b         !: Half step time integrated fluxes  
     
    119120         &        hu_e(jpi,jpj),   hur_e(jpi,jpj),   hv_e(jpi,jpj),   hvr_e(jpi,jpj), STAT=ierr(5) ) 
    120121         ! 
    121       ALLOCATE( ub2_b(jpi,jpj), vb2_b(jpi,jpj)                                      , STAT=ierr(6) ) 
     122      ALLOCATE( ub2_b(jpi,jpj), vb2_b(jpi,jpj), un_bf(jpi,jpj), vn_bf(jpi,jpj)      , STAT=ierr(6) ) 
    122123#if defined key_agrif 
    123124      ALLOCATE( ub2_i_b(jpi,jpj), vb2_i_b(jpi,jpj)                                  , STAT=ierr(7) ) 
Note: See TracChangeset for help on using the changeset viewer.