Ignore:
Timestamp:
2017-05-09T17:36:25+02:00 (3 years ago)
Author:
jchanut
Message:

AGRIF vvl add on

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2017/dev_r7963_nemo_v3_6_AGRIF-3_AGRIFVVL/NEMOGCM/NEMO/NST_SRC/agrif_opa_update.F90

    r7988 r8010  
    3131CONTAINS 
    3232 
    33    RECURSIVE SUBROUTINE Agrif_Update_Tra( ) 
     33   SUBROUTINE Agrif_Update_Tra( ) 
    3434      !!--------------------------------------------- 
    3535      !!   *** ROUTINE Agrif_Update_Tra *** 
     
    6060      Agrif_UseSpecialValueInUpdate = .FALSE. 
    6161      ! 
    62       IF ( lk_agrif_doupd ) THEN ! Initialisation: do recursive update: 
    63          CALL Agrif_ChildGrid_To_ParentGrid() 
    64          CALL Agrif_Update_Tra() 
    65          CALL Agrif_ParentGrid_To_ChildGrid() 
    66       ENDIF 
    67       ! 
    6862#endif 
    6963      ! 
    7064   END SUBROUTINE Agrif_Update_Tra 
    7165 
    72    RECURSIVE SUBROUTINE Agrif_Update_Dyn( ) 
     66   SUBROUTINE Agrif_Update_Dyn( ) 
    7367      !!--------------------------------------------- 
    7468      !!   *** ROUTINE Agrif_Update_Dyn *** 
     
    145139#endif 
    146140      ! 
    147       ! Do recursive update: 
    148       IF ( lk_agrif_doupd ) THEN ! Initialisation: do recursive update: 
    149          CALL Agrif_ChildGrid_To_ParentGrid() 
    150          CALL Agrif_Update_Dyn() 
    151          CALL Agrif_ParentGrid_To_ChildGrid() 
    152       ENDIF 
    153       ! 
    154141   END SUBROUTINE Agrif_Update_Dyn 
    155142 
     
    179166# endif /* key_zdftke */ 
    180167 
    181    RECURSIVE SUBROUTINE Agrif_Update_vvl( ) 
     168   SUBROUTINE Agrif_Update_vvl( ) 
    182169      !!--------------------------------------------- 
    183170      !!   *** ROUTINE Agrif_Update_vvl *** 
     
    205192      CALL Agrif_ParentGrid_To_ChildGrid() 
    206193      ! 
    207       IF ( lk_agrif_doupd ) THEN ! Initialisation: do recursive update: 
    208          CALL Agrif_ChildGrid_To_ParentGrid() 
    209          CALL Agrif_Update_vvl() 
    210          CALL Agrif_ParentGrid_To_ChildGrid() 
    211       ENDIF 
    212       ! 
    213194#endif 
    214195      ! 
     
    232213      fse3u_a(:,:,:) = fse3u_n(:,:,:) 
    233214      fse3v_a(:,:,:) = fse3v_n(:,:,:) 
     215!      ua(:,:,:) = fse3u_b(:,:,:) 
     216!      va(:,:,:) = fse3v_b(:,:,:) 
    234217      hu_a(:,:) = hu(:,:) 
    235218      hv_a(:,:) = hv(:,:) 
     
    290273      !! 
    291274      INTEGER :: ji,jj,jk,jn 
     275      REAL(wp) :: ztb, ztnu, ztno 
    292276      !!--------------------------------------------- 
    293277      ! 
     
    320304                     DO ji=i1,i2 
    321305                        IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN 
    322                            tsb(ji,jj,jk,jn) = ( tsb(ji,jj,jk,jn)*fse3t_b(ji,jj,jk)   & ! jc: should be fse3t_b prior update 
    323                                 & + atfp * ( tabres(ji,jj,jk,jn)                     & 
    324                                 &             - tsn(ji,jj,jk,jn)*fse3t_a(ji,jj,jk) ) &  
    325                                 & * tmask(ji,jj,jk) ) / fse3t_b(ji,jj,jk) 
     306                           ztb  = tsb(ji,jj,jk,jn) * fse3t_b(ji,jj,jk) ! fse3t_b prior update should be used 
     307                           ztnu = tabres(ji,jj,jk,jn) 
     308                           ztno = tsn(ji,jj,jk,jn) * fse3t_a(ji,jj,jk) 
     309                           tsb(ji,jj,jk,jn) = ( ztb + atfp * ( ztnu - ztno) )  &  
     310                                     &        * tmask(ji,jj,jk) / fse3t_b(ji,jj,jk) 
    326311                        ENDIF 
    327312                     ENDDO 
     
    330315            ENDDO 
    331316         ENDIF 
     317 
    332318         DO jn = n1,n2 
    333319            DO jk=k1,k2 
     
    341327            END DO 
    342328         END DO 
     329         ! 
     330         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 
     331            tsb(i1:i2,j1:j2,k1:k2,n1:n2)  = tsn(i1:i2,j1:j2,k1:k2,n1:n2) 
     332         ENDIF 
     333         ! 
    343334      ENDIF 
    344335      !  
    345336   END SUBROUTINE updateTS 
    346337 
    347    SUBROUTINE updateu( tabres, i1, i2, j1, j2, k1, k2, before ) 
     338   SUBROUTINE updateu( tabres, i1, i2, j1, j2, k1, k2, before, nb, ndir ) 
    348339      !!--------------------------------------------- 
    349340      !!           *** ROUTINE updateu *** 
     
    354345      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres 
    355346      LOGICAL, INTENT(in) :: before 
     347      INTEGER, INTENT(in) :: nb , ndir 
    356348      !!  
     349      LOGICAL western_side, eastern_side 
    357350      INTEGER :: ji, jj, jk 
    358351      REAL(wp) :: zrhoy 
     352      REAL(wp) :: zub, zunu, zuno 
    359353      !!--------------------------------------------- 
    360354      !  
     
    371365         tabres = zrhoy * tabres 
    372366      ELSE 
     367         western_side = (nb == 1).AND.(ndir == 1) 
     368         eastern_side = (nb == 1).AND.(ndir == 2) 
    373369         DO jk=k1,k2 
    374370            DO jj=j1,j2 
     
    377373                  ! 
    378374                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 
    379                      ub(ji,jj,jk) = (       ub(ji,jj,jk)*fse3u_b(ji,jj,jk)   & ! jc: should be fse3u_b prior update 
    380                            & + atfp * ( tabres(ji,jj,jk)                     & 
    381                            &              - un(ji,jj,jk)*fse3u_a(ji,jj,jk) ) &  
    382                            & * umask(ji,jj,jk) ) / fse3u_b(ji,jj,jk) 
     375                     zub  = ub(ji,jj,jk) * fse3u_b(ji,jj,jk) 
     376                     zuno = un(ji,jj,jk) * fse3u_a(ji,jj,jk) 
     377                     zunu = tabres(ji,jj,jk) 
     378                     ub(ji,jj,jk) = ( zub + atfp * ( zunu - zuno) ) &       
     379                                    & * umask(ji,jj,jk) / fse3u_b(ji,jj,jk) 
    383380                  ENDIF 
    384381                  ! 
     
    387384            END DO 
    388385         END DO 
     386         ! 
     387!         IF (western_side) THEN 
     388!            DO jk=k1,k2 
     389!               DO jj=j1,j2 
     390!                 un(i1-1,jj,jk) = un(i1-1,jj,jk) * fse3u_a(i1-1,jj,jk) / fse3u_n(i1-1,jj,jk) 
     391!               END DO 
     392!            ENDDO 
     393!         ENDIF 
     394!         IF (eastern_side) THEN 
     395!            DO jk=k1,k2 
     396!               DO jj=j1,j2 
     397!                 un(i2+1,jj,jk) = un(i2+1,jj,jk) * fse3u_a(i2+1,jj,jk) / fse3u_n(i2+1,jj,jk) 
     398!               END DO 
     399!            ENDDO 
     400!         ENDIF 
     401         ! 
     402         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 
     403            ub(i1:i2,j1:j2,k1:k2)  = un(i1:i2,j1:j2,k1:k2) 
     404         ENDIF 
     405         ! 
    389406      ENDIF 
    390407      !  
    391408   END SUBROUTINE updateu 
    392409 
    393    SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, before ) 
     410   SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, before, nb, ndir ) 
    394411      !!--------------------------------------------- 
    395412      !!           *** ROUTINE updatev *** 
     
    401418      REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: tabres 
    402419      LOGICAL :: before 
    403       !! 
     420      INTEGER, INTENT(in) :: nb , ndir 
     421      !! 
     422      LOGICAL :: northern_side, southern_side 
    404423      REAL(wp) :: zrhox 
     424      REAL(wp) :: zvb, zvnu, zvno 
    405425      !!---------------------------------------------       
    406426      ! 
     
    417437         tabres = zrhox * tabres 
    418438      ELSE 
     439         southern_side = (nb == 2).AND.(ndir == 1) 
     440         northern_side = (nb == 2).AND.(ndir == 2) 
    419441         DO jk=k1,k2 
    420442            DO jj=j1,j2 
     
    423445                  ! 
    424446                  IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part 
    425                      vb(ji,jj,jk) = (       vb(ji,jj,jk)*fse3v_b(ji,jj,jk)   & ! jc: should be fse3v_b prior update 
    426                            & + atfp * ( tabres(ji,jj,jk)                     &  
    427                            &              - vn(ji,jj,jk)*fse3v_a(ji,jj,jk) ) &  
    428                            & * vmask(ji,jj,jk) ) / fse3v_b(ji,jj,jk) 
     447                     zvb  = vb(ji,jj,jk) * fse3v_b(ji,jj,jk) 
     448                     zvno = vn(ji,jj,jk) * fse3v_a(ji,jj,jk) 
     449                     zvnu = tabres(ji,jj,jk) 
     450                     vb(ji,jj,jk) = ( zvb + atfp * ( zvnu - zvno) ) &       
     451                                    & * vmask(ji,jj,jk) / fse3v_b(ji,jj,jk) 
    429452                  ENDIF 
    430453                  ! 
     
    433456            END DO 
    434457         END DO 
     458         ! 
     459!         IF (southern_side) THEN 
     460!            DO jk=k1,k2 
     461!               DO ji=i1,i2 
     462!                 vn(ji,j1-1,jk) = vn(ji,j1-1,jk) * fse3v_a(ji,j1-1,jk) / fse3v_n(ji,j1-1,jk) 
     463!               END DO 
     464!            ENDDO 
     465!         ENDIF 
     466!         IF (northern_side) THEN 
     467!            DO jk=k1,k2 
     468!               DO ji=i1,i2 
     469!                 vn(ji,j2+1,jk) = vn(ji,j2+1,jk) * fse3v_a(ji,j2+1,jk) / fse3v_n(ji,j2+1,jk) 
     470!               END DO 
     471!            ENDDO 
     472!         ENDIF 
     473         ! 
     474         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 
     475            vb(i1:i2,j1:j2,k1:k2)  = vn(i1:i2,j1:j2,k1:k2) 
     476         ENDIF 
     477         ! 
    435478      ENDIF 
    436479      !  
    437480   END SUBROUTINE updatev 
    438481 
    439    SUBROUTINE updateu2d( tabres, i1, i2, j1, j2, before ) 
     482   SUBROUTINE updateu2d( tabres, i1, i2, j1, j2, before, nb, ndir ) 
    440483      !!--------------------------------------------- 
    441484      !!          *** ROUTINE updateu2d *** 
     
    446489      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres 
    447490      LOGICAL, INTENT(in) :: before 
    448       !!  
     491      INTEGER, INTENT(in) :: nb , ndir 
     492      !! 
     493      LOGICAL western_side, eastern_side  
    449494      INTEGER :: ji, jj, jk 
    450495      REAL(wp) :: zrhoy 
     
    461506         tabres = zrhoy * tabres 
    462507      ELSE 
     508         western_side = (nb == 1).AND.(ndir == 1) 
     509         eastern_side = (nb == 1).AND.(ndir == 2) 
    463510         DO jj=j1,j2 
    464511            DO ji=i1,i2 
     
    500547            END DO 
    501548         END DO 
     549!         IF (western_side) THEN 
     550!            DO jj=j1,j2 
     551!              un_b(i1-1,jj) = un_b(i1-1,jj) * hu_a(i1-1,jj) * hur(i1-1,jj) 
     552!            END DO 
     553!         ENDIF 
     554!         IF (eastern_side) THEN 
     555!            DO jj=j1,j2 
     556!              un_b(i2+1,jj) = un_b(i2+1,jj) * hu_a(i2+1,jj) * hur(i2+1,jj) 
     557!            ENDDO 
     558!         ENDIF 
     559         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 
     560            ub_b(i1:i2,j1:j2)  = un_b(i1:i2,j1:j2) 
     561         ENDIF 
    502562      ENDIF 
    503563      ! 
    504564   END SUBROUTINE updateu2d 
    505565 
    506    SUBROUTINE updatev2d( tabres, i1, i2, j1, j2, before ) 
     566   SUBROUTINE updatev2d( tabres, i1, i2, j1, j2, before, nb, ndir ) 
    507567      !!--------------------------------------------- 
    508568      !!          *** ROUTINE updatev2d *** 
     
    511571      REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres 
    512572      LOGICAL, INTENT(in) :: before 
     573      INTEGER, INTENT(in) :: nb , ndir 
    513574      !!  
     575      LOGICAL :: northern_side, southern_side 
    514576      INTEGER :: ji, jj, jk 
    515577      REAL(wp) :: zrhox 
     
    526588         tabres = zrhox * tabres 
    527589      ELSE 
     590         southern_side = (nb == 2).AND.(ndir == 1) 
     591         northern_side = (nb == 2).AND.(ndir == 2) 
    528592         DO jj=j1,j2 
    529593            DO ji=i1,i2 
     
    565629            END DO 
    566630         END DO 
     631         ! 
     632!         IF (southern_side) THEN 
     633!            DO ji=i1,i2 
     634!              vn_b(ji,j1-1) = vn_b(ji,j1-1) * hv_a(ji,j1-1) * hvr(ji,j1-1) 
     635!            END DO 
     636!         ENDIF 
     637!         IF (northern_side) THEN 
     638!            DO ji=i1,i2 
     639!              vn_b(ji,j2+1) = vn_b(ji,j2+1) * hv_a(ji,j2+1) * hvr(ji,j2+1) 
     640!            END DO 
     641!         ENDIF 
     642         ! 
     643         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 
     644            vb_b(i1:i2,j1:j2)  = vn_b(i1:i2,j1:j2) 
     645         ENDIF 
     646         ! 
    567647      ENDIF 
    568648      !  
     
    606686            END DO 
    607687         END DO 
     688         ! 
     689         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 
     690            sshb(i1:i2,j1:j2)  = sshn(i1:i2,j1:j2) 
     691         ENDIF 
     692         ! 
    608693      ENDIF 
    609694      ! 
     
    619704      !! 
    620705      INTEGER :: ji, jj 
    621       REAL(wp) :: zrhoy 
     706      REAL(wp) :: zrhoy, za1 
    622707      !!--------------------------------------------- 
    623708      ! 
     
    631716         tabres = zrhoy * tabres 
    632717      ELSE 
     718         za1 = 1._wp / REAL(Agrif_rhot(), wp) 
     719         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) / e2u(i1:i2,j1:j2) 
    633720         DO jj=j1,j2 
    634             DO ji=i1,i2 
    635                ub2_b(ji,jj) = tabres(ji,jj) / e2u(ji,jj) 
     721            DO ji=i1,i2   
     722               ub2_i_b(ji,jj) = ub2_i_b(ji,jj) &  
     723                & + za1 * (tabres(ji,jj) - ub2_b(ji,jj)) 
     724               ub2_b(ji,jj) = tabres(ji,jj) 
    636725            END DO 
    637726         END DO 
     
    649738      !! 
    650739      INTEGER :: ji, jj 
    651       REAL(wp) :: zrhox 
     740      REAL(wp) :: zrhox, za1 
    652741      !!--------------------------------------------- 
    653742      ! 
     
    661750         tabres = zrhox * tabres 
    662751      ELSE 
     752         za1 = 1._wp / REAL(Agrif_rhot(), wp) 
     753         tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) / e1v(i1:i2,j1:j2) 
    663754         DO jj=j1,j2 
    664755            DO ji=i1,i2 
    665                vb2_b(ji,jj) = tabres(ji,jj) / e1v(ji,jj) 
     756               vb2_i_b(ji,jj) = vb2_i_b(ji,jj) &  
     757                & + za1 * (tabres(ji,jj) - vb2_b(ji,jj)) 
     758               vb2_b(ji,jj) = tabres(ji,jj) 
    666759            END DO 
    667760         END DO 
     
    800893         ptab(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2) * e3t_0(i1:i2,j1:j2,k1:k2) 
    801894!< jc tmp: 
     895 
     896         ! Save "old" scale factor (prior update) for subsequent asselin correction 
     897         ! of prognostic variables (needed to update initial state only) 
     898         fse3t_a(i1:i2,j1:j2,k1:k2) = fse3t_n(i1:i2,j1:j2,k1:k2) 
     899!         hdivb(i1:i2,j1:j2,k1:k2)   = fse3t_b(i1:i2,j1:j2,k1:k2) 
    802900 
    803901#if ! defined key_dynspg_ts 
     
    839937         ! ---------------------------- 
    840938         ! 
    841          ! Save "old" scale factor (prior update) for subsequent asselin correction 
    842          ! of prognostic variables (needed to update initial state only) 
    843          fse3t_a(i1:i2,j1:j2,k1:k2) = fse3t_n(i1:i2,j1:j2,k1:k2) 
    844          ! 
    845939         ! Update vertical scale factor at T-points: 
    846940         fse3t_n(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2) 
     
    872966         END DO 
    873967         ! 
     968         IF  ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN 
     969            fse3t_b (i1:i2,j1:j2,1:jpk)  = fse3t_n (i1:i2,j1:j2,1:jpk) 
     970            fse3w_b (i1:i2,j1:j2,1:jpk)  = fse3w_n (i1:i2,j1:j2,1:jpk) 
     971            fsdepw_b(i1:i2,j1:j2,1:jpk)  = fsdepw_n(i1:i2,j1:j2,1:jpk) 
     972            fsdept_b(i1:i2,j1:j2,1:jpk)  = fsdept_n(i1:i2,j1:j2,1:jpk) 
     973         ENDIF 
     974         ! 
    874975      ENDIF 
    875976      ! 
Note: See TracChangeset for help on using the changeset viewer.