#define TWO_WAY /* TWO WAY NESTING */ #undef DECAL_FEEDBACK /* SEPARATION of INTERFACES*/ MODULE agrif_opa_update #if defined key_agrif && ! defined key_offline USE par_oce USE oce USE dom_oce USE agrif_oce USE in_out_manager ! I/O manager USE lib_mpp USE wrk_nemo USE dynspg_oce USE zdf_oce ! vertical physics: ocean variables USE domvvl ! Need interpolation routines IMPLICIT NONE PRIVATE PUBLIC Agrif_Update_Tra, Agrif_Update_Dyn, Update_Scales, Agrif_Update_vvl # if defined key_zdftke PUBLIC Agrif_Update_Tke # endif !!---------------------------------------------------------------------- !! NEMO/NST 3.6 , NEMO Consortium (2010) !! $Id$ !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) !!---------------------------------------------------------------------- CONTAINS SUBROUTINE Agrif_Update_Tra( ) !!--------------------------------------------- !! *** ROUTINE Agrif_Update_Tra *** !!--------------------------------------------- ! IF (Agrif_Root()) RETURN ! #if defined TWO_WAY IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update tracers from grid Number',Agrif_Fixed(), 'nbcline', nbcline Agrif_UseSpecialValueInUpdate = .TRUE. Agrif_SpecialValueFineGrid = 0. ! IF (MOD(nbcline,nbclineupdate) == 0) THEN # if ! defined DECAL_FEEDBACK CALL Agrif_Update_Variable(tsn_id, procname=updateTS) # else CALL Agrif_Update_Variable(tsn_id, locupdate=(/1,0/),procname=updateTS) # endif ELSE # if ! defined DECAL_FEEDBACK CALL Agrif_Update_Variable(tsn_id,locupdate=(/0,2/), procname=updateTS) # else CALL Agrif_Update_Variable(tsn_id,locupdate=(/1,2/), procname=updateTS) # endif ENDIF ! Agrif_UseSpecialValueInUpdate = .FALSE. ! #endif ! END SUBROUTINE Agrif_Update_Tra SUBROUTINE Agrif_Update_Dyn( ) !!--------------------------------------------- !! *** ROUTINE Agrif_Update_Dyn *** !!--------------------------------------------- ! IF (Agrif_Root()) RETURN ! #if defined TWO_WAY IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update momentum from grid Number',Agrif_Fixed(), 'nbcline', nbcline Agrif_UseSpecialValueInUpdate = .FALSE. Agrif_SpecialValueFineGrid = 0. ! IF (mod(nbcline,nbclineupdate) == 0) THEN # if ! defined DECAL_FEEDBACK CALL Agrif_Update_Variable(un_update_id,procname = updateU) CALL Agrif_Update_Variable(vn_update_id,procname = updateV) # else CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU) CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV) # endif ELSE # if ! defined DECAL_FEEDBACK CALL Agrif_Update_Variable(un_update_id,locupdate=(/0,1/),procname = updateU) CALL Agrif_Update_Variable(vn_update_id,locupdate=(/0,1/),procname = updateV) # else CALL Agrif_Update_Variable(un_update_id,locupdate1=(/0,1/),locupdate2=(/1,1/),procname = updateU) CALL Agrif_Update_Variable(vn_update_id,locupdate1=(/1,1/),locupdate2=(/0,1/),procname = updateV) # endif ENDIF # if ! defined DECAL_FEEDBACK CALL Agrif_Update_Variable(e1u_id,procname = updateU2d) CALL Agrif_Update_Variable(e2v_id,procname = updateV2d) # else CALL Agrif_Update_Variable(e1u_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateU2d) CALL Agrif_Update_Variable(e2v_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updateV2d) # endif # if defined key_dynspg_ts IF (ln_bt_fw) THEN ! Update time integrated transports IF (mod(nbcline,nbclineupdate) == 0) THEN # if ! defined DECAL_FEEDBACK CALL Agrif_Update_Variable(ub2b_update_id,procname = updateub2b) CALL Agrif_Update_Variable(vb2b_update_id,procname = updatevb2b) # else CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0,-1/),locupdate2=(/1,-2/),procname = updateub2b) CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1,-2/),locupdate2=(/0,-1/),procname = updatevb2b) # endif ELSE # if ! defined DECAL_FEEDBACK CALL Agrif_Update_Variable(ub2b_update_id,locupdate=(/0,1/),procname = updateub2b) CALL Agrif_Update_Variable(vb2b_update_id,locupdate=(/0,1/),procname = updatevb2b) # else CALL Agrif_Update_Variable(ub2b_update_id,locupdate1=(/0,1/),locupdate2=(/1,1/),procname = updateub2b) CALL Agrif_Update_Variable(vb2b_update_id,locupdate1=(/1,1/),locupdate2=(/0,1/),procname = updatevb2b) # endif ENDIF END IF # endif ! nbcline = nbcline + 1 ! Agrif_UseSpecialValueInUpdate = .TRUE. Agrif_SpecialValueFineGrid = 0. # if ! defined DECAL_FEEDBACK CALL Agrif_Update_Variable(sshn_id,procname = updateSSH) # else CALL Agrif_Update_Variable(sshn_id,locupdate=(/1,0/),procname = updateSSH) # endif Agrif_UseSpecialValueInUpdate = .FALSE. ! #endif ! END SUBROUTINE Agrif_Update_Dyn # if defined key_zdftke SUBROUTINE Agrif_Update_Tke( kt ) !!--------------------------------------------- !! *** ROUTINE Agrif_Update_Tke *** !!--------------------------------------------- !! INTEGER, INTENT(in) :: kt ! IF( ( Agrif_NbStepint() .NE. 0 ) .AND. (kt /= 0) ) RETURN # if defined TWO_WAY Agrif_UseSpecialValueInUpdate = .TRUE. Agrif_SpecialValueFineGrid = 0. CALL Agrif_Update_Variable( en_id, locupdate=(/0,0/), procname=updateEN ) CALL Agrif_Update_Variable(avt_id, locupdate=(/0,0/), procname=updateAVT ) CALL Agrif_Update_Variable(avm_id, locupdate=(/0,0/), procname=updateAVM ) Agrif_UseSpecialValueInUpdate = .FALSE. # endif END SUBROUTINE Agrif_Update_Tke # endif /* key_zdftke */ SUBROUTINE Agrif_Update_vvl( ) !!--------------------------------------------- !! *** ROUTINE Agrif_Update_vvl *** !!--------------------------------------------- ! IF (Agrif_Root()) RETURN ! #if defined TWO_WAY ! IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Update e3 from grid Number',Agrif_Fixed(), 'Step', Agrif_Nb_Step() ! Agrif_UseSpecialValueInUpdate = .FALSE. Agrif_SpecialValueFineGrid = 0. ! # if ! defined DECAL_FEEDBACK CALL Agrif_Update_Variable(e3t_id, procname=updatee3t) # else CALL Agrif_Update_Variable(e3t_id, locupdate=(/1,0/), procname=updatee3t) # endif ! Agrif_UseSpecialValueInUpdate = .FALSE. ! CALL Agrif_ChildGrid_To_ParentGrid() CALL dom_vvl_update_UVF CALL Agrif_ParentGrid_To_ChildGrid() ! #endif ! END SUBROUTINE Agrif_Update_vvl SUBROUTINE dom_vvl_update_UVF !!--------------------------------------------- !! *** ROUTINE dom_vvl_update_UVF *** !!--------------------------------------------- # include "domzgr_substitute.h90" !! INTEGER :: jk REAL(wp):: zcoef !!--------------------------------------------- IF (lwp.AND.lk_agrif_debug) Write(*,*) 'Finalize e3 on grid Number', & & Agrif_Fixed(), 'Step', Agrif_Nb_Step() ! Save "old" scale factor (prior update) for subsequent asselin correction ! of prognostic variables (needed to update initial state only) fse3u_a(:,:,:) = fse3u_n(:,:,:) fse3v_a(:,:,:) = fse3v_n(:,:,:) ! ua(:,:,:) = fse3u_b(:,:,:) ! va(:,:,:) = fse3v_b(:,:,:) hu_a(:,:) = hu(:,:) hv_a(:,:) = hv(:,:) ! Vertical scale factor interpolations ! ------------------------------------ ! CALL dom_vvl_interpol( fse3t_n(:,:,:), fse3u_n(:,:,:), 'U' ) CALL dom_vvl_interpol( fse3t_n(:,:,:), fse3v_n(:,:,:), 'V' ) CALL dom_vvl_interpol( fse3u_n(:,:,:), fse3f_n(:,:,:), 'F' ) CALL dom_vvl_interpol( fse3u_n(:,:,:), fse3uw_n(:,:,:), 'UW' ) CALL dom_vvl_interpol( fse3v_n(:,:,:), fse3vw_n(:,:,:), 'VW' ) ! Update total depths: hu(:,:) = 0._wp ! Ocean depth at U-points hv(:,:) = 0._wp ! Ocean depth at V-points DO jk = 1, jpkm1 hu(:,:) = hu(:,:) + fse3u_n(:,:,jk) * umask(:,:,jk) hv(:,:) = hv(:,:) + fse3v_n(:,:,jk) * vmask(:,:,jk) END DO ! ! Inverse of the local depth hur(:,:) = 1._wp / ( hu(:,:) + 1._wp - umask_i(:,:) ) * umask_i(:,:) hvr(:,:) = 1._wp / ( hv(:,:) + 1._wp - vmask_i(:,:) ) * vmask_i(:,:) #if ! defined key_dynspg_ts IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN #else IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN #endif CALL dom_vvl_interpol( fse3t_b(:,:,:), fse3u_b(:,:,:), 'U' ) CALL dom_vvl_interpol( fse3t_b(:,:,:), fse3v_b(:,:,:), 'V' ) CALL dom_vvl_interpol( fse3u_b(:,:,:), fse3uw_b(:,:,:), 'UW' ) CALL dom_vvl_interpol( fse3v_b(:,:,:), fse3vw_b(:,:,:), 'VW' ) hu_b(:,:) = 0._wp ! Ocean depth at U-points hv_b(:,:) = 0._wp ! Ocean depth at V-points DO jk = 1, jpkm1 hu_b(:,:) = hu_b(:,:) + fse3u_b(:,:,jk) * umask(:,:,jk) hv_b(:,:) = hv_b(:,:) + fse3v_b(:,:,jk) * vmask(:,:,jk) END DO ! ! Inverse of the local depth hur_b(:,:) = 1._wp / ( hu_b(:,:) + 1._wp - umask_i(:,:) ) * umask_i(:,:) hvr_b(:,:) = 1._wp / ( hv_b(:,:) + 1._wp - vmask_i(:,:) ) * vmask_i(:,:) ENDIF ! END SUBROUTINE dom_vvl_update_UVF SUBROUTINE updateTS( tabres, i1, i2, j1, j2, k1, k2, n1, n2, before ) !!--------------------------------------------- !! *** ROUTINE updateT *** !!--------------------------------------------- INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2 REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres LOGICAL, INTENT(in) :: before !! INTEGER :: ji,jj,jk,jn REAL(wp) :: ztb, ztnu, ztno !!--------------------------------------------- ! ! IF (before) THEN DO jn = n1,n2 DO jk=k1,k2 DO jj=j1,j2 DO ji=i1,i2 !> jc tmp tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn) * fse3t_n(ji,jj,jk) / e3t_0(ji,jj,jk) ! tabres(ji,jj,jk,jn) = tsn(ji,jj,jk,jn) * fse3t_n(ji,jj,jk) !< jc tmp END DO END DO END DO END DO ELSE !> jc tmp DO jn = n1,n2 tabres(i1:i2,j1:j2,k1:k2,jn) = tabres(i1:i2,j1:j2,k1:k2,jn) * e3t_0(i1:i2,j1:j2,k1:k2) ENDDO !< jc tmp IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part DO jn = n1,n2 DO jk=k1,k2 DO jj=j1,j2 DO ji=i1,i2 IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN ztb = tsb(ji,jj,jk,jn) * fse3t_b(ji,jj,jk) ! fse3t_b prior update should be used ztnu = tabres(ji,jj,jk,jn) ztno = tsn(ji,jj,jk,jn) * fse3t_a(ji,jj,jk) tsb(ji,jj,jk,jn) = ( ztb + atfp * ( ztnu - ztno) ) & & * tmask(ji,jj,jk) / fse3t_b(ji,jj,jk) ENDIF ENDDO ENDDO ENDDO ENDDO ENDIF DO jn = n1,n2 DO jk=k1,k2 DO jj=j1,j2 DO ji=i1,i2 IF( tabres(ji,jj,jk,jn) .NE. 0. ) THEN tsn(ji,jj,jk,jn) = tabres(ji,jj,jk,jn) / fse3t_n(ji,jj,jk) END IF END DO END DO END DO END DO ! IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN tsb(i1:i2,j1:j2,k1:k2,n1:n2) = tsn(i1:i2,j1:j2,k1:k2,n1:n2) ENDIF ! ENDIF ! END SUBROUTINE updateTS SUBROUTINE updateu( tabres, i1, i2, j1, j2, k1, k2, before, nb, ndir ) !!--------------------------------------------- !! *** ROUTINE updateu *** !!--------------------------------------------- # include "domzgr_substitute.h90" !! INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: tabres LOGICAL, INTENT(in) :: before INTEGER, INTENT(in) :: nb , ndir !! LOGICAL western_side, eastern_side INTEGER :: ji, jj, jk REAL(wp) :: zrhoy REAL(wp) :: zub, zunu, zuno !!--------------------------------------------- ! IF (before) THEN zrhoy = Agrif_Rhoy() DO jk=k1,k2 DO jj=j1,j2 DO ji=i1,i2 tabres(ji,jj,jk) = e2u(ji,jj) * un(ji,jj,jk) tabres(ji,jj,jk) = tabres(ji,jj,jk) * fse3u_n(ji,jj,jk) END DO END DO END DO tabres = zrhoy * tabres ELSE western_side = (nb == 1).AND.(ndir == 1) eastern_side = (nb == 1).AND.(ndir == 2) DO jk=k1,k2 DO jj=j1,j2 DO ji=i1,i2 tabres(ji,jj,jk) = tabres(ji,jj,jk) / e2u(ji,jj) ! IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part zub = ub(ji,jj,jk) * fse3u_b(ji,jj,jk) zuno = un(ji,jj,jk) * fse3u_a(ji,jj,jk) zunu = tabres(ji,jj,jk) ub(ji,jj,jk) = ( zub + atfp * ( zunu - zuno) ) & & * umask(ji,jj,jk) / fse3u_b(ji,jj,jk) ENDIF ! un(ji,jj,jk) = tabres(ji,jj,jk) * umask(ji,jj,jk) / fse3u_n(ji,jj,jk) END DO END DO END DO ! ! IF (western_side) THEN ! DO jk=k1,k2 ! DO jj=j1,j2 ! un(i1-1,jj,jk) = un(i1-1,jj,jk) * fse3u_a(i1-1,jj,jk) / fse3u_n(i1-1,jj,jk) ! END DO ! ENDDO ! ENDIF ! IF (eastern_side) THEN ! DO jk=k1,k2 ! DO jj=j1,j2 ! un(i2+1,jj,jk) = un(i2+1,jj,jk) * fse3u_a(i2+1,jj,jk) / fse3u_n(i2+1,jj,jk) ! END DO ! ENDDO ! ENDIF ! IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN ub(i1:i2,j1:j2,k1:k2) = un(i1:i2,j1:j2,k1:k2) ENDIF ! ENDIF ! END SUBROUTINE updateu SUBROUTINE updatev( tabres, i1, i2, j1, j2, k1, k2, before, nb, ndir ) !!--------------------------------------------- !! *** ROUTINE updatev *** !!--------------------------------------------- # include "domzgr_substitute.h90" !! INTEGER :: i1,i2,j1,j2,k1,k2 INTEGER :: ji,jj,jk REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2) :: tabres LOGICAL :: before INTEGER, INTENT(in) :: nb , ndir !! LOGICAL :: northern_side, southern_side REAL(wp) :: zrhox REAL(wp) :: zvb, zvnu, zvno !!--------------------------------------------- ! IF (before) THEN zrhox = Agrif_Rhox() DO jk=k1,k2 DO jj=j1,j2 DO ji=i1,i2 tabres(ji,jj,jk) = e1v(ji,jj) * vn(ji,jj,jk) tabres(ji,jj,jk) = tabres(ji,jj,jk) * fse3v_n(ji,jj,jk) END DO END DO END DO tabres = zrhox * tabres ELSE southern_side = (nb == 2).AND.(ndir == 1) northern_side = (nb == 2).AND.(ndir == 2) DO jk=k1,k2 DO jj=j1,j2 DO ji=i1,i2 tabres(ji,jj,jk) = tabres(ji,jj,jk) / e1v(ji,jj) ! IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN ! Add asselin part zvb = vb(ji,jj,jk) * fse3v_b(ji,jj,jk) zvno = vn(ji,jj,jk) * fse3v_a(ji,jj,jk) zvnu = tabres(ji,jj,jk) vb(ji,jj,jk) = ( zvb + atfp * ( zvnu - zvno) ) & & * vmask(ji,jj,jk) / fse3v_b(ji,jj,jk) ENDIF ! vn(ji,jj,jk) = tabres(ji,jj,jk) * vmask(ji,jj,jk) / fse3v_n(ji,jj,jk) END DO END DO END DO ! ! IF (southern_side) THEN ! DO jk=k1,k2 ! DO ji=i1,i2 ! vn(ji,j1-1,jk) = vn(ji,j1-1,jk) * fse3v_a(ji,j1-1,jk) / fse3v_n(ji,j1-1,jk) ! END DO ! ENDDO ! ENDIF ! IF (northern_side) THEN ! DO jk=k1,k2 ! DO ji=i1,i2 ! vn(ji,j2+1,jk) = vn(ji,j2+1,jk) * fse3v_a(ji,j2+1,jk) / fse3v_n(ji,j2+1,jk) ! END DO ! ENDDO ! ENDIF ! IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN vb(i1:i2,j1:j2,k1:k2) = vn(i1:i2,j1:j2,k1:k2) ENDIF ! ENDIF ! END SUBROUTINE updatev SUBROUTINE updateu2d( tabres, i1, i2, j1, j2, before, nb, ndir ) !!--------------------------------------------- !! *** ROUTINE updateu2d *** !!--------------------------------------------- # include "domzgr_substitute.h90" !! INTEGER, INTENT(in) :: i1, i2, j1, j2 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres LOGICAL, INTENT(in) :: before INTEGER, INTENT(in) :: nb , ndir !! LOGICAL western_side, eastern_side INTEGER :: ji, jj, jk REAL(wp) :: zrhoy REAL(wp) :: zcorr !!--------------------------------------------- ! IF (before) THEN zrhoy = Agrif_Rhoy() DO jj=j1,j2 DO ji=i1,i2 tabres(ji,jj) = un_b(ji,jj) * hu(ji,jj) * e2u(ji,jj) END DO END DO tabres = zrhoy * tabres ELSE western_side = (nb == 1).AND.(ndir == 1) eastern_side = (nb == 1).AND.(ndir == 2) DO jj=j1,j2 DO ji=i1,i2 tabres(ji,jj) = tabres(ji,jj) / e2u(ji,jj) ! ! Update "now" 3d velocities: spgu(ji,jj) = 0.e0 DO jk=1,jpkm1 spgu(ji,jj) = spgu(ji,jj) + fse3u_n(ji,jj,jk) * un(ji,jj,jk) END DO ! zcorr = (tabres(ji,jj) - spgu(ji,jj)) * hur(ji,jj) DO jk=1,jpkm1 un(ji,jj,jk) = un(ji,jj,jk) + zcorr * umask(ji,jj,jk) END DO ! ! Update barotropic velocities: #if ! defined key_dynspg_ts IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN #else IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN #endif zcorr = (tabres(ji,jj) - un_b(ji,jj) * hu_a(ji,jj)) * hur_b(ji,jj) ub_b(ji,jj) = ub_b(ji,jj) + atfp * zcorr * umask(ji,jj,1) END IF un_b(ji,jj) = tabres(ji,jj) * hur(ji,jj) * umask(ji,jj,1) ! ! Correct "before" velocities to hold correct bt component: spgu(ji,jj) = 0.e0 DO jk=1,jpkm1 spgu(ji,jj) = spgu(ji,jj) + fse3u_b(ji,jj,jk) * ub(ji,jj,jk) END DO ! zcorr = ub_b(ji,jj) - spgu(ji,jj) * hur_b(ji,jj) DO jk=1,jpkm1 ub(ji,jj,jk) = ub(ji,jj,jk) + zcorr * umask(ji,jj,jk) END DO ! END DO END DO ! IF (western_side) THEN ! DO jj=j1,j2 ! un_b(i1-1,jj) = un_b(i1-1,jj) * hu_a(i1-1,jj) * hur(i1-1,jj) ! END DO ! ENDIF ! IF (eastern_side) THEN ! DO jj=j1,j2 ! un_b(i2+1,jj) = un_b(i2+1,jj) * hu_a(i2+1,jj) * hur(i2+1,jj) ! ENDDO ! ENDIF IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN ub_b(i1:i2,j1:j2) = un_b(i1:i2,j1:j2) ENDIF ENDIF ! END SUBROUTINE updateu2d SUBROUTINE updatev2d( tabres, i1, i2, j1, j2, before, nb, ndir ) !!--------------------------------------------- !! *** ROUTINE updatev2d *** !!--------------------------------------------- INTEGER, INTENT(in) :: i1, i2, j1, j2 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres LOGICAL, INTENT(in) :: before INTEGER, INTENT(in) :: nb , ndir !! LOGICAL :: northern_side, southern_side INTEGER :: ji, jj, jk REAL(wp) :: zrhox REAL(wp) :: zcorr !!--------------------------------------------- ! IF (before) THEN zrhox = Agrif_Rhox() DO jj=j1,j2 DO ji=i1,i2 tabres(ji,jj) = vn_b(ji,jj) * hv(ji,jj) * e1v(ji,jj) END DO END DO tabres = zrhox * tabres ELSE southern_side = (nb == 2).AND.(ndir == 1) northern_side = (nb == 2).AND.(ndir == 2) DO jj=j1,j2 DO ji=i1,i2 tabres(ji,jj) = tabres(ji,jj) / e1v(ji,jj) ! ! Update "now" 3d velocities: spgv(ji,jj) = 0.e0 DO jk=1,jpkm1 spgv(ji,jj) = spgv(ji,jj) + fse3v_n(ji,jj,jk) * vn(ji,jj,jk) END DO ! zcorr = (tabres(ji,jj) - spgv(ji,jj)) * hvr(ji,jj) DO jk=1,jpkm1 vn(ji,jj,jk) = vn(ji,jj,jk) + zcorr * vmask(ji,jj,jk) END DO ! ! Update barotropic velocities: #if ! defined key_dynspg_ts IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN #else IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN #endif zcorr = (tabres(ji,jj) - vn_b(ji,jj)*hv_a(ji,jj)) * hvr_b(ji,jj) vb_b(ji,jj) = vb_b(ji,jj) + atfp * zcorr * vmask(ji,jj,1) END IF vn_b(ji,jj) = tabres(ji,jj) * hvr(ji,jj) * vmask(ji,jj,1) ! ! Correct "before" velocities to hold correct bt component: spgv(ji,jj) = 0.e0 DO jk=1,jpkm1 spgv(ji,jj) = spgv(ji,jj) + fse3v_b(ji,jj,jk) * vb(ji,jj,jk) END DO ! zcorr = vb_b(ji,jj) - spgv(ji,jj) * hvr_b(ji,jj) DO jk=1,jpkm1 vb(ji,jj,jk) = vb(ji,jj,jk) + zcorr * vmask(ji,jj,jk) END DO ! END DO END DO ! ! IF (southern_side) THEN ! DO ji=i1,i2 ! vn_b(ji,j1-1) = vn_b(ji,j1-1) * hv_a(ji,j1-1) * hvr(ji,j1-1) ! END DO ! ENDIF ! IF (northern_side) THEN ! DO ji=i1,i2 ! vn_b(ji,j2+1) = vn_b(ji,j2+1) * hv_a(ji,j2+1) * hvr(ji,j2+1) ! END DO ! ENDIF ! IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN vb_b(i1:i2,j1:j2) = vn_b(i1:i2,j1:j2) ENDIF ! ENDIF ! END SUBROUTINE updatev2d SUBROUTINE updateSSH( tabres, i1, i2, j1, j2, before ) !!--------------------------------------------- !! *** ROUTINE updateSSH *** !!--------------------------------------------- INTEGER, INTENT(in) :: i1, i2, j1, j2 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres LOGICAL, INTENT(in) :: before !! INTEGER :: ji, jj !!--------------------------------------------- ! IF (before) THEN DO jj=j1,j2 DO ji=i1,i2 tabres(ji,jj) = sshn(ji,jj) END DO END DO ELSE #if ! defined key_dynspg_ts IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN #else IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN #endif DO jj=j1,j2 DO ji=i1,i2 sshb(ji,jj) = sshb(ji,jj) & & + atfp * ( tabres(ji,jj) - sshn(ji,jj) ) * tmask(ji,jj,1) END DO END DO ENDIF ! DO jj=j1,j2 DO ji=i1,i2 sshn(ji,jj) = tabres(ji,jj) * tmask(ji,jj,1) END DO END DO ! IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN sshb(i1:i2,j1:j2) = sshn(i1:i2,j1:j2) ENDIF ! ENDIF ! END SUBROUTINE updateSSH SUBROUTINE updateub2b( tabres, i1, i2, j1, j2, before ) !!--------------------------------------------- !! *** ROUTINE updateub2b *** !!--------------------------------------------- INTEGER, INTENT(in) :: i1, i2, j1, j2 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres LOGICAL, INTENT(in) :: before !! INTEGER :: ji, jj REAL(wp) :: zrhoy, za1 !!--------------------------------------------- ! IF (before) THEN zrhoy = Agrif_Rhoy() DO jj=j1,j2 DO ji=i1,i2 tabres(ji,jj) = ub2_i_b(ji,jj) * e2u(ji,jj) END DO END DO tabres = zrhoy * tabres ELSE za1 = 1._wp / REAL(Agrif_rhot(), wp) tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) / e2u(i1:i2,j1:j2) DO jj=j1,j2 DO ji=i1,i2 ub2_i_b(ji,jj) = ub2_i_b(ji,jj) & & + za1 * (tabres(ji,jj) - ub2_b(ji,jj)) ub2_b(ji,jj) = tabres(ji,jj) END DO END DO ENDIF ! END SUBROUTINE updateub2b SUBROUTINE updatevb2b( tabres, i1, i2, j1, j2, before ) !!--------------------------------------------- !! *** ROUTINE updatevb2b *** !!--------------------------------------------- INTEGER, INTENT(in) :: i1, i2, j1, j2 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres LOGICAL, INTENT(in) :: before !! INTEGER :: ji, jj REAL(wp) :: zrhox, za1 !!--------------------------------------------- ! IF (before) THEN zrhox = Agrif_Rhox() DO jj=j1,j2 DO ji=i1,i2 tabres(ji,jj) = vb2_i_b(ji,jj) * e1v(ji,jj) END DO END DO tabres = zrhox * tabres ELSE za1 = 1._wp / REAL(Agrif_rhot(), wp) tabres(i1:i2,j1:j2) = tabres(i1:i2,j1:j2) / e1v(i1:i2,j1:j2) DO jj=j1,j2 DO ji=i1,i2 vb2_i_b(ji,jj) = vb2_i_b(ji,jj) & & + za1 * (tabres(ji,jj) - vb2_b(ji,jj)) vb2_b(ji,jj) = tabres(ji,jj) END DO END DO ENDIF ! END SUBROUTINE updatevb2b SUBROUTINE update_scales( tabres, i1, i2, j1, j2, k1, k2, n1,n2, before ) ! currently not used !!--------------------------------------------- !! *** ROUTINE updateT *** !!--------------------------------------------- # include "domzgr_substitute.h90" INTEGER, INTENT(in) :: i1,i2,j1,j2,k1,k2,n1,n2 REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,n1:n2), INTENT(inout) :: tabres LOGICAL, iNTENT(in) :: before INTEGER :: ji,jj,jk REAL(wp) :: ztemp IF (before) THEN DO jk=k1,k2 DO jj=j1,j2 DO ji=i1,i2 tabres(ji,jj,jk,1) = e1t(ji,jj)*e2t(ji,jj)*tmask(ji,jj,jk) tabres(ji,jj,jk,2) = e1t(ji,jj)*tmask(ji,jj,jk) tabres(ji,jj,jk,3) = e2t(ji,jj)*tmask(ji,jj,jk) END DO END DO END DO tabres(:,:,:,1)=tabres(:,:,:,1)*Agrif_Rhox()*Agrif_Rhoy() tabres(:,:,:,2)=tabres(:,:,:,2)*Agrif_Rhox() tabres(:,:,:,3)=tabres(:,:,:,3)*Agrif_Rhoy() ELSE DO jk=k1,k2 DO jj=j1,j2 DO ji=i1,i2 IF( tabres(ji,jj,jk,1) .NE. 0. ) THEN print *,'VAL = ',ji,jj,jk,tabres(ji,jj,jk,1),e1t(ji,jj)*e2t(ji,jj)*tmask(ji,jj,jk) print *,'VAL2 = ',ji,jj,jk,tabres(ji,jj,jk,2),e1t(ji,jj)*tmask(ji,jj,jk) print *,'VAL3 = ',ji,jj,jk,tabres(ji,jj,jk,3),e2t(ji,jj)*tmask(ji,jj,jk) ztemp = sqrt(tabres(ji,jj,jk,1)/(tabres(ji,jj,jk,2)*tabres(ji,jj,jk,3))) print *,'CORR = ',ztemp-1. print *,'NEW VALS = ',tabres(ji,jj,jk,2)*ztemp,tabres(ji,jj,jk,3)*ztemp, & tabres(ji,jj,jk,2)*ztemp*tabres(ji,jj,jk,3)*ztemp e1t(ji,jj) = tabres(ji,jj,jk,2)*ztemp e2t(ji,jj) = tabres(ji,jj,jk,3)*ztemp END IF END DO END DO END DO ENDIF ! END SUBROUTINE update_scales # if defined key_zdftke SUBROUTINE updateEN( ptab, i1, i2, j1, j2, k1, k2, before ) !!--------------------------------------------- !! *** ROUTINE updateen *** !!--------------------------------------------- INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab LOGICAL, INTENT(in) :: before !!--------------------------------------------- ! IF (before) THEN ptab (i1:i2,j1:j2,k1:k2) = en(i1:i2,j1:j2,k1:k2) ELSE en(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) ENDIF ! END SUBROUTINE updateEN SUBROUTINE updateAVT( ptab, i1, i2, j1, j2, k1, k2, before ) !!--------------------------------------------- !! *** ROUTINE updateavt *** !!--------------------------------------------- INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab LOGICAL, INTENT(in) :: before !!--------------------------------------------- ! IF (before) THEN ptab (i1:i2,j1:j2,k1:k2) = avt_k(i1:i2,j1:j2,k1:k2) ELSE avt_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) ENDIF ! END SUBROUTINE updateAVT SUBROUTINE updateAVM( ptab, i1, i2, j1, j2, k1, k2, before ) !!--------------------------------------------- !! *** ROUTINE updateavm *** !!--------------------------------------------- INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab LOGICAL, INTENT(in) :: before !!--------------------------------------------- ! IF (before) THEN ptab (i1:i2,j1:j2,k1:k2) = avm_k(i1:i2,j1:j2,k1:k2) ELSE avm_k(i1:i2,j1:j2,k1:k2) = ptab (i1:i2,j1:j2,k1:k2) ENDIF ! END SUBROUTINE updateAVM # endif /* key_zdftke */ SUBROUTINE updatee3t( ptab, i1, i2, j1, j2, k1, k2, before ) !!--------------------------------------------- !! *** ROUTINE updatee3t *** !!--------------------------------------------- INTEGER, INTENT(in) :: i1, i2, j1, j2, k1, k2 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab LOGICAL, INTENT(in) :: before INTEGER :: ji,jj,jk REAL(wp) :: zcoef !!--------------------------------------------- ! IF (before) THEN !> jc tmp: ! ptab(i1:i2,j1:j2,k1:k2) = fse3t_n(i1:i2,j1:j2,k1:k2) ptab(i1:i2,j1:j2,k1:k2) = fse3t_n(i1:i2,j1:j2,k1:k2) / e3t_0(i1:i2,j1:j2,k1:k2) !< jc tmp: ELSE ! ! 1) Updates at before time step: ! ------------------------------- ! !> jc tmp: ptab(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2) * e3t_0(i1:i2,j1:j2,k1:k2) !< jc tmp: ! Save "old" scale factor (prior update) for subsequent asselin correction ! of prognostic variables (needed to update initial state only) fse3t_a(i1:i2,j1:j2,k1:k2) = fse3t_n(i1:i2,j1:j2,k1:k2) ! hdivb(i1:i2,j1:j2,k1:k2) = fse3t_b(i1:i2,j1:j2,k1:k2) #if ! defined key_dynspg_ts IF (.NOT.(lk_agrif_fstep.AND.(neuler==0))) THEN #else IF ((.NOT.(lk_agrif_fstep.AND.(neuler==0))).AND.(.NOT.ln_bt_fw)) THEN #endif DO jk = 1, jpkm1 DO jj=j1,j2 DO ji=i1,i2 fse3t_b(ji,jj,jk) = fse3t_b(ji,jj,jk) & & + atfp * ( ptab(ji,jj,jk) - fse3t_n(ji,jj,jk) ) END DO END DO END DO ! fse3w_b (i1:i2,j1:j2,1) = e3w_0(i1:i2,j1:j2,1) + fse3t_b(i1:i2,j1:j2,1) - e3t_0(i1:i2,j1:j2,1) fsdepw_b(i1:i2,j1:j2,1) = 0.0_wp fsdept_b(i1:i2,j1:j2,1) = 0.5_wp * fse3w_b(i1:i2,j1:j2,1) ! DO jk = 2, jpk DO jj = j1,j2 DO ji = i1,i2 zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) fse3w_b(ji,jj,jk) = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) * & & ( fse3t_b(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) ) & & + 0.5_wp * tmask(ji,jj,jk) * & & ( fse3t_b(ji,jj,jk ) - e3t_0(ji,jj,jk ) ) fsdepw_b(ji,jj,jk) = fsdepw_b(ji,jj,jk-1) + fse3t_b(ji,jj,jk-1) fsdept_b(ji,jj,jk) = zcoef * ( fsdepw_b(ji,jj,jk ) + 0.5 * fse3w_b(ji,jj,jk)) & & + (1-zcoef) * ( fsdept_b(ji,jj,jk-1) + fse3w_b(ji,jj,jk)) END DO END DO END DO ! ENDIF ! ! 2) Updates at now time step: ! ---------------------------- ! ! Update vertical scale factor at T-points: fse3t_n(i1:i2,j1:j2,k1:k2) = ptab(i1:i2,j1:j2,k1:k2) ! ! Update total depth: ht(i1:i2,j1:j2) = 0._wp DO jk = 1, jpkm1 ht(i1:i2,j1:j2) = ht(i1:i2,j1:j2) + fse3t_n(i1:i2,j1:j2,jk) * tmask(i1:i2,j1:j2,jk) END DO ! ! Update vertical scale factor at W-points and depths: fse3w_n (i1:i2,j1:j2,1) = e3w_0(i1:i2,j1:j2,1) + fse3t_n(i1:i2,j1:j2,1) - e3t_0(i1:i2,j1:j2,1) fsdept_n(i1:i2,j1:j2,1) = 0.5_wp * fse3w_n(i1:i2,j1:j2,1) fsdepw_n(i1:i2,j1:j2,1) = 0.0_wp fsde3w_n(i1:i2,j1:j2,1) = fsdept_n(i1:i2,j1:j2,1) - (ht(i1:i2,j1:j2)-ht_0(i1:i2,j1:j2)) ! Last term in the rhs is ssh ! DO jk = 2, jpk DO jj = j1,j2 DO ji = i1,i2 zcoef = (tmask(ji,jj,jk) - wmask(ji,jj,jk)) fse3w_n(ji,jj,jk) = e3w_0(ji,jj,jk) + ( 1.0_wp - 0.5_wp * tmask(ji,jj,jk) ) * ( fse3t_n(ji,jj,jk-1) - e3t_0(ji,jj,jk-1) ) & & + 0.5_wp * tmask(ji,jj,jk) * ( fse3t_n(ji,jj,jk ) - e3t_0(ji,jj,jk ) ) fsdepw_n(ji,jj,jk) = fsdepw_n(ji,jj,jk-1) + fse3t_n(ji,jj,jk-1) fsdept_n(ji,jj,jk) = zcoef * ( fsdepw_n(ji,jj,jk ) + 0.5 * fse3w_n(ji,jj,jk)) & & + (1-zcoef) * ( fsdept_n(ji,jj,jk-1) + fse3w_n(ji,jj,jk)) fsde3w_n(ji,jj,jk) = fsdept_n(ji,jj,jk) - (ht(ji,jj)-ht_0(ji,jj)) ! Last term in the rhs is ssh END DO END DO END DO ! IF ((neuler==0).AND.(Agrif_Nb_Step()==0) ) THEN fse3t_b (i1:i2,j1:j2,1:jpk) = fse3t_n (i1:i2,j1:j2,1:jpk) fse3w_b (i1:i2,j1:j2,1:jpk) = fse3w_n (i1:i2,j1:j2,1:jpk) fsdepw_b(i1:i2,j1:j2,1:jpk) = fsdepw_n(i1:i2,j1:j2,1:jpk) fsdept_b(i1:i2,j1:j2,1:jpk) = fsdept_n(i1:i2,j1:j2,1:jpk) ENDIF ! ENDIF ! END SUBROUTINE updatee3t #else CONTAINS SUBROUTINE agrif_opa_update_empty !!--------------------------------------------- !! *** ROUTINE agrif_opa_update_empty *** !!--------------------------------------------- WRITE(*,*) 'agrif_opa_update : You should not have seen this print! error?' END SUBROUTINE agrif_opa_update_empty #endif END MODULE agrif_opa_update