- Timestamp:
- 2019-03-26T09:50:57+01:00 (5 years ago)
- Location:
- NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src
- Files:
-
- 12 edited
-
OCE/DYN/dynadv_ubs.F90 (modified) (1 diff)
-
OCE/DYN/dynvor.F90 (modified) (1 diff)
-
OCE/TRA/traadv.F90 (modified) (6 diffs)
-
OCE/TRA/traadv_cen.F90 (modified) (11 diffs)
-
OCE/TRA/traadv_fct.F90 (modified) (21 diffs)
-
OCE/TRA/traadv_mus.F90 (modified) (11 diffs)
-
OCE/TRA/traadv_qck.F90 (modified) (20 diffs)
-
OCE/TRA/traadv_ubs.F90 (modified) (20 diffs)
-
OCE/oce.F90 (modified) (3 diffs)
-
OCE/step.F90 (modified) (2 diffs)
-
TOP/TRP/trcadv.F90 (modified) (3 diffs)
-
TOP/TRP/trctrp.F90 (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynadv_ubs.F90
r10789 r10802 234 234 ENDIF 235 235 ! ! Control print 236 IF(ln_ctl) CALL prt_ctl( tab3d_1=ua, clinfo1=' u u(:,:,:,ktlev1)s2 adv - Ua: ', mask1=umask, &236 IF(ln_ctl) CALL prt_ctl( tab3d_1=ua, clinfo1=' ubs2 adv - Ua: ', mask1=umask, & 237 237 & tab3d_2=va, clinfo2= ' Va: ', mask2=vmask, clinfo3='dyn' ) 238 238 ! -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynvor.F90
r10789 r10802 878 878 IF( ln_dynvor_msk ) CALL ctl_stop( 'dyn_vor_init: masked vorticity is not currently not available') 879 879 880 !!gm this should be removed when choosing a u u(:,:,:,ktlev)ique strategy for fmask at the coast880 !!gm this should be removed when choosing a unique strategy for fmask at the coast 881 881 ! If energy, enstrophy or mixed advection of momentum in vector form change the value for masks 882 882 ! at angles with three ocean points and one land point -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv.F90
r10068 r10802 75 75 CONTAINS 76 76 77 SUBROUTINE tra_adv( kt )77 SUBROUTINE tra_adv( kt, ktlev1, ktlev2, ktlev3, pts_rhs ) 78 78 !!---------------------------------------------------------------------- 79 79 !! *** ROUTINE tra_adv *** … … 81 81 !! ** Purpose : compute the ocean tracer advection trend. 82 82 !! 83 !! ** Method : - Update (ua,va) with the advection term following nadv 84 !!---------------------------------------------------------------------- 85 INTEGER, INTENT(in) :: kt ! ocean time-step index 83 !! ** Method : - Update (pu_rhs,pv_rhs) with the advection term following nadv 84 !!---------------------------------------------------------------------- 85 INTEGER, INTENT(in) :: kt ! ocean time-step index 86 INTEGER, INTENT(in) :: ktlev1, ktlev2, ktlev3 ! time level indices for source terms 87 REAL(wp), INTENT( inout), DIMENSION(jpi,jpj,jpk,jpts) :: pts_rhs ! temperature and salinity trends 86 88 ! 87 89 INTEGER :: jk ! dummy loop index … … 103 105 IF( ln_wave .AND. ln_sdw ) THEN 104 106 DO jk = 1, jpkm1 ! eulerian transport + Stokes Drift 105 zun(:,:,jk) = e2u (:,:) * e3u _n(:,:,jk) * ( un(:,:,jk) + usd(:,:,jk) )106 zvn(:,:,jk) = e1v (:,:) * e3v _n(:,:,jk) * ( vn(:,:,jk) + vsd(:,:,jk) )107 zwn(:,:,jk) = e1e2t(:,:) * ( w n(:,:,jk) + wsd(:,:,jk) )107 zun(:,:,jk) = e2u (:,:) * e3u(:,:,jk,ktlev2) * ( uu(:,:,jk,ktlev2) + usd(:,:,jk) ) 108 zvn(:,:,jk) = e1v (:,:) * e3v(:,:,jk,ktlev2) * ( vv(:,:,jk,ktlev2) + vsd(:,:,jk) ) 109 zwn(:,:,jk) = e1e2t(:,:) * ( ww(:,:,jk) + wsd(:,:,jk) ) 108 110 END DO 109 111 ELSE 110 112 DO jk = 1, jpkm1 111 zun(:,:,jk) = e2u (:,:) * e3u _n(:,:,jk) * un(:,:,jk) ! eulerian transport only112 zvn(:,:,jk) = e1v (:,:) * e3v _n(:,:,jk) * vn(:,:,jk)113 zwn(:,:,jk) = e1e2t(:,:) * w n(:,:,jk)113 zun(:,:,jk) = e2u (:,:) * e3u(:,:,jk,ktlev2) * uu(:,:,jk,ktlev2) ! eulerian transport only 114 zvn(:,:,jk) = e1v (:,:) * e3v(:,:,jk,ktlev2) * vv(:,:,jk,ktlev2) 115 zwn(:,:,jk) = e1e2t(:,:) * ww(:,:,jk) 114 116 END DO 115 117 ENDIF … … 139 141 IF( l_trdtra ) THEN !* Save ta and sa trends 140 142 ALLOCATE( ztrdt(jpi,jpj,jpk), ztrds(jpi,jpj,jpk) ) 141 ztrdt(:,:,:) = tsa(:,:,:,jp_tem)142 ztrds(:,:,:) = tsa(:,:,:,jp_sal)143 ztrdt(:,:,:) = pts_rhs(:,:,:,jp_tem) 144 ztrds(:,:,:) = pts_rhs(:,:,:,jp_sal) 143 145 ENDIF 144 146 ! … … 146 148 ! 147 149 CASE ( np_CEN ) ! Centered scheme : 2nd / 4th order 148 CALL tra_adv_cen ( kt, nit000, 'TRA', zun, zvn, zwn , tsn, tsa, jpts, nn_cen_h, nn_cen_v )150 CALL tra_adv_cen ( kt, nit000, ktlev2, 'TRA', zun, zvn, zwn , ts(:,:,:,:,ktlev2), pts_rhs, jpts, nn_cen_h, nn_cen_v ) 149 151 CASE ( np_FCT ) ! FCT scheme : 2nd / 4th order 150 CALL tra_adv_fct ( kt, nit000, 'TRA', r2dt, zun, zvn, zwn, tsb, tsn, tsa, jpts, nn_fct_h, nn_fct_v ) 152 CALL tra_adv_fct ( kt, nit000, ktlev1, ktlev2, ktlev3, 'TRA', r2dt, zun, zvn, zwn, & 153 & ts(:,:,:,:,ktlev1), ts(:,:,:,:,ktlev2), pts_rhs, jpts, nn_fct_h, nn_fct_v ) 151 154 CASE ( np_MUS ) ! MUSCL 152 CALL tra_adv_mus ( kt, nit000, 'TRA', r2dt, zun, zvn, zwn, tsb, tsa, jpts , ln_mus_ups )155 CALL tra_adv_mus ( kt, nit000, ktlev2, 'TRA', r2dt, zun, zvn, zwn, ts(:,:,:,:,ktlev1), pts_rhs, jpts , ln_mus_ups ) 153 156 CASE ( np_UBS ) ! UBS 154 CALL tra_adv_ubs ( kt, nit000, 'TRA', r2dt, zun, zvn, zwn, tsb, tsn, tsa, jpts , nn_ubs_v )157 CALL tra_adv_ubs ( kt, nit000, ktlev2, 'TRA', r2dt, zun, zvn, zwn, ts(:,:,:,:,ktlev1), ts(:,:,:,:,ktlev2), pts_rhs, jpts , nn_ubs_v ) 155 158 CASE ( np_QCK ) ! QUICKEST 156 CALL tra_adv_qck ( kt, nit000, 'TRA', r2dt, zun, zvn, zwn, tsb, tsn, tsa, jpts )159 CALL tra_adv_qck ( kt, nit000, ktlev2, 'TRA', r2dt, zun, zvn, zwn, ts(:,:,:,:,ktlev1), ts(:,:,:,:,ktlev2), pts_rhs, jpts ) 157 160 ! 158 161 END SELECT … … 160 163 IF( l_trdtra ) THEN ! save the advective trends for further diagnostics 161 164 DO jk = 1, jpkm1 162 ztrdt(:,:,jk) = tsa(:,:,jk,jp_tem) - ztrdt(:,:,jk)163 ztrds(:,:,jk) = tsa(:,:,jk,jp_sal) - ztrds(:,:,jk)165 ztrdt(:,:,jk) = pts_rhs(:,:,jk,jp_tem) - ztrdt(:,:,jk) 166 ztrds(:,:,jk) = pts_rhs(:,:,jk,jp_sal) - ztrds(:,:,jk) 164 167 END DO 165 168 CALL trd_tra( kt, 'TRA', jp_tem, jptra_totad, ztrdt ) -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_cen.F90
r10425 r10802 44 44 CONTAINS 45 45 46 SUBROUTINE tra_adv_cen( kt, kit000, cdtype, pun, pvn, pwn, &47 & ptn, pta, kjpt, kn_cen_h, kn_cen_v )46 SUBROUTINE tra_adv_cen( kt, kit000, ktlev, cdtype, pu, pv, pwn, & 47 & pt, pt_rhs, kjpt, kn_cen_h, kn_cen_v ) 48 48 !!---------------------------------------------------------------------- 49 49 !! *** ROUTINE tra_adv_cen *** … … 59 59 !! = 4 ==>> 4th order COMPACT scheme - - 60 60 !! 61 !! ** Action : - update pt awith the now advective tracer trends61 !! ** Action : - update pt_rhs with the now advective tracer trends 62 62 !! - send trends to trdtra module for further diagnostcs (l_trdtra=T) 63 63 !! - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) … … 65 65 INTEGER , INTENT(in ) :: kt ! ocean time-step index 66 66 INTEGER , INTENT(in ) :: kit000 ! first time step index 67 INTEGER , INTENT(in ) :: ktlev ! time level index for source terms 67 68 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 68 69 INTEGER , INTENT(in ) :: kjpt ! number of tracers 69 70 INTEGER , INTENT(in ) :: kn_cen_h ! =2/4 (2nd or 4th order scheme) 70 71 INTEGER , INTENT(in ) :: kn_cen_v ! =2/4 (2nd or 4th order scheme) 71 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu n, pvn, pwn ! 3 ocean velocity components72 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt n! now tracer fields73 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt a! tracer trend72 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu, pv, pwn ! 3 ocean velocity components 73 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt ! now tracer fields 74 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 74 75 ! 75 76 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 106 107 DO jj = 1, jpjm1 107 108 DO ji = 1, fs_jpim1 ! vector opt. 108 zwx(ji,jj,jk) = 0.5_wp * pu n(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji+1,jj ,jk,jn) )109 zwy(ji,jj,jk) = 0.5_wp * pv n(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji ,jj+1,jk,jn) )109 zwx(ji,jj,jk) = 0.5_wp * pu(ji,jj,jk) * ( pt(ji,jj,jk,jn) + pt(ji+1,jj ,jk,jn) ) 110 zwy(ji,jj,jk) = 0.5_wp * pv(ji,jj,jk) * ( pt(ji,jj,jk,jn) + pt(ji ,jj+1,jk,jn) ) 110 111 END DO 111 112 END DO … … 118 119 DO jj = 2, jpjm1 119 120 DO ji = fs_2, fs_jpim1 ! vector opt. 120 ztu(ji,jj,jk) = ( pt n(ji+1,jj ,jk,jn) - ptn(ji,jj,jk,jn) ) * umask(ji,jj,jk)121 ztv(ji,jj,jk) = ( pt n(ji ,jj+1,jk,jn) - ptn(ji,jj,jk,jn) ) * vmask(ji,jj,jk)121 ztu(ji,jj,jk) = ( pt(ji+1,jj ,jk,jn) - pt(ji,jj,jk,jn) ) * umask(ji,jj,jk) 122 ztv(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 122 123 END DO 123 124 END DO … … 128 129 DO jj = 2, jpjm1 129 130 DO ji = 1, fs_jpim1 ! vector opt. 130 zC2t_u = pt n(ji,jj,jk,jn) + ptn(ji+1,jj ,jk,jn) ! C2 interpolation of T at u- & v-points (x2)131 zC2t_v = pt n(ji,jj,jk,jn) + ptn(ji ,jj+1,jk,jn)131 zC2t_u = pt(ji,jj,jk,jn) + pt(ji+1,jj ,jk,jn) ! C2 interpolation of T at u- & v-points (x2) 132 zC2t_v = pt(ji,jj,jk,jn) + pt(ji ,jj+1,jk,jn) 132 133 ! ! C4 interpolation of T at u- & v-points (x2) 133 134 zC4t_u = zC2t_u + r1_6 * ( ztu(ji-1,jj,jk) - ztu(ji+1,jj,jk) ) 134 135 zC4t_v = zC2t_v + r1_6 * ( ztv(ji,jj-1,jk) - ztv(ji,jj+1,jk) ) 135 136 ! ! C4 fluxes 136 zwx(ji,jj,jk) = 0.5_wp * pu n(ji,jj,jk) * zC4t_u137 zwy(ji,jj,jk) = 0.5_wp * pv n(ji,jj,jk) * zC4t_v137 zwx(ji,jj,jk) = 0.5_wp * pu(ji,jj,jk) * zC4t_u 138 zwy(ji,jj,jk) = 0.5_wp * pv(ji,jj,jk) * zC4t_v 138 139 END DO 139 140 END DO … … 150 151 DO jj = 2, jpjm1 151 152 DO ji = fs_2, fs_jpim1 ! vector opt. 152 zwz(ji,jj,jk) = 0.5 * pwn(ji,jj,jk) * ( pt n(ji,jj,jk,jn) + ptn(ji,jj,jk-1,jn) ) * wmask(ji,jj,jk)153 zwz(ji,jj,jk) = 0.5 * pwn(ji,jj,jk) * ( pt(ji,jj,jk,jn) + pt(ji,jj,jk-1,jn) ) * wmask(ji,jj,jk) 153 154 END DO 154 155 END DO … … 156 157 ! 157 158 CASE( 4 ) !* 4th order compact 158 CALL interp_4th_cpt( pt n(:,:,:,jn) , ztw ) ! ztw = interpolated value of T at w-point159 CALL interp_4th_cpt( pt(:,:,:,jn) , ztw ) ! ztw = interpolated value of T at w-point 159 160 DO jk = 2, jpkm1 160 161 DO jj = 2, jpjm1 … … 171 172 DO jj = 1, jpj 172 173 DO ji = 1, jpi 173 zwz(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * pt n(ji,jj,mikt(ji,jj),jn)174 zwz(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn) 174 175 END DO 175 176 END DO 176 177 ELSE ! no ice-shelf cavities (only ocean surface) 177 zwz(:,:,1) = pwn(:,:,1) * pt n(:,:,1,jn)178 zwz(:,:,1) = pwn(:,:,1) * pt(:,:,1,jn) 178 179 ENDIF 179 180 ENDIF … … 182 183 DO jj = 2, jpjm1 183 184 DO ji = fs_2, fs_jpim1 ! vector opt. 184 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) &185 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) & 185 186 & - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & 186 187 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) & 187 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) * r1_e1e2t(ji,jj) / e3t _n(ji,jj,jk)188 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev) 188 189 END DO 189 190 END DO … … 191 192 ! ! trend diagnostics 192 193 IF( l_trd ) THEN 193 CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pu n, ptn(:,:,:,jn) )194 CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pv n, ptn(:,:,:,jn) )195 CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pwn, pt n(:,:,:,jn) )194 CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pu, pt(:,:,:,jn) ) 195 CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pv, pt(:,:,:,jn) ) 196 CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pwn, pt(:,:,:,jn) ) 196 197 END IF 197 198 ! ! "Poleward" heat and salt transports -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_fct.F90
r10425 r10802 52 52 CONTAINS 53 53 54 SUBROUTINE tra_adv_fct( kt, kit000, cdtype, p2dt, pun, pvn, pwn, &55 & ptb, ptn, pta, kjpt, kn_fct_h, kn_fct_v )54 SUBROUTINE tra_adv_fct( kt, kit000, ktlev1, ktlev2, ktlev3, cdtype, p2dt, pu, pv, pwn, & 55 & pt_lev1, pt_lev2, pt_rhs, kjpt, kn_fct_h, kn_fct_v ) 56 56 !!---------------------------------------------------------------------- 57 57 !! *** ROUTINE tra_adv_fct *** … … 65 65 !! - corrected flux (monotonic correction) 66 66 !! 67 !! ** Action : - update pt awith the now advective tracer trends67 !! ** Action : - update pt_rhs with the now advective tracer trends 68 68 !! - send trends to trdtra module for further diagnostics (l_trdtra=T) 69 69 !! - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) 70 70 !!---------------------------------------------------------------------- 71 71 INTEGER , INTENT(in ) :: kt ! ocean time-step index 72 INTEGER , INTENT(in ) :: ktlev1, ktlev2, ktlev3 ! time level indices for source terms 72 73 INTEGER , INTENT(in ) :: kit000 ! first time step index 73 74 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) … … 76 77 INTEGER , INTENT(in ) :: kn_fct_v ! order of the FCT scheme (=2 or 4) 77 78 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 78 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu n, pvn, pwn ! 3 ocean velocity components79 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt b, ptn! before and now tracer fields80 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt a! tracer trend79 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu, pv, pwn ! 3 ocean velocity components 80 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt_lev1, pt_lev2 ! before and now tracer fields 81 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 81 82 ! 82 83 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 125 126 DO ji = 1, fs_jpim1 ! vector opt. 126 127 ! upstream scheme 127 zfp_ui = pu n(ji,jj,jk) + ABS( pun(ji,jj,jk) )128 zfm_ui = pu n(ji,jj,jk) - ABS( pun(ji,jj,jk) )129 zfp_vj = pv n(ji,jj,jk) + ABS( pvn(ji,jj,jk) )130 zfm_vj = pv n(ji,jj,jk) - ABS( pvn(ji,jj,jk) )131 zwx(ji,jj,jk) = 0.5 * ( zfp_ui * pt b(ji,jj,jk,jn) + zfm_ui * ptb(ji+1,jj ,jk,jn) )132 zwy(ji,jj,jk) = 0.5 * ( zfp_vj * pt b(ji,jj,jk,jn) + zfm_vj * ptb(ji ,jj+1,jk,jn) )128 zfp_ui = pu(ji,jj,jk) + ABS( pu(ji,jj,jk) ) 129 zfm_ui = pu(ji,jj,jk) - ABS( pu(ji,jj,jk) ) 130 zfp_vj = pv(ji,jj,jk) + ABS( pv(ji,jj,jk) ) 131 zfm_vj = pv(ji,jj,jk) - ABS( pv(ji,jj,jk) ) 132 zwx(ji,jj,jk) = 0.5 * ( zfp_ui * pt_lev1(ji,jj,jk,jn) + zfm_ui * pt_lev1(ji+1,jj ,jk,jn) ) 133 zwy(ji,jj,jk) = 0.5 * ( zfp_vj * pt_lev1(ji,jj,jk,jn) + zfm_vj * pt_lev1(ji ,jj+1,jk,jn) ) 133 134 END DO 134 135 END DO … … 140 141 zfp_wk = pwn(ji,jj,jk) + ABS( pwn(ji,jj,jk) ) 141 142 zfm_wk = pwn(ji,jj,jk) - ABS( pwn(ji,jj,jk) ) 142 zwz(ji,jj,jk) = 0.5 * ( zfp_wk * pt b(ji,jj,jk,jn) + zfm_wk * ptb(ji,jj,jk-1,jn) ) * wmask(ji,jj,jk)143 zwz(ji,jj,jk) = 0.5 * ( zfp_wk * pt_lev1(ji,jj,jk,jn) + zfm_wk * pt_lev1(ji,jj,jk-1,jn) ) * wmask(ji,jj,jk) 143 144 END DO 144 145 END DO … … 148 149 DO jj = 1, jpj 149 150 DO ji = 1, jpi 150 zwz(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * pt b(ji,jj,mikt(ji,jj),jn) ! linear free surface151 zwz(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * pt_lev1(ji,jj,mikt(ji,jj),jn) ! linear free surface 151 152 END DO 152 153 END DO 153 154 ELSE ! no cavities: only at the ocean surface 154 zwz(:,:,1) = pwn(:,:,1) * pt b(:,:,1,jn)155 zwz(:,:,1) = pwn(:,:,1) * pt_lev1(:,:,1,jn) 155 156 ENDIF 156 157 ENDIF … … 164 165 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) * r1_e1e2t(ji,jj) 165 166 ! ! update and guess with monotonic sheme 166 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra / e3t_n(ji,jj,jk) * tmask(ji,jj,jk)167 zwi(ji,jj,jk) = ( e3t _b(ji,jj,jk) * ptb(ji,jj,jk,jn) + p2dt * ztra ) / e3t_a(ji,jj,jk) * tmask(ji,jj,jk)167 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + ztra / e3t(ji,jj,jk,ktlev2) * tmask(ji,jj,jk) 168 zwi(ji,jj,jk) = ( e3t(ji,jj,jk,ktlev1) * pt_lev1(ji,jj,jk,jn) + p2dt * ztra ) / e3t(ji,jj,jk,ktlev3) * tmask(ji,jj,jk) 168 169 END DO 169 170 END DO … … 184 185 DO jj = 1, jpjm1 185 186 DO ji = 1, fs_jpim1 ! vector opt. 186 zwx(ji,jj,jk) = 0.5_wp * pu n(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji+1,jj,jk,jn) ) - zwx(ji,jj,jk)187 zwy(ji,jj,jk) = 0.5_wp * pv n(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji,jj+1,jk,jn) ) - zwy(ji,jj,jk)187 zwx(ji,jj,jk) = 0.5_wp * pu(ji,jj,jk) * ( pt_lev2(ji,jj,jk,jn) + pt_lev2(ji+1,jj,jk,jn) ) - zwx(ji,jj,jk) 188 zwy(ji,jj,jk) = 0.5_wp * pv(ji,jj,jk) * ( pt_lev2(ji,jj,jk,jn) + pt_lev2(ji,jj+1,jk,jn) ) - zwy(ji,jj,jk) 188 189 END DO 189 190 END DO … … 196 197 DO jj = 1, jpjm1 ! 1st derivative (gradient) 197 198 DO ji = 1, fs_jpim1 ! vector opt. 198 ztu(ji,jj,jk) = ( pt n(ji+1,jj ,jk,jn) - ptn(ji,jj,jk,jn) ) * umask(ji,jj,jk)199 ztv(ji,jj,jk) = ( pt n(ji ,jj+1,jk,jn) - ptn(ji,jj,jk,jn) ) * vmask(ji,jj,jk)199 ztu(ji,jj,jk) = ( pt_lev2(ji+1,jj ,jk,jn) - pt_lev2(ji,jj,jk,jn) ) * umask(ji,jj,jk) 200 ztv(ji,jj,jk) = ( pt_lev2(ji ,jj+1,jk,jn) - pt_lev2(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 200 201 END DO 201 202 END DO … … 212 213 DO jj = 1, jpjm1 213 214 DO ji = 1, fs_jpim1 ! vector opt. 214 zC2t_u = pt n(ji,jj,jk,jn) + ptn(ji+1,jj ,jk,jn) ! 2 x C2 interpolation of T at u- & v-points215 zC2t_v = pt n(ji,jj,jk,jn) + ptn(ji ,jj+1,jk,jn)215 zC2t_u = pt_lev2(ji,jj,jk,jn) + pt_lev2(ji+1,jj ,jk,jn) ! 2 x C2 interpolation of T at u- & v-points 216 zC2t_v = pt_lev2(ji,jj,jk,jn) + pt_lev2(ji ,jj+1,jk,jn) 216 217 ! ! C4 minus upstream advective fluxes 217 zwx(ji,jj,jk) = 0.5_wp * pu n(ji,jj,jk) * ( zC2t_u + zltu(ji,jj,jk) - zltu(ji+1,jj,jk) ) - zwx(ji,jj,jk)218 zwy(ji,jj,jk) = 0.5_wp * pv n(ji,jj,jk) * ( zC2t_v + zltv(ji,jj,jk) - zltv(ji,jj+1,jk) ) - zwy(ji,jj,jk)218 zwx(ji,jj,jk) = 0.5_wp * pu(ji,jj,jk) * ( zC2t_u + zltu(ji,jj,jk) - zltu(ji+1,jj,jk) ) - zwx(ji,jj,jk) 219 zwy(ji,jj,jk) = 0.5_wp * pv(ji,jj,jk) * ( zC2t_v + zltv(ji,jj,jk) - zltv(ji,jj+1,jk) ) - zwy(ji,jj,jk) 219 220 END DO 220 221 END DO … … 227 228 DO jj = 1, jpjm1 228 229 DO ji = 1, fs_jpim1 ! vector opt. 229 ztu(ji,jj,jk) = ( pt n(ji+1,jj ,jk,jn) - ptn(ji,jj,jk,jn) ) * umask(ji,jj,jk)230 ztv(ji,jj,jk) = ( pt n(ji ,jj+1,jk,jn) - ptn(ji,jj,jk,jn) ) * vmask(ji,jj,jk)230 ztu(ji,jj,jk) = ( pt_lev2(ji+1,jj ,jk,jn) - pt_lev2(ji,jj,jk,jn) ) * umask(ji,jj,jk) 231 ztv(ji,jj,jk) = ( pt_lev2(ji ,jj+1,jk,jn) - pt_lev2(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 231 232 END DO 232 233 END DO … … 237 238 DO jj = 2, jpjm1 238 239 DO ji = 2, fs_jpim1 ! vector opt. 239 zC2t_u = pt n(ji,jj,jk,jn) + ptn(ji+1,jj ,jk,jn) ! 2 x C2 interpolation of T at u- & v-points (x2)240 zC2t_v = pt n(ji,jj,jk,jn) + ptn(ji ,jj+1,jk,jn)240 zC2t_u = pt_lev2(ji,jj,jk,jn) + pt_lev2(ji+1,jj ,jk,jn) ! 2 x C2 interpolation of T at u- & v-points (x2) 241 zC2t_v = pt_lev2(ji,jj,jk,jn) + pt_lev2(ji ,jj+1,jk,jn) 241 242 ! ! C4 interpolation of T at u- & v-points (x2) 242 243 zC4t_u = zC2t_u + r1_6 * ( ztu(ji-1,jj ,jk) - ztu(ji+1,jj ,jk) ) 243 244 zC4t_v = zC2t_v + r1_6 * ( ztv(ji ,jj-1,jk) - ztv(ji ,jj+1,jk) ) 244 245 ! ! C4 minus upstream advective fluxes 245 zwx(ji,jj,jk) = 0.5_wp * pu n(ji,jj,jk) * zC4t_u - zwx(ji,jj,jk)246 zwy(ji,jj,jk) = 0.5_wp * pv n(ji,jj,jk) * zC4t_v - zwy(ji,jj,jk)246 zwx(ji,jj,jk) = 0.5_wp * pu(ji,jj,jk) * zC4t_u - zwx(ji,jj,jk) 247 zwy(ji,jj,jk) = 0.5_wp * pv(ji,jj,jk) * zC4t_v - zwy(ji,jj,jk) 247 248 END DO 248 249 END DO … … 257 258 DO jj = 2, jpjm1 258 259 DO ji = fs_2, fs_jpim1 259 zwz(ji,jj,jk) = ( pwn(ji,jj,jk) * 0.5_wp * ( pt n(ji,jj,jk,jn) + ptn(ji,jj,jk-1,jn) ) &260 zwz(ji,jj,jk) = ( pwn(ji,jj,jk) * 0.5_wp * ( pt_lev2(ji,jj,jk,jn) + pt_lev2(ji,jj,jk-1,jn) ) & 260 261 & - zwz(ji,jj,jk) ) * wmask(ji,jj,jk) 261 262 END DO … … 264 265 ! 265 266 CASE( 4 ) !- 4th order COMPACT 266 CALL interp_4th_cpt( pt n(:,:,:,jn) , ztw ) ! zwt = COMPACT interpolation of T at w-point267 CALL interp_4th_cpt( pt_lev2(:,:,:,jn) , ztw ) ! zwt = COMPACT interpolation of T at w-point 267 268 DO jk = 2, jpkm1 268 269 DO jj = 2, jpjm1 … … 282 283 ! !== monotonicity algorithm ==! 283 284 ! 284 CALL nonosc( pt b(:,:,:,jn), zwx, zwy, zwz, zwi, p2dt)285 CALL nonosc( pt_lev1(:,:,:,jn), zwx, zwy, zwz, zwi, p2dt, e3t(:,:,:,ktlev2) ) 285 286 ! 286 287 ! !== final trend with corrected fluxes ==! … … 289 290 DO jj = 2, jpjm1 290 291 DO ji = fs_2, fs_jpim1 ! vector opt. 291 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) &292 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & 292 293 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) & 293 294 & + zwz(ji,jj,jk) - zwz(ji ,jj ,jk+1) ) & 294 & * r1_e1e2t(ji,jj) / e3t _n(ji,jj,jk)295 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev2) 295 296 END DO 296 297 END DO … … 303 304 ! 304 305 IF( l_trd ) THEN ! trend diagnostics 305 CALL trd_tra( kt, cdtype, jn, jptra_xad, ztrdx, pu n, ptn(:,:,:,jn) )306 CALL trd_tra( kt, cdtype, jn, jptra_yad, ztrdy, pv n, ptn(:,:,:,jn) )307 CALL trd_tra( kt, cdtype, jn, jptra_zad, ztrdz, pwn, pt n(:,:,:,jn) )306 CALL trd_tra( kt, cdtype, jn, jptra_xad, ztrdx, pu, pt_lev2(:,:,:,jn) ) 307 CALL trd_tra( kt, cdtype, jn, jptra_yad, ztrdy, pv, pt_lev2(:,:,:,jn) ) 308 CALL trd_tra( kt, cdtype, jn, jptra_zad, ztrdz, pwn, pt_lev2(:,:,:,jn) ) 308 309 ENDIF 309 310 ! ! heat/salt transport … … 328 329 329 330 330 SUBROUTINE nonosc( pbef, paa, pbb, pcc, paft, p2dt )331 SUBROUTINE nonosc( pbef, paa, pbb, pcc, paft, p2dt, pe3t ) 331 332 !!--------------------------------------------------------------------- 332 333 !! *** ROUTINE nonosc *** … … 341 342 !! in-space based differencing for fluid 342 343 !!---------------------------------------------------------------------- 343 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step344 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(in ) :: pbef, paft ! before & afterfield345 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(inout) :: paa, pbb, pcc ! monotonic fluxes in the 3 directions344 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 345 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(in ) :: pbef, paft, pe3t ! before & after field, now e3t field 346 REAL(wp), DIMENSION (jpi,jpj,jpk), INTENT(inout) :: paa, pbb, pcc ! monotonic fluxes in the 3 directions 346 347 ! 347 348 INTEGER :: ji, jj, jk ! dummy loop indices … … 392 393 393 394 ! up & down beta terms 394 zbt = e1e2t(ji,jj) * e3t_n(ji,jj,jk) / p2dt395 zbt = e1e2t(ji,jj) * pe3t(ji,jj,jk) / p2dt 395 396 zbetup(ji,jj,jk) = ( zup - paft(ji,jj,jk) ) / ( zpos + zrtrn ) * zbt 396 397 zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - zdo ) / ( zneg + zrtrn ) * zbt … … 634 635 !! The tri-diagonals matrix is given as input 3D arrays: pD, pU, pL 635 636 !! (i.e. the Diagonal, the Upper diagonal, and the Lower diagonal). 636 !! The solution is pt a.637 !! The solution is pt_rhs. 637 638 !! The 3d array zwt is used as a work space array. 638 639 !!---------------------------------------------------------------------- -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_mus.F90
r10425 r10802 54 54 CONTAINS 55 55 56 SUBROUTINE tra_adv_mus( kt, kit000, cdtype, p2dt, pun, pvn, pwn, &57 & ptb, pta, kjpt, ld_msc_ups )56 SUBROUTINE tra_adv_mus( kt, kit000, ktlev, cdtype, p2dt, pu, pv, pwn, & 57 & pt, pt_rhs, kjpt, ld_msc_ups ) 58 58 !!---------------------------------------------------------------------- 59 59 !! *** ROUTINE tra_adv_mus *** … … 66 66 !! ld_msc_ups=T : 67 67 !! 68 !! ** Action : - update pt awith the now advective tracer trends68 !! ** Action : - update pt_rhs with the now advective tracer trends 69 69 !! - send trends to trdtra module for further diagnostcs (l_trdtra=T) 70 70 !! - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) … … 75 75 INTEGER , INTENT(in ) :: kt ! ocean time-step index 76 76 INTEGER , INTENT(in ) :: kit000 ! first time step index 77 INTEGER , INTENT(in ) :: ktlev ! time level index for source terms 77 78 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 78 79 INTEGER , INTENT(in ) :: kjpt ! number of tracers 79 80 LOGICAL , INTENT(in ) :: ld_msc_ups ! use upstream scheme within muscl 80 81 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 81 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu n, pvn, pwn ! 3 ocean velocity components82 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt b! before tracer field83 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt a! tracer trend82 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu, pv, pwn ! 3 ocean velocity components 83 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt ! before tracer field 84 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 84 85 ! 85 86 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 134 135 DO jj = 1, jpjm1 135 136 DO ji = 1, fs_jpim1 ! vector opt. 136 zwx(ji,jj,jk) = umask(ji,jj,jk) * ( pt b(ji+1,jj,jk,jn) - ptb(ji,jj,jk,jn) )137 zwy(ji,jj,jk) = vmask(ji,jj,jk) * ( pt b(ji,jj+1,jk,jn) - ptb(ji,jj,jk,jn) )137 zwx(ji,jj,jk) = umask(ji,jj,jk) * ( pt(ji+1,jj,jk,jn) - pt(ji,jj,jk,jn) ) 138 zwy(ji,jj,jk) = vmask(ji,jj,jk) * ( pt(ji,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) 138 139 END DO 139 140 END DO … … 172 173 DO ji = fs_2, fs_jpim1 ! vector opt. 173 174 ! MUSCL fluxes 174 z0u = SIGN( 0.5, pu n(ji,jj,jk) )175 z0u = SIGN( 0.5, pu(ji,jj,jk) ) 175 176 zalpha = 0.5 - z0u 176 zu = z0u - 0.5 * pu n(ji,jj,jk) * p2dt * r1_e1e2u(ji,jj) / e3u_n(ji,jj,jk)177 zzwx = pt b(ji+1,jj,jk,jn) + xind(ji,jj,jk) * zu * zslpx(ji+1,jj,jk)178 zzwy = pt b(ji ,jj,jk,jn) + xind(ji,jj,jk) * zu * zslpx(ji ,jj,jk)179 zwx(ji,jj,jk) = pu n(ji,jj,jk) * ( zalpha * zzwx + (1.-zalpha) * zzwy )177 zu = z0u - 0.5 * pu(ji,jj,jk) * p2dt * r1_e1e2u(ji,jj) / e3u(ji,jj,jk,ktlev) 178 zzwx = pt(ji+1,jj,jk,jn) + xind(ji,jj,jk) * zu * zslpx(ji+1,jj,jk) 179 zzwy = pt(ji ,jj,jk,jn) + xind(ji,jj,jk) * zu * zslpx(ji ,jj,jk) 180 zwx(ji,jj,jk) = pu(ji,jj,jk) * ( zalpha * zzwx + (1.-zalpha) * zzwy ) 180 181 ! 181 z0v = SIGN( 0.5, pv n(ji,jj,jk) )182 z0v = SIGN( 0.5, pv(ji,jj,jk) ) 182 183 zalpha = 0.5 - z0v 183 zv = z0v - 0.5 * pv n(ji,jj,jk) * p2dt * r1_e1e2v(ji,jj) / e3v_n(ji,jj,jk)184 zzwx = pt b(ji,jj+1,jk,jn) + xind(ji,jj,jk) * zv * zslpy(ji,jj+1,jk)185 zzwy = pt b(ji,jj ,jk,jn) + xind(ji,jj,jk) * zv * zslpy(ji,jj ,jk)186 zwy(ji,jj,jk) = pv n(ji,jj,jk) * ( zalpha * zzwx + (1.-zalpha) * zzwy )184 zv = z0v - 0.5 * pv(ji,jj,jk) * p2dt * r1_e1e2v(ji,jj) / e3v(ji,jj,jk,ktlev) 185 zzwx = pt(ji,jj+1,jk,jn) + xind(ji,jj,jk) * zv * zslpy(ji,jj+1,jk) 186 zzwy = pt(ji,jj ,jk,jn) + xind(ji,jj,jk) * zv * zslpy(ji,jj ,jk) 187 zwy(ji,jj,jk) = pv(ji,jj,jk) * ( zalpha * zzwx + (1.-zalpha) * zzwy ) 187 188 END DO 188 189 END DO … … 193 194 DO jj = 2, jpjm1 194 195 DO ji = fs_2, fs_jpim1 ! vector opt. 195 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) &196 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) - ( zwx(ji,jj,jk) - zwx(ji-1,jj ,jk ) & 196 197 & + zwy(ji,jj,jk) - zwy(ji ,jj-1,jk ) ) & 197 & * r1_e1e2t(ji,jj) / e3t _n(ji,jj,jk)198 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev) 198 199 END DO 199 200 END DO … … 201 202 ! ! trend diagnostics 202 203 IF( l_trd ) THEN 203 CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pu n, ptb(:,:,:,jn) )204 CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pv n, ptb(:,:,:,jn) )204 CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pu, pt(:,:,:,jn) ) 205 CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pv, pt(:,:,:,jn) ) 205 206 END IF 206 207 ! ! "Poleward" heat and salt transports … … 215 216 zwx(:,:,jpk) = 0._wp 216 217 DO jk = 2, jpkm1 ! interior values 217 zwx(:,:,jk) = tmask(:,:,jk) * ( pt b(:,:,jk-1,jn) - ptb(:,:,jk,jn) )218 zwx(:,:,jk) = tmask(:,:,jk) * ( pt(:,:,jk-1,jn) - pt(:,:,jk,jn) ) 218 219 END DO 219 220 ! !-- Slopes of tracer … … 242 243 zalpha = 0.5 + z0w 243 244 zw = z0w - 0.5 * pwn(ji,jj,jk+1) * p2dt * r1_e1e2t(ji,jj) / e3w_n(ji,jj,jk+1) 244 zzwx = pt b(ji,jj,jk+1,jn) + xind(ji,jj,jk) * zw * zslpx(ji,jj,jk+1)245 zzwy = pt b(ji,jj,jk ,jn) + xind(ji,jj,jk) * zw * zslpx(ji,jj,jk )245 zzwx = pt(ji,jj,jk+1,jn) + xind(ji,jj,jk) * zw * zslpx(ji,jj,jk+1) 246 zzwy = pt(ji,jj,jk ,jn) + xind(ji,jj,jk) * zw * zslpx(ji,jj,jk ) 246 247 zwx(ji,jj,jk+1) = pwn(ji,jj,jk+1) * ( zalpha * zzwx + (1.-zalpha) * zzwy ) * wmask(ji,jj,jk) 247 248 END DO … … 252 253 DO jj = 1, jpj 253 254 DO ji = 1, jpi 254 zwx(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * pt b(ji,jj,mikt(ji,jj),jn)255 zwx(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn) 255 256 END DO 256 257 END DO 257 258 ELSE ! no cavities: only at the ocean surface 258 zwx(:,:,1) = pwn(:,:,1) * pt b(:,:,1,jn)259 zwx(:,:,1) = pwn(:,:,1) * pt(:,:,1,jn) 259 260 ENDIF 260 261 ENDIF … … 263 264 DO jj = 2, jpjm1 264 265 DO ji = fs_2, fs_jpim1 ! vector opt. 265 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - ( zwx(ji,jj,jk) - zwx(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk)266 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) - ( zwx(ji,jj,jk) - zwx(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev) 266 267 END DO 267 268 END DO 268 269 END DO 269 270 ! ! send trends for diagnostic 270 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_zad, zwx, pwn, pt b(:,:,:,jn) )271 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_zad, zwx, pwn, pt(:,:,:,jn) ) 271 272 ! 272 273 END DO ! end of tracer loop -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_qck.F90
r10425 r10802 47 47 CONTAINS 48 48 49 SUBROUTINE tra_adv_qck ( kt, kit000, cdtype, p2dt, pun, pvn, pwn, &50 & ptb, ptn, pta, kjpt )49 SUBROUTINE tra_adv_qck ( kt, kit000, ktlev, cdtype, p2dt, pu, pv, pwn, & 50 & pt_lev1, pt_lev2, pt_rhs, kjpt ) 51 51 !!---------------------------------------------------------------------- 52 52 !! *** ROUTINE tra_adv_qck *** … … 72 72 !! dt = 2*rdtra and the scalar values are tb and sb 73 73 !! 74 !! On the vertical, the simple centered scheme used pt n74 !! On the vertical, the simple centered scheme used pt_lev2 75 75 !! 76 76 !! The fluxes are bounded by the ULTIMATE limiter to … … 78 78 !! prevent the appearance of spurious numerical oscillations 79 79 !! 80 !! ** Action : - update pt awith the now advective tracer trends80 !! ** Action : - update pt_rhs with the now advective tracer trends 81 81 !! - send trends to trdtra module for further diagnostcs (l_trdtra=T) 82 82 !! - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) … … 86 86 INTEGER , INTENT(in ) :: kt ! ocean time-step index 87 87 INTEGER , INTENT(in ) :: kit000 ! first time step index 88 INTEGER , INTENT(in ) :: ktlev ! time level index for source terms 88 89 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 89 90 INTEGER , INTENT(in ) :: kjpt ! number of tracers 90 91 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 91 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu n, pvn, pwn ! 3 ocean velocity components92 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt b, ptn! before and now tracer fields93 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt a! tracer trend92 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu, pv, pwn ! 3 ocean velocity components 93 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt_lev1, pt_lev2 ! before and now tracer fields 94 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 94 95 !!---------------------------------------------------------------------- 95 96 ! … … 108 109 ! 109 110 ! ! horizontal fluxes are computed with the QUICKEST + ULTIMATE scheme 110 CALL tra_adv_qck_i( kt, cdtype, p2dt, pun, ptb, ptn, pta, kjpt )111 CALL tra_adv_qck_j( kt, cdtype, p2dt, pvn, ptb, ptn, pta, kjpt )111 CALL tra_adv_qck_i( kt, ktlev, cdtype, p2dt, pu, pt_lev1, pt_lev2, pt_rhs, kjpt ) 112 CALL tra_adv_qck_j( kt, ktlev, cdtype, p2dt, pv, pt_lev1, pt_lev2, pt_rhs, kjpt ) 112 113 113 114 ! ! vertical fluxes are computed with the 2nd order centered scheme 114 CALL tra_adv_cen2_k( kt, cdtype, pwn, ptn, pta, kjpt )115 CALL tra_adv_cen2_k( kt, ktlev, cdtype, pwn, pt_lev2, pt_rhs, kjpt ) 115 116 ! 116 117 END SUBROUTINE tra_adv_qck 117 118 118 119 119 SUBROUTINE tra_adv_qck_i( kt, cdtype, p2dt, pun, &120 & pt b, ptn, pta, kjpt )120 SUBROUTINE tra_adv_qck_i( kt, ktlev, cdtype, p2dt, pu, & 121 & pt_lev1, pt_lev2, pt_rhs, kjpt ) 121 122 !!---------------------------------------------------------------------- 122 123 !! 123 124 !!---------------------------------------------------------------------- 124 125 INTEGER , INTENT(in ) :: kt ! ocean time-step index 126 INTEGER , INTENT(in ) :: ktlev ! time level index for source terms 125 127 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 126 128 INTEGER , INTENT(in ) :: kjpt ! number of tracers 127 129 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 128 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu n! i-velocity components129 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt b, ptn! before and now tracer fields130 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt a! tracer trend130 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu ! i-velocity components 131 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt_lev1, pt_lev2 ! before and now tracer fields 132 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 131 133 !! 132 134 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 145 147 DO jj = 2, jpjm1 146 148 DO ji = fs_2, fs_jpim1 ! vector opt. 147 zfc(ji,jj,jk) = pt b(ji-1,jj,jk,jn) ! Upstream in the x-direction for the tracer148 zfd(ji,jj,jk) = pt b(ji+1,jj,jk,jn) ! Downstream in the x-direction for the tracer149 zfc(ji,jj,jk) = pt_lev1(ji-1,jj,jk,jn) ! Upstream in the x-direction for the tracer 150 zfd(ji,jj,jk) = pt_lev1(ji+1,jj,jk,jn) ! Downstream in the x-direction for the tracer 149 151 END DO 150 152 END DO … … 158 160 DO jj = 2, jpjm1 159 161 DO ji = fs_2, fs_jpim1 ! vector opt. 160 zdir = 0.5 + SIGN( 0.5, pu n(ji,jj,jk) ) ! if pun> 0 : zdir = 1 otherwise zdir = 0162 zdir = 0.5 + SIGN( 0.5, pu(ji,jj,jk) ) ! if pu > 0 : zdir = 1 otherwise zdir = 0 161 163 zfu(ji,jj,jk) = zdir * zfc(ji,jj,jk ) + ( 1. - zdir ) * zfd(ji+1,jj,jk) ! FU in the x-direction for T 162 164 END DO … … 167 169 DO jj = 2, jpjm1 168 170 DO ji = fs_2, fs_jpim1 ! vector opt. 169 zdir = 0.5 + SIGN( 0.5, pu n(ji,jj,jk) ) ! if pun> 0 : zdir = 1 otherwise zdir = 0170 zdx = ( zdir * e1t(ji,jj) + ( 1. - zdir ) * e1t(ji+1,jj) ) * e2u(ji,jj) * e3u _n(ji,jj,jk)171 zwx(ji,jj,jk) = ABS( pu n(ji,jj,jk) ) * p2dt / zdx ! (0<zc_cfl<1 : Courant number on x-direction)172 zfc(ji,jj,jk) = zdir * pt b(ji ,jj,jk,jn) + ( 1. - zdir ) * ptb(ji+1,jj,jk,jn) ! FC in the x-direction for T173 zfd(ji,jj,jk) = zdir * pt b(ji+1,jj,jk,jn) + ( 1. - zdir ) * ptb(ji ,jj,jk,jn) ! FD in the x-direction for T171 zdir = 0.5 + SIGN( 0.5, pu(ji,jj,jk) ) ! if pu > 0 : zdir = 1 otherwise zdir = 0 172 zdx = ( zdir * e1t(ji,jj) + ( 1. - zdir ) * e1t(ji+1,jj) ) * e2u(ji,jj) * e3u(ji,jj,jk,ktlev) 173 zwx(ji,jj,jk) = ABS( pu(ji,jj,jk) ) * p2dt / zdx ! (0<zc_cfl<1 : Courant number on x-direction) 174 zfc(ji,jj,jk) = zdir * pt_lev1(ji ,jj,jk,jn) + ( 1. - zdir ) * pt_lev1(ji+1,jj,jk,jn) ! FC in the x-direction for T 175 zfd(ji,jj,jk) = zdir * pt_lev1(ji+1,jj,jk,jn) + ( 1. - zdir ) * pt_lev1(ji ,jj,jk,jn) ! FD in the x-direction for T 174 176 END DO 175 177 END DO … … 197 199 DO jj = 2, jpjm1 198 200 DO ji = fs_2, fs_jpim1 ! vector opt. 199 zdir = 0.5 + SIGN( 0.5, pu n(ji,jj,jk) ) ! if pun> 0 : zdir = 1 otherwise zdir = 0201 zdir = 0.5 + SIGN( 0.5, pu(ji,jj,jk) ) ! if pu > 0 : zdir = 1 otherwise zdir = 0 200 202 !--- If the second ustream point is a land point 201 203 !--- the flux is computed by the 1st order UPWIND scheme 202 204 zmsk = zdir * zfu(ji,jj,jk) + ( 1. - zdir ) * zfu(ji+1,jj,jk) 203 205 zwx(ji,jj,jk) = zmsk * zwx(ji,jj,jk) + ( 1. - zmsk ) * zfc(ji,jj,jk) 204 zwx(ji,jj,jk) = zwx(ji,jj,jk) * pu n(ji,jj,jk)206 zwx(ji,jj,jk) = zwx(ji,jj,jk) * pu(ji,jj,jk) 205 207 END DO 206 208 END DO … … 213 215 DO jj = 2, jpjm1 214 216 DO ji = fs_2, fs_jpim1 ! vector opt. 215 zbtr = r1_e1e2t(ji,jj) / e3t _n(ji,jj,jk)217 zbtr = r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev) 216 218 ! horizontal advective trends 217 219 ztra = - zbtr * ( zwx(ji,jj,jk) - zwx(ji-1,jj,jk) ) 218 220 !--- add it to the general tracer trends 219 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra221 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + ztra 220 222 END DO 221 223 END DO 222 224 END DO 223 225 ! ! trend diagnostics 224 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pu n, ptn(:,:,:,jn) )226 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pu, pt_lev2(:,:,:,jn) ) 225 227 ! 226 228 END DO … … 229 231 230 232 231 SUBROUTINE tra_adv_qck_j( kt, cdtype, p2dt, pvn, &232 & pt b, ptn, pta, kjpt )233 SUBROUTINE tra_adv_qck_j( kt, ktlev, cdtype, p2dt, pv, & 234 & pt_lev1, pt_lev2, pt_rhs, kjpt ) 233 235 !!---------------------------------------------------------------------- 234 236 !! 235 237 !!---------------------------------------------------------------------- 236 238 INTEGER , INTENT(in ) :: kt ! ocean time-step index 239 INTEGER , INTENT(in ) :: ktlev ! time level index for source terms 237 240 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 238 241 INTEGER , INTENT(in ) :: kjpt ! number of tracers 239 242 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 240 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pv n! j-velocity components241 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt b, ptn! before and now tracer fields242 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt a! tracer trend243 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pv ! j-velocity components 244 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt_lev1, pt_lev2 ! before and now tracer fields 245 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 243 246 !! 244 247 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 259 262 DO ji = fs_2, fs_jpim1 ! vector opt. 260 263 ! Upstream in the x-direction for the tracer 261 zfc(ji,jj,jk) = pt b(ji,jj-1,jk,jn)264 zfc(ji,jj,jk) = pt_lev1(ji,jj-1,jk,jn) 262 265 ! Downstream in the x-direction for the tracer 263 zfd(ji,jj,jk) = pt b(ji,jj+1,jk,jn)266 zfd(ji,jj,jk) = pt_lev1(ji,jj+1,jk,jn) 264 267 END DO 265 268 END DO … … 275 278 DO jj = 2, jpjm1 276 279 DO ji = fs_2, fs_jpim1 ! vector opt. 277 zdir = 0.5 + SIGN( 0.5, pv n(ji,jj,jk) ) ! if pun> 0 : zdir = 1 otherwise zdir = 0280 zdir = 0.5 + SIGN( 0.5, pv(ji,jj,jk) ) ! if pu > 0 : zdir = 1 otherwise zdir = 0 278 281 zfu(ji,jj,jk) = zdir * zfc(ji,jj,jk ) + ( 1. - zdir ) * zfd(ji,jj+1,jk) ! FU in the x-direction for T 279 282 END DO … … 284 287 DO jj = 2, jpjm1 285 288 DO ji = fs_2, fs_jpim1 ! vector opt. 286 zdir = 0.5 + SIGN( 0.5, pv n(ji,jj,jk) ) ! if pun> 0 : zdir = 1 otherwise zdir = 0287 zdx = ( zdir * e2t(ji,jj) + ( 1. - zdir ) * e2t(ji,jj+1) ) * e1v(ji,jj) * e3v _n(ji,jj,jk)288 zwy(ji,jj,jk) = ABS( pv n(ji,jj,jk) ) * p2dt / zdx ! (0<zc_cfl<1 : Courant number on x-direction)289 zfc(ji,jj,jk) = zdir * pt b(ji,jj ,jk,jn) + ( 1. - zdir ) * ptb(ji,jj+1,jk,jn) ! FC in the x-direction for T290 zfd(ji,jj,jk) = zdir * pt b(ji,jj+1,jk,jn) + ( 1. - zdir ) * ptb(ji,jj ,jk,jn) ! FD in the x-direction for T289 zdir = 0.5 + SIGN( 0.5, pv(ji,jj,jk) ) ! if pu > 0 : zdir = 1 otherwise zdir = 0 290 zdx = ( zdir * e2t(ji,jj) + ( 1. - zdir ) * e2t(ji,jj+1) ) * e1v(ji,jj) * e3v(ji,jj,jk,ktlev) 291 zwy(ji,jj,jk) = ABS( pv(ji,jj,jk) ) * p2dt / zdx ! (0<zc_cfl<1 : Courant number on x-direction) 292 zfc(ji,jj,jk) = zdir * pt_lev1(ji,jj ,jk,jn) + ( 1. - zdir ) * pt_lev1(ji,jj+1,jk,jn) ! FC in the x-direction for T 293 zfd(ji,jj,jk) = zdir * pt_lev1(ji,jj+1,jk,jn) + ( 1. - zdir ) * pt_lev1(ji,jj ,jk,jn) ! FD in the x-direction for T 291 294 END DO 292 295 END DO … … 314 317 DO jj = 2, jpjm1 315 318 DO ji = fs_2, fs_jpim1 ! vector opt. 316 zdir = 0.5 + SIGN( 0.5, pv n(ji,jj,jk) ) ! if pun> 0 : zdir = 1 otherwise zdir = 0319 zdir = 0.5 + SIGN( 0.5, pv(ji,jj,jk) ) ! if pu > 0 : zdir = 1 otherwise zdir = 0 317 320 !--- If the second ustream point is a land point 318 321 !--- the flux is computed by the 1st order UPWIND scheme 319 322 zmsk = zdir * zfu(ji,jj,jk) + ( 1. - zdir ) * zfu(ji,jj+1,jk) 320 323 zwy(ji,jj,jk) = zmsk * zwy(ji,jj,jk) + ( 1. - zmsk ) * zfc(ji,jj,jk) 321 zwy(ji,jj,jk) = zwy(ji,jj,jk) * pv n(ji,jj,jk)324 zwy(ji,jj,jk) = zwy(ji,jj,jk) * pv(ji,jj,jk) 322 325 END DO 323 326 END DO … … 330 333 DO jj = 2, jpjm1 331 334 DO ji = fs_2, fs_jpim1 ! vector opt. 332 zbtr = r1_e1e2t(ji,jj) / e3t _n(ji,jj,jk)335 zbtr = r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev) 333 336 ! horizontal advective trends 334 337 ztra = - zbtr * ( zwy(ji,jj,jk) - zwy(ji,jj-1,jk) ) 335 338 !--- add it to the general tracer trends 336 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztra339 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + ztra 337 340 END DO 338 341 END DO 339 342 END DO 340 343 ! ! trend diagnostics 341 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pv n, ptn(:,:,:,jn) )344 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pv, pt_lev2(:,:,:,jn) ) 342 345 ! ! "Poleward" heat and salt transports (contribution of upstream fluxes) 343 346 IF( l_ptr ) CALL dia_ptr_hst( jn, 'adv', zwy(:,:,:) ) … … 348 351 349 352 350 SUBROUTINE tra_adv_cen2_k( kt, cdtype, pwn, &351 & pt n, pta, kjpt )353 SUBROUTINE tra_adv_cen2_k( kt, ktlev, cdtype, pwn, & 354 & pt_lev2, pt_rhs, kjpt ) 352 355 !!---------------------------------------------------------------------- 353 356 !! 354 357 !!---------------------------------------------------------------------- 355 358 INTEGER , INTENT(in ) :: kt ! ocean time-step index 359 INTEGER , INTENT(in ) :: ktlev ! time level index for source terms 356 360 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 357 361 INTEGER , INTENT(in ) :: kjpt ! number of tracers 358 362 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pwn ! vertical velocity 359 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt n! before and now tracer fields360 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt a! tracer trend363 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt_lev2 ! before and now tracer fields 364 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 361 365 ! 362 366 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 374 378 DO jj = 2, jpjm1 375 379 DO ji = fs_2, fs_jpim1 ! vector opt. 376 zwz(ji,jj,jk) = 0.5 * pwn(ji,jj,jk) * ( pt n(ji,jj,jk-1,jn) + ptn(ji,jj,jk,jn) ) * wmask(ji,jj,jk)380 zwz(ji,jj,jk) = 0.5 * pwn(ji,jj,jk) * ( pt_lev2(ji,jj,jk-1,jn) + pt_lev2(ji,jj,jk,jn) ) * wmask(ji,jj,jk) 377 381 END DO 378 382 END DO … … 382 386 DO jj = 1, jpj 383 387 DO ji = 1, jpi 384 zwz(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * pt n(ji,jj,mikt(ji,jj),jn) ! linear free surface388 zwz(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * pt_lev2(ji,jj,mikt(ji,jj),jn) ! linear free surface 385 389 END DO 386 390 END DO 387 391 ELSE ! no ocean cavities (only ocean surface) 388 zwz(:,:,1) = pwn(:,:,1) * pt n(:,:,1,jn)392 zwz(:,:,1) = pwn(:,:,1) * pt_lev2(:,:,1,jn) 389 393 ENDIF 390 394 ENDIF … … 393 397 DO jj = 2, jpjm1 394 398 DO ji = fs_2, fs_jpim1 ! vector opt. 395 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - ( zwz(ji,jj,jk) - zwz(ji,jj,jk+1) ) &396 & * r1_e1e2t(ji,jj) / e3t _n(ji,jj,jk)399 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) - ( zwz(ji,jj,jk) - zwz(ji,jj,jk+1) ) & 400 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev) 397 401 END DO 398 402 END DO 399 403 END DO 400 404 ! ! Send trends for diagnostic 401 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pwn, pt n(:,:,:,jn) )405 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pwn, pt_lev2(:,:,:,jn) ) 402 406 ! 403 407 END DO -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_ubs.F90
r10425 r10802 46 46 CONTAINS 47 47 48 SUBROUTINE tra_adv_ubs( kt, kit000, cdtype, p2dt, pun, pvn, pwn, &49 & pt b, ptn, pta, kjpt, kn_ubs_v )48 SUBROUTINE tra_adv_ubs( kt, kit000, ktlev, cdtype, p2dt, pu, pv, pwn, & 49 & pt_lev1, pt_lev2, pt_rhs, kjpt, kn_ubs_v ) 50 50 !!---------------------------------------------------------------------- 51 51 !! *** ROUTINE tra_adv_ubs *** … … 58 58 !! It is only used in the horizontal direction. 59 59 !! For example the i-component of the advective fluxes are given by : 60 !! ! e2u e3u u n ( mi(Tn) - zltu(i ) ) if un(i) >= 060 !! ! e2u e3u uu ( mi(Tn) - zltu(i ) ,ktlev) if uu(i,ktlev) >= 0 61 61 !! ztu = ! or 62 !! ! e2u e3u u n ( mi(Tn) - zltu(i+1) ) if un(i) < 062 !! ! e2u e3u uu ( mi(Tn) - zltu(i+1) ,ktlev) if uu(i,ktlev) < 0 63 63 !! where zltu is the second derivative of the before temperature field: 64 64 !! zltu = 1/e3t di[ e2u e3u / e1u di[Tb] ] … … 77 77 !! scheme (kn_ubs_v=4). 78 78 !! 79 !! ** Action : - update pt awith the now advective tracer trends79 !! ** Action : - update pt_rhs with the now advective tracer trends 80 80 !! - send trends to trdtra module for further diagnostcs (l_trdtra=T) 81 81 !! - htr_adv, str_adv : poleward advective heat and salt transport (ln_diaptr=T) … … 86 86 INTEGER , INTENT(in ) :: kt ! ocean time-step index 87 87 INTEGER , INTENT(in ) :: kit000 ! first time step index 88 INTEGER , INTENT(in ) :: ktlev ! time level index for source terms 88 89 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 89 90 INTEGER , INTENT(in ) :: kjpt ! number of tracers 90 91 INTEGER , INTENT(in ) :: kn_ubs_v ! number of tracers 91 92 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 92 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu n, pvn, pwn ! 3 ocean transport components93 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt b, ptn! before and now tracer fields94 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt a! tracer trend93 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu, pv, pwn ! 3 ocean transport components 94 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt_lev1, pt_lev2 ! before and now tracer fields 95 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 95 96 ! 96 97 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 126 127 DO jj = 1, jpjm1 ! First derivative (masked gradient) 127 128 DO ji = 1, fs_jpim1 ! vector opt. 128 zeeu = e2_e1u(ji,jj) * e3u _n(ji,jj,jk) * umask(ji,jj,jk)129 zeev = e1_e2v(ji,jj) * e3v _n(ji,jj,jk) * vmask(ji,jj,jk)130 ztu(ji,jj,jk) = zeeu * ( pt b(ji+1,jj ,jk,jn) - ptb(ji,jj,jk,jn) )131 ztv(ji,jj,jk) = zeev * ( pt b(ji ,jj+1,jk,jn) - ptb(ji,jj,jk,jn) )129 zeeu = e2_e1u(ji,jj) * e3u(ji,jj,jk,ktlev) * umask(ji,jj,jk) 130 zeev = e1_e2v(ji,jj) * e3v(ji,jj,jk,ktlev) * vmask(ji,jj,jk) 131 ztu(ji,jj,jk) = zeeu * ( pt_lev1(ji+1,jj ,jk,jn) - pt_lev1(ji,jj,jk,jn) ) 132 ztv(ji,jj,jk) = zeev * ( pt_lev1(ji ,jj+1,jk,jn) - pt_lev1(ji,jj,jk,jn) ) 132 133 END DO 133 134 END DO 134 135 DO jj = 2, jpjm1 ! Second derivative (divergence) 135 136 DO ji = fs_2, fs_jpim1 ! vector opt. 136 zcoef = 1._wp / ( 6._wp * e3t _n(ji,jj,jk) )137 zcoef = 1._wp / ( 6._wp * e3t(ji,jj,jk,ktlev) ) 137 138 zltu(ji,jj,jk) = ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) ) * zcoef 138 139 zltv(ji,jj,jk) = ( ztv(ji,jj,jk) - ztv(ji,jj-1,jk) ) * zcoef … … 146 147 DO jj = 1, jpjm1 147 148 DO ji = 1, fs_jpim1 ! vector opt. 148 zfp_ui = pu n(ji,jj,jk) + ABS( pun(ji,jj,jk) ) ! upstream transport (x2)149 zfm_ui = pu n(ji,jj,jk) - ABS( pun(ji,jj,jk) )150 zfp_vj = pv n(ji,jj,jk) + ABS( pvn(ji,jj,jk) )151 zfm_vj = pv n(ji,jj,jk) - ABS( pvn(ji,jj,jk) )149 zfp_ui = pu(ji,jj,jk) + ABS( pu(ji,jj,jk) ) ! upstream transport (x2) 150 zfm_ui = pu(ji,jj,jk) - ABS( pu(ji,jj,jk) ) 151 zfp_vj = pv(ji,jj,jk) + ABS( pv(ji,jj,jk) ) 152 zfm_vj = pv(ji,jj,jk) - ABS( pv(ji,jj,jk) ) 152 153 ! ! 2nd order centered advective fluxes (x2) 153 zcenut = pu n(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji+1,jj ,jk,jn) )154 zcenvt = pv n(ji,jj,jk) * ( ptn(ji,jj,jk,jn) + ptn(ji ,jj+1,jk,jn) )154 zcenut = pu(ji,jj,jk) * ( pt_lev2(ji,jj,jk,jn) + pt_lev2(ji+1,jj ,jk,jn) ) 155 zcenvt = pv(ji,jj,jk) * ( pt_lev2(ji,jj,jk,jn) + pt_lev2(ji ,jj+1,jk,jn) ) 155 156 ! ! UBS advective fluxes 156 157 ztu(ji,jj,jk) = 0.5 * ( zcenut - zfp_ui * zltu(ji,jj,jk) - zfm_ui * zltu(ji+1,jj,jk) ) … … 160 161 END DO 161 162 ! 162 zltu(:,:,:) = pt a(:,:,:,jn) ! store the initial trends before its update163 zltu(:,:,:) = pt_rhs(:,:,:,jn) ! store the initial trends before its update 163 164 ! 164 165 DO jk = 1, jpkm1 !== add the horizontal advective trend ==! 165 166 DO jj = 2, jpjm1 166 167 DO ji = fs_2, fs_jpim1 ! vector opt. 167 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) &168 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) & 168 169 & - ( ztu(ji,jj,jk) - ztu(ji-1,jj ,jk) & 169 & + ztv(ji,jj,jk) - ztv(ji ,jj-1,jk) ) * r1_e1e2t(ji,jj) / e3t _n(ji,jj,jk)170 & + ztv(ji,jj,jk) - ztv(ji ,jj-1,jk) ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev) 170 171 END DO 171 172 END DO … … 173 174 END DO 174 175 ! 175 zltu(:,:,:) = pt a(:,:,:,jn) - zltu(:,:,:) ! Horizontal advective trend used in vertical 2nd order FCT case176 zltu(:,:,:) = pt_rhs(:,:,:,jn) - zltu(:,:,:) ! Horizontal advective trend used in vertical 2nd order FCT case 176 177 ! ! and/or in trend diagnostic (l_trd=T) 177 178 ! 178 179 IF( l_trd ) THEN ! trend diagnostics 179 CALL trd_tra( kt, cdtype, jn, jptra_xad, ztu, pu n, ptn(:,:,:,jn) )180 CALL trd_tra( kt, cdtype, jn, jptra_yad, ztv, pv n, ptn(:,:,:,jn) )180 CALL trd_tra( kt, cdtype, jn, jptra_xad, ztu, pu, pt_lev2(:,:,:,jn) ) 181 CALL trd_tra( kt, cdtype, jn, jptra_yad, ztv, pv, pt_lev2(:,:,:,jn) ) 181 182 END IF 182 183 ! … … 193 194 CASE( 2 ) ! 2nd order FCT 194 195 ! 195 IF( l_trd ) zltv(:,:,:) = pt a(:,:,:,jn) ! store ptaif trend diag.196 IF( l_trd ) zltv(:,:,:) = pt_rhs(:,:,:,jn) ! store pt_rhs if trend diag. 196 197 ! 197 198 ! !* upstream advection with initial mass fluxes & intermediate update ==! … … 201 202 zfp_wk = pwn(ji,jj,jk) + ABS( pwn(ji,jj,jk) ) 202 203 zfm_wk = pwn(ji,jj,jk) - ABS( pwn(ji,jj,jk) ) 203 ztw(ji,jj,jk) = 0.5_wp * ( zfp_wk * pt b(ji,jj,jk,jn) + zfm_wk * ptb(ji,jj,jk-1,jn) ) * wmask(ji,jj,jk)204 ztw(ji,jj,jk) = 0.5_wp * ( zfp_wk * pt_lev1(ji,jj,jk,jn) + zfm_wk * pt_lev1(ji,jj,jk-1,jn) ) * wmask(ji,jj,jk) 204 205 END DO 205 206 END DO … … 209 210 DO jj = 1, jpj 210 211 DO ji = 1, jpi 211 ztw(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * pt b(ji,jj,mikt(ji,jj),jn) ! linear free surface212 ztw(ji,jj, mikt(ji,jj) ) = pwn(ji,jj,mikt(ji,jj)) * pt_lev1(ji,jj,mikt(ji,jj),jn) ! linear free surface 212 213 END DO 213 214 END DO 214 215 ELSE ! no cavities: only at the ocean surface 215 ztw(:,:,1) = pwn(:,:,1) * pt b(:,:,1,jn)216 ztw(:,:,1) = pwn(:,:,1) * pt_lev1(:,:,1,jn) 216 217 ENDIF 217 218 ENDIF … … 220 221 DO jj = 2, jpjm1 221 222 DO ji = fs_2, fs_jpim1 ! vector opt. 222 ztak = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t _n(ji,jj,jk)223 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ztak224 zti(ji,jj,jk) = ( pt b(ji,jj,jk,jn) + p2dt * ( ztak + zltu(ji,jj,jk) ) ) * tmask(ji,jj,jk)223 ztak = - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev) 224 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + ztak 225 zti(ji,jj,jk) = ( pt_lev1(ji,jj,jk,jn) + p2dt * ( ztak + zltu(ji,jj,jk) ) ) * tmask(ji,jj,jk) 225 226 END DO 226 227 END DO … … 232 233 DO jj = 1, jpj 233 234 DO ji = 1, jpi 234 ztw(ji,jj,jk) = ( 0.5_wp * pwn(ji,jj,jk) * ( pt n(ji,jj,jk,jn) + ptn(ji,jj,jk-1,jn) ) &235 ztw(ji,jj,jk) = ( 0.5_wp * pwn(ji,jj,jk) * ( pt_lev2(ji,jj,jk,jn) + pt_lev2(ji,jj,jk-1,jn) ) & 235 236 & - ztw(ji,jj,jk) ) * wmask(ji,jj,jk) 236 237 END DO … … 240 241 IF( ln_linssh ) ztw(:,:, 1 ) = 0._wp ! only ocean surface as interior zwz values have been w-masked 241 242 ! 242 CALL nonosc_z( pt b(:,:,:,jn), ztw, zti, p2dt) ! monotonicity algorithm243 CALL nonosc_z( pt_lev1(:,:,:,jn), ztw, zti, p2dt, e3t(:,:,:,ktlev) ) ! monotonicity algorithm 243 244 ! 244 245 CASE( 4 ) ! 4th order COMPACT 245 CALL interp_4th_cpt( pt n(:,:,:,jn) , ztw ) ! 4th order compact interpolation of T at w-point246 CALL interp_4th_cpt( pt_lev2(:,:,:,jn) , ztw ) ! 4th order compact interpolation of T at w-point 246 247 DO jk = 2, jpkm1 247 248 DO jj = 2, jpjm1 … … 251 252 END DO 252 253 END DO 253 IF( ln_linssh ) ztw(:,:, 1 ) = pwn(:,:,1) * pt n(:,:,1,jn) !!gm ISF & 4th COMPACT doesn't work254 IF( ln_linssh ) ztw(:,:, 1 ) = pwn(:,:,1) * pt_lev2(:,:,1,jn) !!gm ISF & 4th COMPACT doesn't work 254 255 ! 255 256 END SELECT … … 258 259 DO jj = 2, jpjm1 259 260 DO ji = fs_2, fs_jpim1 ! vector opt. 260 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t_n(ji,jj,jk)261 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) - ( ztw(ji,jj,jk) - ztw(ji,jj,jk+1) ) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev) 261 262 END DO 262 263 END DO … … 264 265 ! 265 266 IF( l_trd ) THEN ! vertical advective trend diagnostics 266 DO jk = 1, jpkm1 ! (compute -w.dk[ptn]= -dk[w.ptn] + pt n.dk[w])267 DO jk = 1, jpkm1 ! (compute -w.dk[ptn]= -dk[w.ptn] + pt_lev2.dk[w]) 267 268 DO jj = 2, jpjm1 268 269 DO ji = fs_2, fs_jpim1 ! vector opt. 269 zltv(ji,jj,jk) = pt a(ji,jj,jk,jn) - zltv(ji,jj,jk) &270 & + pt n(ji,jj,jk,jn) * ( pwn(ji,jj,jk) - pwn(ji,jj,jk+1) ) &271 & * r1_e1e2t(ji,jj) / e3t _n(ji,jj,jk)270 zltv(ji,jj,jk) = pt_rhs(ji,jj,jk,jn) - zltv(ji,jj,jk) & 271 & + pt_lev2(ji,jj,jk,jn) * ( pwn(ji,jj,jk) - pwn(ji,jj,jk+1) ) & 272 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev) 272 273 END DO 273 274 END DO … … 281 282 282 283 283 SUBROUTINE nonosc_z( pbef, pcc, paft, p2dt )284 SUBROUTINE nonosc_z( pbef, pcc, paft, p2dt, pe3t ) 284 285 !!--------------------------------------------------------------------- 285 286 !! *** ROUTINE nonosc_z *** … … 296 297 REAL(wp), INTENT(in ) :: p2dt ! tracer time-step 297 298 REAL(wp), DIMENSION (jpi,jpj,jpk) :: pbef ! before field 299 REAL(wp), INTENT(in ), DIMENSION (jpi,jpj,jpk) :: pe3t ! now cell thickness field 298 300 REAL(wp), INTENT(inout), DIMENSION (jpi,jpj,jpk) :: paft ! after field 299 301 REAL(wp), INTENT(inout), DIMENSION (jpi,jpj,jpk) :: pcc ! monotonic flux in the k direction … … 352 354 zneg = MAX( 0., pcc(ji ,jj ,jk ) ) - MIN( 0., pcc(ji ,jj ,jk+1) ) 353 355 ! up & down beta terms 354 zbt = e1e2t(ji,jj) * e3t_n(ji,jj,jk) / p2dt356 zbt = e1e2t(ji,jj) * pe3t(ji,jj,jk) / p2dt 355 357 zbetup(ji,jj,jk) = ( zbetup(ji,jj,jk) - paft(ji,jj,jk) ) / (zpos+zrtrn) * zbt 356 358 zbetdo(ji,jj,jk) = ( paft(ji,jj,jk) - zbetdo(ji,jj,jk) ) / (zneg+zrtrn) * zbt -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/oce.F90
r10789 r10802 23 23 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: wi !: vertical vel. (adaptive-implicit) [m/s] 24 24 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hdivn !: horizontal divergence [s-1] 25 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,: ) :: tsb , tsn , tsa!: 4D T-S fields [Celsius,psu]25 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:,:), TARGET :: ts !: 4D T-S fields [Celsius,psu] 26 26 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: rab_b, rab_n !: thermal/haline expansion coef. [Celsius-1,psu-1] 27 27 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: rn2b , rn2 !: brunt-vaisala frequency**2 [s-2] … … 72 72 REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:,:) :: vb , vn , va !: j-horizontal velocity [m/s] 73 73 REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:,:) :: wn !: k-vertical velocity [m/s] 74 REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:,:,:) :: tsb , tsn , tsa !: 4D T-S fields [Celsius,psu] 75 !! TEMPORARY POINTERS FOR DEVELOPMENT ONLY 74 76 75 77 !!---------------------------------------------------------------------- … … 90 92 ALLOCATE( uu (jpi,jpj,jpk, jpt) , vv (jpi,jpj,jpk,jpt) , & 91 93 & ww (jpi,jpj,jpk) , hdivn(jpi,jpj,jpk) , & 92 & ts b (jpi,jpj,jpk,jpts) , tsn (jpi,jpj,jpk,jpts) , tsa(jpi,jpj,jpk,jpts) ,&94 & ts (jpi,jpj,jpk,jpts,jpt) , & 93 95 & rab_b(jpi,jpj,jpk,jpts) , rab_n(jpi,jpj,jpk,jpts) , & 94 96 & rn2b (jpi,jpj,jpk) , rn2 (jpi,jpj,jpk) , & -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/step.F90
r10799 r10802 247 247 & CALL Agrif_Sponge_tra ! tracers sponge 248 248 #endif 249 CALL tra_adv ( kstp ) ! horizontal & vertical advection249 CALL tra_adv ( kstp, Nm1, Nnn, Np1, ts(:,:,:,:,Nrhs) ) ! horizontal & vertical advection 250 250 IF( ln_zdfosm ) CALL tra_osm ( kstp ) ! OSMOSIS non-local tracer fluxes 251 251 IF( lrst_oce .AND. ln_zdfosm ) & … … 344 344 wn => ww(:,:,:) 345 345 346 tsb => ts(:,:,:,:,Nm1); tsn => ts(:,:,:,:,Nnn); tsa => ts(:,:,:,:,Np1) 347 346 348 e3t_b => e3t(:,:,:,Nm1); e3t_n => e3t(:,:,:,Nnn); e3t_a => e3t(:,:,:,Np1) 347 349 e3u_b => e3u(:,:,:,Nm1); e3u_n => e3u(:,:,:,Nnn); e3u_a => e3u(:,:,:,Np1) -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcadv.F90
r10068 r10802 68 68 CONTAINS 69 69 70 SUBROUTINE trc_adv( kt )70 SUBROUTINE trc_adv( kt, ktlev1, ktlev2, ktlev3 ) 71 71 !!---------------------------------------------------------------------- 72 72 !! *** ROUTINE trc_adv *** … … 77 77 !!---------------------------------------------------------------------- 78 78 INTEGER, INTENT(in) :: kt ! ocean time-step index 79 INTEGER, INTENT(in) :: ktlev1, ktlev2, ktlev3 ! time level indices for source terms 79 80 ! 80 81 INTEGER :: jk ! dummy loop index … … 123 124 ! 124 125 CASE ( np_CEN ) ! Centered : 2nd / 4th order 125 CALL tra_adv_cen( kt, nittrc000, 'TRC', zun, zvn, zwn , trn, tra, jptra, nn_cen_h, nn_cen_v )126 CALL tra_adv_cen( kt, nittrc000, ktlev2, 'TRC', zun, zvn, zwn , trn, tra, jptra, nn_cen_h, nn_cen_v ) 126 127 CASE ( np_FCT ) ! FCT : 2nd / 4th order 127 CALL tra_adv_fct( kt, nittrc000, 'TRC', r2dttrc, zun, zvn, zwn, trb, trn, tra, jptra, nn_fct_h, nn_fct_v )128 CALL tra_adv_fct( kt, nittrc000, ktlev1, ktlev2, ktlev3, 'TRC', r2dttrc, zun, zvn, zwn, trb, trn, tra, jptra, nn_fct_h, nn_fct_v ) 128 129 CASE ( np_MUS ) ! MUSCL 129 CALL tra_adv_mus( kt, nittrc000, 'TRC', r2dttrc, zun, zvn, zwn, trb, tra, jptra , ln_mus_ups )130 CALL tra_adv_mus( kt, nittrc000, ktlev2, 'TRC', r2dttrc, zun, zvn, zwn, trb, tra, jptra , ln_mus_ups ) 130 131 CASE ( np_UBS ) ! UBS 131 CALL tra_adv_ubs( kt, nittrc000, 'TRC', r2dttrc, zun, zvn, zwn, trb, trn, tra, jptra , nn_ubs_v )132 CALL tra_adv_ubs( kt, nittrc000, ktlev2, 'TRC', r2dttrc, zun, zvn, zwn, trb, trn, tra, jptra , nn_ubs_v ) 132 133 CASE ( np_QCK ) ! QUICKEST 133 CALL tra_adv_qck( kt, nittrc000, 'TRC', r2dttrc, zun, zvn, zwn, trb, trn, tra, jptra )134 CALL tra_adv_qck( kt, nittrc000, ktlev2, 'TRC', r2dttrc, zun, zvn, zwn, trb, trn, tra, jptra ) 134 135 ! 135 136 END SELECT -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trctrp.F90
r10068 r10802 64 64 IF( ln_trcdmp ) CALL trc_dmp ( kt ) ! internal damping trends 65 65 IF( ln_bdy ) CALL trc_bdy_dmp( kt ) ! BDY damping trends 66 CALL trc_adv ( kt ) ! horizontal & vertical advection66 CALL trc_adv ( kt, Nm1, Nnn, Np1 ) ! horizontal & vertical advection 67 67 ! ! Partial top/bottom cell: GRADh( trb ) 68 68 IF( ln_zps ) THEN
Note: See TracChangeset
for help on using the changeset viewer.
