Changeset 10806
- Timestamp:
- 2019-03-27T17:55:22+01:00 (6 years ago)
- Location:
- NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps
- Files:
-
- 24 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/cfgs/SHARED/field_def_nemo-oce.xml
r10645 r10806 62 62 <field id="mldr10_1max" long_name="Max of Mixed Layer Depth (dsigma = 0.01 wrt 10m)" field_ref="mldr10_1" operation="maximum" /> 63 63 <field id="mldr10_1min" long_name="Min of Mixed Layer Depth (dsigma = 0.01 wrt 10m)" field_ref="mldr10_1" operation="minimum" /> 64 <field id="mldzint_1" long_name="Mixed Layer Depth interpolated" standard_name="ocean_mixed_layer_thickness" unit="m" /> 65 <field id="mldzint_2" long_name="Mixed Layer Depth interpolated" standard_name="ocean_mixed_layer_thickness" unit="m" /> 66 <field id="mldzint_3" long_name="Mixed Layer Depth interpolated" standard_name="ocean_mixed_layer_thickness" unit="m" /> 67 <field id="mldzint_4" long_name="Mixed Layer Depth interpolated" standard_name="ocean_mixed_layer_thickness" unit="m" /> 68 <field id="mldzint_5" long_name="Mixed Layer Depth interpolated" standard_name="ocean_mixed_layer_thickness" unit="m" /> 69 <field id="mldhtc_1" long_name="Mixed Layer Depth integrated heat content" standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content" unit="J/m2" /> 70 <field id="mldhtc_2" long_name="Mixed Layer Depth integrated heat content" standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content" unit="J/m2" /> 71 <field id="mldhtc_3" long_name="Mixed Layer Depth integrated heat content" standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content" unit="J/m2" /> 72 <field id="mldhtc_4" long_name="Mixed Layer Depth integrated heat content" standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content" unit="J/m2" /> 73 <field id="mldhtc_5" long_name="Mixed Layer Depth integrated heat content" standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content" unit="J/m2" /> 64 74 <field id="heatc" long_name="Heat content vertically integrated" standard_name="integral_of_sea_water_potential_temperature_wrt_depth_expressed_as_heat_content" unit="J/m2" /> 65 75 <field id="saltc" long_name="Salt content vertically integrated" unit="1e-3*kg/m2" /> … … 77 87 <!-- variables available with MLE --> 78 88 <field id="Lf_NHpf" long_name="MLE: Lf = N H / f" unit="m" /> 79 80 <!-- next variables available with ln_zad_Aimp=.true. -->81 <field id="Courant" long_name="Courant number" unit="#" grid_ref="grid_T_3D" />82 <field id="wimp" long_name="Implicit vertical velocity" unit="m/s" grid_ref="grid_T_3D" />83 <field id="wexp" long_name="Explicit vertical velocity" unit="m/s" grid_ref="grid_T_3D" />84 <field id="wi_cff" long_name="Fraction of implicit vertical velocity" unit="#" grid_ref="grid_T_3D" />85 89 86 90 <!-- next variables available with key_diahth --> … … 351 355 <field id="uoce" long_name="ocean current along i-axis" standard_name="sea_water_x_velocity" unit="m/s" grid_ref="grid_U_3D" /> 352 356 <field id="uoce_e3u" long_name="ocean current along i-axis (thickness weighted)" unit="m/s" grid_ref="grid_U_3D" > uoce * e3u </field> 357 <field id="uoce2_e3u" long_name="ocean current along i-axis squared (thickness weighted)" unit="m3/s2" grid_ref="grid_U_3D" > uoce * uoce * e3u </field> 353 358 <field id="ssu" long_name="ocean surface current along i-axis" unit="m/s" /> 354 359 <field id="sbu" long_name="ocean bottom current along i-axis" unit="m/s" /> … … 405 410 <field id="voce" long_name="ocean current along j-axis" standard_name="sea_water_y_velocity" unit="m/s" grid_ref="grid_V_3D" /> 406 411 <field id="voce_e3v" long_name="ocean current along j-axis (thickness weighted)" unit="m/s" grid_ref="grid_V_3D" > voce * e3v </field> 412 <field id="voce2_e3v" long_name="ocean current along j-axis squared (thickness weighted)" unit="m3/s2" grid_ref="grid_V_3D" > voce * voce * e3v </field> 407 413 <field id="ssv" long_name="ocean surface current along j-axis" unit="m/s" /> 408 414 <field id="sbv" long_name="ocean bottom current along j-axis" unit="m/s" /> … … 975 981 976 982 <field_group id="25h_grid_W" grid_ref="grid_W_3D" operation="instant"> 977 <field id="vo vecrtz25h" name="k current 25h mean" unit="m/s" />983 <field id="vomecrtz25h" name="k current 25h mean" unit="m/s" /> 978 984 <field id="avt25h" name="vertical diffusivity25h mean" unit="m2/s" /> 979 985 <field id="avm25h" name="vertical viscosity 25h mean" unit="m2/s" /> -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DOM/dom_oce.F90
r10787 r10806 154 154 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:,:) :: gde3w 155 155 156 ! ! reference depths of water column 157 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ht_0 !: t-depth [m] 158 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hu_0 !: u-depth [m] 159 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hv_0 !: v-depth [m] 160 ! ! time-dependent depths of water column 161 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:,:) :: hu, hv, r1_hu, r1_hv 162 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, TARGET, DIMENSION(:,:) :: ht 156 ! ! ref. ! before ! now ! after ! 157 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ht_0 , ht_n !: t-depth [m] 158 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hu_0 , hu_b , hu_n , hu_a !: u-depth [m] 159 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hv_0 , hv_b , hv_n , hv_a !: v-depth [m] 160 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: r1_hu_b , r1_hu_n , r1_hu_a !: inverse of u-depth [1/m] 161 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: r1_hv_b , r1_hv_n , r1_hv_a !: inverse of v-depth [1/m] 163 162 164 163 !! TEMPORARY POINTERS FOR DEVELOPMENT ONLY … … 173 172 REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:,:) :: gdepw_b , gdepw_n !: w- depth [m] 174 173 REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:,:) :: gde3w_n !: w- depth (sum of e3w) [m] 175 REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:) :: ht_n !: t-depth [m]176 REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:) :: hu_b , hu_n , hu_a !: u-depth [m]177 REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:) :: hv_b , hv_n , hv_a !: v-depth [m]178 REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:) :: r1_hu_b , r1_hu_n , r1_hu_a !: inverse of u-depth [1/m]179 REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:) :: r1_hv_b , r1_hv_n , r1_hv_a !: inverse of v-depth [1/m]180 174 !! TEMPORARY POINTERS FOR DEVELOPMENT ONLY 181 175 … … 295 289 & e3f (jpi,jpj,jpk), STAT=ierr(5) ) 296 290 ! 297 ALLOCATE( ht_0(jpi,jpj) , hu_0(jpi,jpj) , hv_0(jpi,jpj) , & 298 & ht (jpi,jpj) , hu (jpi,jpj,jpt) , hv (jpi,jpj,jpt) , & 299 & r1_hu (jpi,jpj,jpt) , r1_hv (jpi,jpj,jpt) , STAT=ierr(6) ) 291 ALLOCATE( ht_0(jpi,jpj) , hu_0(jpi,jpj) , hv_0(jpi,jpj) , & 292 & hu_b(jpi,jpj) , hv_b(jpi,jpj) , r1_hu_b(jpi,jpj) , r1_hv_b(jpi,jpj) , & 293 & ht_n(jpi,jpj) , hu_n(jpi,jpj) , hv_n(jpi,jpj) , r1_hu_n(jpi,jpj) , r1_hv_n(jpi,jpj) , & 294 & hu_a(jpi,jpj) , hv_a(jpi,jpj) , r1_hu_a(jpi,jpj) , r1_hv_a(jpi,jpj) , STAT=ierr(6) ) 300 295 ! 301 296 ! -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynldf.F90
r10068 r10806 43 43 CONTAINS 44 44 45 SUBROUTINE dyn_ldf( kt )45 SUBROUTINE dyn_ldf( kt, ktlev1, ktlev2, pu_rhs, pv_rhs ) 46 46 !!---------------------------------------------------------------------- 47 47 !! *** ROUTINE dyn_ldf *** … … 49 49 !! ** Purpose : compute the lateral ocean dynamics physics. 50 50 !!---------------------------------------------------------------------- 51 INTEGER, INTENT(in) :: kt ! ocean time-step index 51 INTEGER, INTENT(in) :: kt ! ocean time-step index 52 INTEGER, INTENT(in) :: ktlev1, ktlev2 ! time level index for source terms 53 REAL(wp), INTENT(inout), DIMENSION(jpi,jpj,jpk) :: pu_rhs, pv_rhs ! momentum trends 52 54 ! 53 55 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdu, ztrdv … … 58 60 IF( l_trddyn ) THEN ! temporary save of momentum trends 59 61 ALLOCATE( ztrdu(jpi,jpj,jpk) , ztrdv(jpi,jpj,jpk) ) 60 ztrdu(:,:,:) = ua(:,:,:)61 ztrdv(:,:,:) = va(:,:,:)62 ztrdu(:,:,:) = pu_rhs(:,:,:) 63 ztrdv(:,:,:) = pv_rhs(:,:,:) 62 64 ENDIF 63 65 64 66 SELECT CASE ( nldf_dyn ) ! compute lateral mixing trend and add it to the general trend 65 67 ! 66 CASE ( np_lap ) ; CALL dyn_ldf_lap( kt, ub, vb, ua, va, 1 ) ! iso-level laplacian68 CASE ( np_lap ) ; CALL dyn_ldf_lap( kt, ktlev1, ktlev2, uu(:,:,:,ktlev1), vv(:,:,:,ktlev1), pu_rhs, pv_rhs, 1 ) ! iso-level laplacian 67 69 CASE ( np_lap_i ) ; CALL dyn_ldf_iso( kt ) ! rotated laplacian 68 CASE ( np_blp ) ; CALL dyn_ldf_blp( kt, ub, vb, ua, va) ! iso-level bi-laplacian70 CASE ( np_blp ) ; CALL dyn_ldf_blp( kt, ktlev1, ktlev2, uu(:,:,:,ktlev1), vv(:,:,:,ktlev1), pu_rhs, pv_rhs ) ! iso-level bi-laplacian 69 71 ! 70 72 END SELECT 71 73 72 74 IF( l_trddyn ) THEN ! save the horizontal diffusive trends for further diagnostics 73 ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:)74 ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)75 ztrdu(:,:,:) = pu_rhs(:,:,:) - ztrdu(:,:,:) 76 ztrdv(:,:,:) = pv_rhs(:,:,:) - ztrdv(:,:,:) 75 77 CALL trd_dyn( ztrdu, ztrdv, jpdyn_ldf, kt ) 76 78 DEALLOCATE ( ztrdu , ztrdv ) -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynldf_lap_blp.F90
r10425 r10806 35 35 CONTAINS 36 36 37 SUBROUTINE dyn_ldf_lap( kt, pub, pvb, pua, pva, kpass )37 SUBROUTINE dyn_ldf_lap( kt, ktlev1, ktlev2, pu, pv, pu_rhs, pva_rhs, kpass ) 38 38 !!---------------------------------------------------------------------- 39 39 !! *** ROUTINE dyn_ldf_lap *** … … 45 45 !! writen as : grad_h( ahmt div_h(U )) - curl_h( ahmf curl_z(U) ) 46 46 !! 47 !! ** Action : - pu a, pva increased by the harmonic operator applied on pub, pvb.47 !! ** Action : - pu_rhs, pva_rhs increased by the harmonic operator applied on pu, pv. 48 48 !!---------------------------------------------------------------------- 49 INTEGER , INTENT(in ) :: kt ! ocean time-step index 49 INTEGER , INTENT(in ) :: kt ! ocean time-step index 50 INTEGER , INTENT(in ) :: ktlev1, ktlev2 ! time level index for scale factors 50 51 INTEGER , INTENT(in ) :: kpass ! =1/2 first or second passage 51 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pu b, pvb! before velocity [m/s]52 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu a, pva! velocity trend [m/s2]52 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pu, pv ! before velocity [m/s] 53 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu_rhs, pva_rhs ! velocity trend [m/s2] 53 54 ! 54 55 INTEGER :: ji, jj, jk ! dummy loop indices … … 76 77 !!gm open question here : e3f at before or now ? probably now... 77 78 !!gm note that ahmf has already been multiplied by fmask 78 zcur(ji-1,jj-1) = ahmf(ji-1,jj-1,jk) * e3f _n(ji-1,jj-1,jk) * r1_e1e2f(ji-1,jj-1) &79 & * ( e2v(ji ,jj-1) * pv b(ji ,jj-1,jk) - e2v(ji-1,jj-1) * pvb(ji-1,jj-1,jk) &80 & - e1u(ji-1,jj ) * pu b(ji-1,jj ,jk) + e1u(ji-1,jj-1) * pub(ji-1,jj-1,jk) )79 zcur(ji-1,jj-1) = ahmf(ji-1,jj-1,jk) * e3f(ji-1,jj-1,jk) * r1_e1e2f(ji-1,jj-1) & 80 & * ( e2v(ji ,jj-1) * pv(ji ,jj-1,jk) - e2v(ji-1,jj-1) * pv(ji-1,jj-1,jk) & 81 & - e1u(ji-1,jj ) * pu(ji-1,jj ,jk) + e1u(ji-1,jj-1) * pu(ji-1,jj-1,jk) ) 81 82 ! ! ahm * div (computed from 2 to jpi/jpj) 82 83 !!gm note that ahmt has already been multiplied by tmask 83 zdiv(ji,jj) = ahmt(ji,jj,jk) * r1_e1e2t(ji,jj) / e3t _b(ji,jj,jk) &84 & * ( e2u(ji,jj)*e3u _b(ji,jj,jk) * pub(ji,jj,jk) - e2u(ji-1,jj)*e3u_b(ji-1,jj,jk) * pub(ji-1,jj,jk) &85 & + e1v(ji,jj)*e3v _b(ji,jj,jk) * pvb(ji,jj,jk) - e1v(ji,jj-1)*e3v_b(ji,jj-1,jk) * pvb(ji,jj-1,jk) )84 zdiv(ji,jj) = ahmt(ji,jj,jk) * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev1) & 85 & * ( e2u(ji,jj)*e3u(ji,jj,jk,ktlev1) * pu(ji,jj,jk) - e2u(ji-1,jj)*e3u(ji-1,jj,jk,ktlev1) * pu(ji-1,jj,jk) & 86 & + e1v(ji,jj)*e3v(ji,jj,jk,ktlev1) * pv(ji,jj,jk) - e1v(ji,jj-1)*e3v(ji,jj-1,jk,ktlev1) * pv(ji,jj-1,jk) ) 86 87 END DO 87 88 END DO … … 89 90 DO jj = 2, jpjm1 ! - curl( curl) + grad( div ) 90 91 DO ji = fs_2, fs_jpim1 ! vector opt. 91 pu a(ji,jj,jk) = pua(ji,jj,jk) + zsign * ( &92 & - ( zcur(ji ,jj) - zcur(ji,jj-1) ) * r1_e2u(ji,jj) / e3u _n(ji,jj,jk) &92 pu_rhs(ji,jj,jk) = pu_rhs(ji,jj,jk) + zsign * ( & 93 & - ( zcur(ji ,jj) - zcur(ji,jj-1) ) * r1_e2u(ji,jj) / e3u(ji,jj,jk,ktlev2) & 93 94 & + ( zdiv(ji+1,jj) - zdiv(ji,jj ) ) * r1_e1u(ji,jj) ) 94 95 ! 95 pva (ji,jj,jk) = pva(ji,jj,jk) + zsign * ( &96 & ( zcur(ji,jj ) - zcur(ji-1,jj) ) * r1_e1v(ji,jj) / e3v _n(ji,jj,jk) &96 pva_rhs(ji,jj,jk) = pva_rhs(ji,jj,jk) + zsign * ( & 97 & ( zcur(ji,jj ) - zcur(ji-1,jj) ) * r1_e1v(ji,jj) / e3v(ji,jj,jk,ktlev2) & 97 98 & + ( zdiv(ji,jj+1) - zdiv(ji ,jj) ) * r1_e2v(ji,jj) ) 98 99 END DO … … 105 106 106 107 107 SUBROUTINE dyn_ldf_blp( kt, pub, pvb, pua, pva)108 SUBROUTINE dyn_ldf_blp( kt, ktlev1, ktlev2, pu, pv, pu_rhs, pva_rhs ) 108 109 !!---------------------------------------------------------------------- 109 110 !! *** ROUTINE dyn_ldf_blp *** … … 116 117 !! It is computed by two successive calls to dyn_ldf_lap routine 117 118 !! 118 !! ** Action : pt aupdated with the before rotated bilaplacian diffusion119 !! ** Action : pt_rhs updated with the before rotated bilaplacian diffusion 119 120 !!---------------------------------------------------------------------- 120 121 INTEGER , INTENT(in ) :: kt ! ocean time-step index 121 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pub, pvb ! before velocity fields 122 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pua, pva ! momentum trend 122 INTEGER , INTENT(in ) :: ktlev1, ktlev2 ! time level index for scale factors 123 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pu, pv ! before velocity fields 124 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(inout) :: pu_rhs, pva_rhs ! momentum trend 123 125 ! 124 126 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zulap, zvlap ! laplacian at u- and v-point … … 134 136 zvlap(:,:,:) = 0._wp 135 137 ! 136 CALL dyn_ldf_lap( kt, pub, pvb, zulap, zvlap, 1 ) ! rotated laplacian applied to ptb(output in zlap)138 CALL dyn_ldf_lap( kt, ktlev1, ktlev2, pu, pv, zulap, zvlap, 1 ) ! rotated laplacian applied to pt (output in zlap) 137 139 ! 138 140 CALL lbc_lnk_multi( 'dynldf_lap_blp', zulap, 'U', -1., zvlap, 'V', -1. ) ! Lateral boundary conditions 139 141 ! 140 CALL dyn_ldf_lap( kt, zulap, zvlap, pua, pva, 2 ) ! rotated laplacian applied to zlap (output in pta)142 CALL dyn_ldf_lap( kt, ktlev1, ktlev2, zulap, zvlap, pu_rhs, pva_rhs, 2 ) ! rotated laplacian applied to zlap (output in pt_rhs) 141 143 ! 142 144 END SUBROUTINE dyn_ldf_blp -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/DYN/dynvor.F90
r10802 r10806 119 119 ALLOCATE( ztrdu(jpi,jpj,jpk), ztrdv(jpi,jpj,jpk) ) 120 120 ! 121 ztrdu(:,:,:) = ua(:,:,:) !* planetary vorticity trend (including Stokes-Coriolis force)122 ztrdv(:,:,:) = va(:,:,:)121 ztrdu(:,:,:) = pu_rhs(:,:,:) !* planetary vorticity trend (including Stokes-Coriolis force) 122 ztrdv(:,:,:) = pv_rhs(:,:,:) 123 123 SELECT CASE( nvor_scheme ) 124 124 CASE( np_ENS ) ; CALL vor_ens( kt, ktlev, ncor, uu(:,:,:,ktlev) , vv(:,:,:,ktlev) , pu_rhs, pv_rhs ) ! enstrophy conserving scheme … … 133 133 IF( ln_stcor ) CALL vor_een( kt, ktlev, ncor, usd, vsd, pu_rhs, pv_rhs ) ! add the Stokes-Coriolis trend 134 134 END SELECT 135 ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:)136 ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)135 ztrdu(:,:,:) = pu_rhs(:,:,:) - ztrdu(:,:,:) 136 ztrdv(:,:,:) = pv_rhs(:,:,:) - ztrdv(:,:,:) 137 137 CALL trd_dyn( ztrdu, ztrdv, jpdyn_pvo, kt ) 138 138 ! 139 139 IF( n_dynadv /= np_LIN_dyn ) THEN !* relative vorticity or metric trend (only in non-linear case) 140 ztrdu(:,:,:) = ua(:,:,:)141 ztrdv(:,:,:) = va(:,:,:)140 ztrdu(:,:,:) = pu_rhs(:,:,:) 141 ztrdv(:,:,:) = pv_rhs(:,:,:) 142 142 SELECT CASE( nvor_scheme ) 143 143 CASE( np_ENT ) ; CALL vor_enT( kt, ktlev, nrvm, uu(:,:,:,ktlev) , vv(:,:,:,ktlev) , pu_rhs, pv_rhs ) ! energy conserving scheme (T-pts) … … 147 147 CASE( np_EEN ) ; CALL vor_een( kt, ktlev, nrvm, uu(:,:,:,ktlev) , vv(:,:,:,ktlev) , pu_rhs, pv_rhs ) ! energy & enstrophy scheme 148 148 END SELECT 149 ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:)150 ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:)149 ztrdu(:,:,:) = pu_rhs(:,:,:) - ztrdu(:,:,:) 150 ztrdv(:,:,:) = pv_rhs(:,:,:) - ztrdv(:,:,:) 151 151 CALL trd_dyn( ztrdu, ztrdv, jpdyn_rvo, kt ) 152 152 ENDIF … … 410 410 411 411 IF( ln_sco ) THEN 412 zwz(:,:) = zwz(:,:) / e3f _n(:,:,jk)412 zwz(:,:) = zwz(:,:) / e3f(:,:,jk) 413 413 zwx(:,:) = e2u(:,:) * e3u(:,:,jk,ktlev) * pu(:,:,jk) 414 414 zwy(:,:) = e1v(:,:) * e3v(:,:,jk,ktlev) * pv(:,:,jk) … … 518 518 ! 519 519 IF( ln_sco ) THEN !== horizontal fluxes ==! 520 zwz(:,:) = zwz(:,:) / e3f _n(:,:,jk)520 zwz(:,:) = zwz(:,:) / e3f(:,:,jk) 521 521 zwx(:,:) = e2u(:,:) * e3u(:,:,jk,ktlev) * pu(:,:,jk) 522 522 zwy(:,:) = e1v(:,:) * e3v(:,:,jk,ktlev) * pv(:,:,jk) -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv.F90
r10802 r10806 75 75 CONTAINS 76 76 77 SUBROUTINE tra_adv( kt, ktlev1, ktlev2, ktlev3, pts_rhs )77 SUBROUTINE tra_adv( kt, ktlev1, ktlev2, ktlev3, kt2lev, pts_rhs ) 78 78 !!---------------------------------------------------------------------- 79 79 !! *** ROUTINE tra_adv *** … … 84 84 !!---------------------------------------------------------------------- 85 85 INTEGER, INTENT(in) :: kt ! ocean time-step index 86 INTEGER, INTENT(in) :: ktlev1, ktlev2, ktlev3 ! time level indices for source terms 86 INTEGER, INTENT(in) :: ktlev1, ktlev2, ktlev3 ! time level indices for 3-time-level source terms 87 INTEGER, INTENT(in) :: kt2lev ! time level index for 2-time-level source terms 87 88 REAL(wp), INTENT( inout), DIMENSION(jpi,jpj,jpk,jpts) :: pts_rhs ! temperature and salinity trends 88 89 ! … … 153 154 & ts(:,:,:,:,ktlev1), ts(:,:,:,:,ktlev2), pts_rhs, jpts, nn_fct_h, nn_fct_v ) 154 155 CASE ( np_MUS ) ! MUSCL 155 CALL tra_adv_mus ( kt, nit000, ktlev2, 'TRA', r2dt, zun, zvn, zwn, ts(:,:,:,:,ktlev1), pts_rhs, jpts , ln_mus_ups )156 CALL tra_adv_mus ( kt, nit000, ktlev2, kt2lev, 'TRA', r2dt, zun, zvn, zwn, ts(:,:,:,:,ktlev1), pts_rhs, jpts , ln_mus_ups ) 156 157 CASE ( np_UBS ) ! UBS 157 158 CALL tra_adv_ubs ( kt, nit000, ktlev2, 'TRA', r2dt, zun, zvn, zwn, ts(:,:,:,:,ktlev1), ts(:,:,:,:,ktlev2), pts_rhs, jpts , nn_ubs_v ) -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_cen.F90
r10802 r10806 44 44 CONTAINS 45 45 46 SUBROUTINE tra_adv_cen( kt, kit000, ktlev, cdtype, pu, pv, pw n, &46 SUBROUTINE tra_adv_cen( kt, kit000, ktlev, cdtype, pu, pv, pw, & 47 47 & pt, pt_rhs, kjpt, kn_cen_h, kn_cen_v ) 48 48 !!---------------------------------------------------------------------- … … 70 70 INTEGER , INTENT(in ) :: kn_cen_h ! =2/4 (2nd or 4th order scheme) 71 71 INTEGER , INTENT(in ) :: kn_cen_v ! =2/4 (2nd or 4th order scheme) 72 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu, pv, pw n! 3 ocean velocity components72 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu, pv, pw ! 3 ocean velocity components 73 73 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt ! now tracer fields 74 74 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend … … 151 151 DO jj = 2, jpjm1 152 152 DO ji = fs_2, fs_jpim1 ! vector opt. 153 zwz(ji,jj,jk) = 0.5 * pw n(ji,jj,jk) * ( pt(ji,jj,jk,jn) + pt(ji,jj,jk-1,jn) ) * wmask(ji,jj,jk)153 zwz(ji,jj,jk) = 0.5 * pw(ji,jj,jk) * ( pt(ji,jj,jk,jn) + pt(ji,jj,jk-1,jn) ) * wmask(ji,jj,jk) 154 154 END DO 155 155 END DO … … 161 161 DO jj = 2, jpjm1 162 162 DO ji = fs_2, fs_jpim1 163 zwz(ji,jj,jk) = pw n(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk)163 zwz(ji,jj,jk) = pw(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk) 164 164 END DO 165 165 END DO … … 172 172 DO jj = 1, jpj 173 173 DO ji = 1, jpi 174 zwz(ji,jj, mikt(ji,jj) ) = pw n(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn)174 zwz(ji,jj, mikt(ji,jj) ) = pw(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn) 175 175 END DO 176 176 END DO 177 177 ELSE ! no ice-shelf cavities (only ocean surface) 178 zwz(:,:,1) = pw n(:,:,1) * pt(:,:,1,jn)178 zwz(:,:,1) = pw(:,:,1) * pt(:,:,1,jn) 179 179 ENDIF 180 180 ENDIF … … 194 194 CALL trd_tra( kt, cdtype, jn, jptra_xad, zwx, pu, pt(:,:,:,jn) ) 195 195 CALL trd_tra( kt, cdtype, jn, jptra_yad, zwy, pv, pt(:,:,:,jn) ) 196 CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pw n, pt(:,:,:,jn) )196 CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pw, pt(:,:,:,jn) ) 197 197 END IF 198 198 ! ! "Poleward" heat and salt transports -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_fct.F90
r10802 r10806 52 52 CONTAINS 53 53 54 SUBROUTINE tra_adv_fct( kt, kit000, ktlev1, ktlev2, ktlev3, cdtype, p2dt, pu, pv, pw n, &54 SUBROUTINE tra_adv_fct( kt, kit000, ktlev1, ktlev2, ktlev3, cdtype, p2dt, pu, pv, pw, & 55 55 & pt_lev1, pt_lev2, pt_rhs, kjpt, kn_fct_h, kn_fct_v ) 56 56 !!---------------------------------------------------------------------- … … 77 77 INTEGER , INTENT(in ) :: kn_fct_v ! order of the FCT scheme (=2 or 4) 78 78 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 79 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu, pv, pw n! 3 ocean velocity components79 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu, pv, pw ! 3 ocean velocity components 80 80 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt_lev1, pt_lev2 ! before and now tracer fields 81 81 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend … … 139 139 DO jj = 1, jpj 140 140 DO ji = 1, jpi 141 zfp_wk = pw n(ji,jj,jk) + ABS( pwn(ji,jj,jk) )142 zfm_wk = pw n(ji,jj,jk) - ABS( pwn(ji,jj,jk) )141 zfp_wk = pw(ji,jj,jk) + ABS( pw(ji,jj,jk) ) 142 zfm_wk = pw(ji,jj,jk) - ABS( pw(ji,jj,jk) ) 143 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) 144 144 END DO … … 149 149 DO jj = 1, jpj 150 150 DO ji = 1, jpi 151 zwz(ji,jj, mikt(ji,jj) ) = pw n(ji,jj,mikt(ji,jj)) * pt_lev1(ji,jj,mikt(ji,jj),jn) ! linear free surface151 zwz(ji,jj, mikt(ji,jj) ) = pw(ji,jj,mikt(ji,jj)) * pt_lev1(ji,jj,mikt(ji,jj),jn) ! linear free surface 152 152 END DO 153 153 END DO 154 154 ELSE ! no cavities: only at the ocean surface 155 zwz(:,:,1) = pw n(:,:,1) * pt_lev1(:,:,1,jn)155 zwz(:,:,1) = pw(:,:,1) * pt_lev1(:,:,1,jn) 156 156 ENDIF 157 157 ENDIF … … 258 258 DO jj = 2, jpjm1 259 259 DO ji = fs_2, fs_jpim1 260 zwz(ji,jj,jk) = ( pw n(ji,jj,jk) * 0.5_wp * ( pt_lev2(ji,jj,jk,jn) + pt_lev2(ji,jj,jk-1,jn) ) &260 zwz(ji,jj,jk) = ( pw(ji,jj,jk) * 0.5_wp * ( pt_lev2(ji,jj,jk,jn) + pt_lev2(ji,jj,jk-1,jn) ) & 261 261 & - zwz(ji,jj,jk) ) * wmask(ji,jj,jk) 262 262 END DO … … 269 269 DO jj = 2, jpjm1 270 270 DO ji = fs_2, fs_jpim1 271 zwz(ji,jj,jk) = ( pw n(ji,jj,jk) * ztw(ji,jj,jk) - zwz(ji,jj,jk) ) * wmask(ji,jj,jk)271 zwz(ji,jj,jk) = ( pw(ji,jj,jk) * ztw(ji,jj,jk) - zwz(ji,jj,jk) ) * wmask(ji,jj,jk) 272 272 END DO 273 273 END DO … … 306 306 CALL trd_tra( kt, cdtype, jn, jptra_xad, ztrdx, pu, pt_lev2(:,:,:,jn) ) 307 307 CALL trd_tra( kt, cdtype, jn, jptra_yad, ztrdy, pv, pt_lev2(:,:,:,jn) ) 308 CALL trd_tra( kt, cdtype, jn, jptra_zad, ztrdz, pw n, pt_lev2(:,:,:,jn) )308 CALL trd_tra( kt, cdtype, jn, jptra_zad, ztrdz, pw, pt_lev2(:,:,:,jn) ) 309 309 ENDIF 310 310 ! ! heat/salt transport -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_mus.F90
r10802 r10806 54 54 CONTAINS 55 55 56 SUBROUTINE tra_adv_mus( kt, kit000, ktlev, cdtype, p2dt, pu, pv, pwn, &57 & pt, pt_rhs, kjpt, ld_msc_ups )56 SUBROUTINE tra_adv_mus( kt, kit000, ktlev, kt2lev, cdtype, p2dt, pu, pv, pw, & 57 & pt, pt_rhs, kjpt, ld_msc_ups ) 58 58 !!---------------------------------------------------------------------- 59 59 !! *** ROUTINE tra_adv_mus *** … … 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 INTEGER , INTENT(in ) :: ktlev ! time level index for 3-time-level source terms 78 INTEGER , INTENT(in ) :: kt2lev ! time level index for 2-time-level source terms 78 79 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 79 80 INTEGER , INTENT(in ) :: kjpt ! number of tracers 80 81 LOGICAL , INTENT(in ) :: ld_msc_ups ! use upstream scheme within muscl 81 82 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 82 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu, pv, pw n! 3 ocean velocity components83 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu, pv, pw ! 3 ocean velocity components 83 84 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt ! before tracer field 84 85 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend … … 240 241 DO jj = 2, jpjm1 241 242 DO ji = fs_2, fs_jpim1 ! vector opt. 242 z0w = SIGN( 0.5, pw n(ji,jj,jk+1) )243 z0w = SIGN( 0.5, pw(ji,jj,jk+1) ) 243 244 zalpha = 0.5 + z0w 244 zw = z0w - 0.5 * pw n(ji,jj,jk+1) * p2dt * r1_e1e2t(ji,jj) / e3w_n(ji,jj,jk+1)245 zw = z0w - 0.5 * pw(ji,jj,jk+1) * p2dt * r1_e1e2t(ji,jj) / e3w(ji,jj,jk+1,kt2lev) 245 246 zzwx = pt(ji,jj,jk+1,jn) + xind(ji,jj,jk) * zw * zslpx(ji,jj,jk+1) 246 247 zzwy = pt(ji,jj,jk ,jn) + xind(ji,jj,jk) * zw * zslpx(ji,jj,jk ) 247 zwx(ji,jj,jk+1) = pw n(ji,jj,jk+1) * ( zalpha * zzwx + (1.-zalpha) * zzwy ) * wmask(ji,jj,jk)248 zwx(ji,jj,jk+1) = pw(ji,jj,jk+1) * ( zalpha * zzwx + (1.-zalpha) * zzwy ) * wmask(ji,jj,jk) 248 249 END DO 249 250 END DO … … 253 254 DO jj = 1, jpj 254 255 DO ji = 1, jpi 255 zwx(ji,jj, mikt(ji,jj) ) = pw n(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn)256 zwx(ji,jj, mikt(ji,jj) ) = pw(ji,jj,mikt(ji,jj)) * pt(ji,jj,mikt(ji,jj),jn) 256 257 END DO 257 258 END DO 258 259 ELSE ! no cavities: only at the ocean surface 259 zwx(:,:,1) = pw n(:,:,1) * pt(:,:,1,jn)260 zwx(:,:,1) = pw(:,:,1) * pt(:,:,1,jn) 260 261 ENDIF 261 262 ENDIF … … 269 270 END DO 270 271 ! ! send trends for diagnostic 271 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_zad, zwx, pw n, pt(:,:,:,jn) )272 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_zad, zwx, pw, pt(:,:,:,jn) ) 272 273 ! 273 274 END DO ! end of tracer loop -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_qck.F90
r10802 r10806 47 47 CONTAINS 48 48 49 SUBROUTINE tra_adv_qck ( kt, kit000, ktlev, cdtype, p2dt, pu, pv, pw n, &49 SUBROUTINE tra_adv_qck ( kt, kit000, ktlev, cdtype, p2dt, pu, pv, pw, & 50 50 & pt_lev1, pt_lev2, pt_rhs, kjpt ) 51 51 !!---------------------------------------------------------------------- … … 90 90 INTEGER , INTENT(in ) :: kjpt ! number of tracers 91 91 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 92 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu, pv, pw n! 3 ocean velocity components92 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu, pv, pw ! 3 ocean velocity components 93 93 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt_lev1, pt_lev2 ! before and now tracer fields 94 94 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend … … 113 113 114 114 ! ! vertical fluxes are computed with the 2nd order centered scheme 115 CALL tra_adv_cen2_k( kt, ktlev, cdtype, pw n, pt_lev2, pt_rhs, kjpt )115 CALL tra_adv_cen2_k( kt, ktlev, cdtype, pw, pt_lev2, pt_rhs, kjpt ) 116 116 ! 117 117 END SUBROUTINE tra_adv_qck … … 351 351 352 352 353 SUBROUTINE tra_adv_cen2_k( kt, ktlev, cdtype, pw n, &353 SUBROUTINE tra_adv_cen2_k( kt, ktlev, cdtype, pw, & 354 354 & pt_lev2, pt_rhs, kjpt ) 355 355 !!---------------------------------------------------------------------- … … 360 360 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 361 361 INTEGER , INTENT(in ) :: kjpt ! number of tracers 362 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pw n! vertical velocity362 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pw ! vertical velocity 363 363 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt_lev2 ! before and now tracer fields 364 364 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend … … 378 378 DO jj = 2, jpjm1 379 379 DO ji = fs_2, fs_jpim1 ! vector opt. 380 zwz(ji,jj,jk) = 0.5 * pw n(ji,jj,jk) * ( pt_lev2(ji,jj,jk-1,jn) + pt_lev2(ji,jj,jk,jn) ) * wmask(ji,jj,jk)380 zwz(ji,jj,jk) = 0.5 * pw(ji,jj,jk) * ( pt_lev2(ji,jj,jk-1,jn) + pt_lev2(ji,jj,jk,jn) ) * wmask(ji,jj,jk) 381 381 END DO 382 382 END DO … … 386 386 DO jj = 1, jpj 387 387 DO ji = 1, jpi 388 zwz(ji,jj, mikt(ji,jj) ) = pw n(ji,jj,mikt(ji,jj)) * pt_lev2(ji,jj,mikt(ji,jj),jn) ! linear free surface388 zwz(ji,jj, mikt(ji,jj) ) = pw(ji,jj,mikt(ji,jj)) * pt_lev2(ji,jj,mikt(ji,jj),jn) ! linear free surface 389 389 END DO 390 390 END DO 391 391 ELSE ! no ocean cavities (only ocean surface) 392 zwz(:,:,1) = pw n(:,:,1) * pt_lev2(:,:,1,jn)392 zwz(:,:,1) = pw(:,:,1) * pt_lev2(:,:,1,jn) 393 393 ENDIF 394 394 ENDIF … … 403 403 END DO 404 404 ! ! Send trends for diagnostic 405 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pw n, pt_lev2(:,:,:,jn) )405 IF( l_trd ) CALL trd_tra( kt, cdtype, jn, jptra_zad, zwz, pw, pt_lev2(:,:,:,jn) ) 406 406 ! 407 407 END DO -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traadv_ubs.F90
r10802 r10806 46 46 CONTAINS 47 47 48 SUBROUTINE tra_adv_ubs( kt, kit000, ktlev, cdtype, p2dt, pu, pv, pw n, &48 SUBROUTINE tra_adv_ubs( kt, kit000, ktlev, cdtype, p2dt, pu, pv, pw, & 49 49 & pt_lev1, pt_lev2, pt_rhs, kjpt, kn_ubs_v ) 50 50 !!---------------------------------------------------------------------- … … 91 91 INTEGER , INTENT(in ) :: kn_ubs_v ! number of tracers 92 92 REAL(wp) , INTENT(in ) :: p2dt ! tracer time-step 93 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu, pv, pw n! 3 ocean transport components93 REAL(wp), DIMENSION(jpi,jpj,jpk ), INTENT(in ) :: pu, pv, pw ! 3 ocean transport components 94 94 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt_lev1, pt_lev2 ! before and now tracer fields 95 95 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend … … 200 200 DO jj = 1, jpj 201 201 DO ji = 1, jpi 202 zfp_wk = pw n(ji,jj,jk) + ABS( pwn(ji,jj,jk) )203 zfm_wk = pw n(ji,jj,jk) - ABS( pwn(ji,jj,jk) )202 zfp_wk = pw(ji,jj,jk) + ABS( pw(ji,jj,jk) ) 203 zfm_wk = pw(ji,jj,jk) - ABS( pw(ji,jj,jk) ) 204 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) 205 205 END DO … … 210 210 DO jj = 1, jpj 211 211 DO ji = 1, jpi 212 ztw(ji,jj, mikt(ji,jj) ) = pw n(ji,jj,mikt(ji,jj)) * pt_lev1(ji,jj,mikt(ji,jj),jn) ! linear free surface212 ztw(ji,jj, mikt(ji,jj) ) = pw(ji,jj,mikt(ji,jj)) * pt_lev1(ji,jj,mikt(ji,jj),jn) ! linear free surface 213 213 END DO 214 214 END DO 215 215 ELSE ! no cavities: only at the ocean surface 216 ztw(:,:,1) = pw n(:,:,1) * pt_lev1(:,:,1,jn)216 ztw(:,:,1) = pw(:,:,1) * pt_lev1(:,:,1,jn) 217 217 ENDIF 218 218 ENDIF … … 233 233 DO jj = 1, jpj 234 234 DO ji = 1, jpi 235 ztw(ji,jj,jk) = ( 0.5_wp * pw n(ji,jj,jk) * ( pt_lev2(ji,jj,jk,jn) + pt_lev2(ji,jj,jk-1,jn) ) &235 ztw(ji,jj,jk) = ( 0.5_wp * pw(ji,jj,jk) * ( pt_lev2(ji,jj,jk,jn) + pt_lev2(ji,jj,jk-1,jn) ) & 236 236 & - ztw(ji,jj,jk) ) * wmask(ji,jj,jk) 237 237 END DO … … 248 248 DO jj = 2, jpjm1 249 249 DO ji = fs_2, fs_jpim1 250 ztw(ji,jj,jk) = pw n(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk)251 END DO 252 END DO 253 END DO 254 IF( ln_linssh ) ztw(:,:, 1 ) = pw n(:,:,1) * pt_lev2(:,:,1,jn) !!gm ISF & 4th COMPACT doesn't work250 ztw(ji,jj,jk) = pw(ji,jj,jk) * ztw(ji,jj,jk) * wmask(ji,jj,jk) 251 END DO 252 END DO 253 END DO 254 IF( ln_linssh ) ztw(:,:, 1 ) = pw(:,:,1) * pt_lev2(:,:,1,jn) !!gm ISF & 4th COMPACT doesn't work 255 255 ! 256 256 END SELECT … … 269 269 DO ji = fs_2, fs_jpim1 ! vector opt. 270 270 zltv(ji,jj,jk) = pt_rhs(ji,jj,jk,jn) - zltv(ji,jj,jk) & 271 & + pt_lev2(ji,jj,jk,jn) * ( pw n(ji,jj,jk) - pwn(ji,jj,jk+1) ) &271 & + pt_lev2(ji,jj,jk,jn) * ( pw(ji,jj,jk) - pw(ji,jj,jk+1) ) & 272 272 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev) 273 273 END DO -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trabbc.F90
r10425 r10806 51 51 CONTAINS 52 52 53 SUBROUTINE tra_bbc( kt )53 SUBROUTINE tra_bbc( kt, ktlev, pts_rhs ) 54 54 !!---------------------------------------------------------------------- 55 55 !! *** ROUTINE tra_bbc *** … … 74 74 !!---------------------------------------------------------------------- 75 75 INTEGER, INTENT(in) :: kt ! ocean time-step index 76 INTEGER, INTENT(in) :: ktlev ! time level index for source terms 77 REAL(wp), INTENT( inout), DIMENSION(jpi,jpj,jpk,jpts) :: pts_rhs ! temperature and salinity trends 76 78 ! 77 79 INTEGER :: ji, jj ! dummy loop indices … … 83 85 IF( l_trdtra ) THEN ! Save the input temperature trend 84 86 ALLOCATE( ztrdt(jpi,jpj,jpk) ) 85 ztrdt(:,:,:) = tsa(:,:,:,jp_tem)87 ztrdt(:,:,:) = pts_rhs(:,:,:,jp_tem) 86 88 ENDIF 87 89 ! ! Add the geothermal trend on temperature 88 90 DO jj = 2, jpjm1 89 91 DO ji = 2, jpim1 90 tsa(ji,jj,mbkt(ji,jj),jp_tem) = tsa(ji,jj,mbkt(ji,jj),jp_tem) + qgh_trd0(ji,jj) / e3t_n(ji,jj,mbkt(ji,jj))92 pts_rhs(ji,jj,mbkt(ji,jj),jp_tem) = pts_rhs(ji,jj,mbkt(ji,jj),jp_tem) + qgh_trd0(ji,jj) / e3t(ji,jj,mbkt(ji,jj),ktlev) 91 93 END DO 92 94 END DO 93 95 ! 94 CALL lbc_lnk( 'trabbc', tsa(:,:,:,jp_tem) , 'T', 1. )96 CALL lbc_lnk( 'trabbc', pts_rhs(:,:,:,jp_tem) , 'T', 1. ) 95 97 ! 96 98 IF( l_trdtra ) THEN ! Send the trend for diagnostics 97 ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)99 ztrdt(:,:,:) = pts_rhs(:,:,:,jp_tem) - ztrdt(:,:,:) 98 100 CALL trd_tra( kt, 'TRA', jp_tem, jptra_bbc, ztrdt ) 99 101 DEALLOCATE( ztrdt ) -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trabbl.F90
r10425 r10806 89 89 90 90 91 SUBROUTINE tra_bbl( kt )91 SUBROUTINE tra_bbl( kt, ktlev1, ktlev2, kt2lev, pts_rhs ) 92 92 !!---------------------------------------------------------------------- 93 93 !! *** ROUTINE bbl *** … … 101 101 !! is added to the general tracer trend 102 102 !!---------------------------------------------------------------------- 103 INTEGER, INTENT( in ) :: kt ! ocean time-step 103 INTEGER, INTENT( in ) :: kt ! ocean time-step 104 INTEGER, INTENT( in ) :: ktlev1, ktlev2 ! time level indices for 3-time-level source terms 105 INTEGER, INTENT( in ) :: kt2lev ! time level index for 2-time-level source terms 106 REAL(wp), INTENT( inout), DIMENSION(jpi,jpj,jpk,jpts) :: pts_rhs ! temperature and salinity trends 104 107 ! 105 108 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdt, ztrds … … 110 113 IF( l_trdtra ) THEN !* Save the T-S input trends 111 114 ALLOCATE( ztrdt(jpi,jpj,jpk) , ztrds(jpi,jpj,jpk) ) 112 ztrdt(:,:,:) = tsa(:,:,:,jp_tem)113 ztrds(:,:,:) = tsa(:,:,:,jp_sal)114 ENDIF 115 116 IF( l_bbl ) CALL bbl( kt, nit000, 'TRA' ) !* bbl coef. and transport (only if not already done in trcbbl)115 ztrdt(:,:,:) = pts_rhs(:,:,:,jp_tem) 116 ztrds(:,:,:) = pts_rhs(:,:,:,jp_sal) 117 ENDIF 118 119 IF( l_bbl ) CALL bbl( kt, nit000, ktlev1, ktlev2, kt2lev, 'TRA' ) !* bbl coef. and transport (only if not already done in trcbbl) 117 120 118 121 IF( nn_bbl_ldf == 1 ) THEN !* Diffusive bbl 119 122 ! 120 CALL tra_bbl_dif( ts b, tsa, jpts )123 CALL tra_bbl_dif( ts(:,:,:,:,ktlev1), e3t(:,:,:,ktlev2), pts_rhs, jpts ) 121 124 IF( ln_ctl ) & 122 125 CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' bbl_ldf - Ta: ', mask1=tmask, & … … 131 134 IF( nn_bbl_adv /= 0 ) THEN !* Advective bbl 132 135 ! 133 CALL tra_bbl_adv( ts b, tsa, jpts )136 CALL tra_bbl_adv( ts(:,:,:,:,ktlev1), e3t(:,:,:,ktlev2), pts_rhs, jpts ) 134 137 IF(ln_ctl) & 135 138 CALL prt_ctl( tab3d_1=tsa(:,:,:,jp_tem), clinfo1=' bbl_adv - Ta: ', mask1=tmask, & … … 143 146 144 147 IF( l_trdtra ) THEN ! send the trends for further diagnostics 145 ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)146 ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:)148 ztrdt(:,:,:) = pts_rhs(:,:,:,jp_tem) - ztrdt(:,:,:) 149 ztrds(:,:,:) = pts_rhs(:,:,:,jp_sal) - ztrds(:,:,:) 147 150 CALL trd_tra( kt, 'TRA', jp_tem, jptra_bbl, ztrdt ) 148 151 CALL trd_tra( kt, 'TRA', jp_sal, jptra_bbl, ztrds ) … … 155 158 156 159 157 SUBROUTINE tra_bbl_dif( pt b, pta, kjpt )160 SUBROUTINE tra_bbl_dif( pt, pe3t, pt_rhs, kjpt ) 158 161 !!---------------------------------------------------------------------- 159 162 !! *** ROUTINE tra_bbl_dif *** … … 171 174 !! convection is satified) 172 175 !! 173 !! ** Action : pt aincreased by the bbl diffusive trend176 !! ** Action : pt_rhs increased by the bbl diffusive trend 174 177 !! 175 178 !! References : Beckmann, A., and R. Doscher, 1997, J. Phys.Oceanogr., 581-591. … … 177 180 !!---------------------------------------------------------------------- 178 181 INTEGER , INTENT(in ) :: kjpt ! number of tracers 179 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before and now tracer fields 180 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! tracer trend 182 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt ! tracer fields 183 REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT(in ) :: pe3t ! thickness fields 184 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 181 185 ! 182 186 INTEGER :: ji, jj, jn ! dummy loop indices … … 191 195 DO ji = 1, jpi 192 196 ik = mbkt(ji,jj) ! bottom T-level index 193 zptb(ji,jj) = pt b(ji,jj,ik,jn) ! bottom before T and S197 zptb(ji,jj) = pt(ji,jj,ik,jn) ! bottom before T and S 194 198 END DO 195 199 END DO … … 198 202 DO ji = 2, jpim1 199 203 ik = mbkt(ji,jj) ! bottom T-level index 200 pt a(ji,jj,ik,jn) = pta(ji,jj,ik,jn) &204 pt_rhs(ji,jj,ik,jn) = pt_rhs(ji,jj,ik,jn) & 201 205 & + ( ahu_bbl(ji ,jj ) * ( zptb(ji+1,jj ) - zptb(ji ,jj ) ) & 202 206 & - ahu_bbl(ji-1,jj ) * ( zptb(ji ,jj ) - zptb(ji-1,jj ) ) & 203 207 & + ahv_bbl(ji ,jj ) * ( zptb(ji ,jj+1) - zptb(ji ,jj ) ) & 204 208 & - ahv_bbl(ji ,jj-1) * ( zptb(ji ,jj ) - zptb(ji ,jj-1) ) ) & 205 & * r1_e1e2t(ji,jj) / e3t_n(ji,jj,ik)209 & * r1_e1e2t(ji,jj) / pe3t(ji,jj,ik) 206 210 END DO 207 211 END DO … … 212 216 213 217 214 SUBROUTINE tra_bbl_adv( pt b, pta, kjpt )218 SUBROUTINE tra_bbl_adv( pt, pe3t, pt_rhs, kjpt ) 215 219 !!---------------------------------------------------------------------- 216 220 !! *** ROUTINE trc_bbl *** … … 228 232 !!---------------------------------------------------------------------- 229 233 INTEGER , INTENT(in ) :: kjpt ! number of tracers 230 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: ptb ! before and now tracer fields 231 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pta ! tracer trend 234 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt ! before and now tracer fields 235 REAL(wp), DIMENSION(jpi,jpj,jpk) , INTENT(in ) :: pe3t ! thickness fields 236 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 232 237 ! 233 238 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 250 255 ! 251 256 ! ! up -slope T-point (shelf bottom point) 252 zbtr = r1_e1e2t(iis,jj) / e3t_n(iis,jj,ikus)253 ztra = zu_bbl * ( pt b(iid,jj,ikus,jn) - ptb(iis,jj,ikus,jn) ) * zbtr254 pt a(iis,jj,ikus,jn) = pta(iis,jj,ikus,jn) + ztra257 zbtr = r1_e1e2t(iis,jj) / pe3t(iis,jj,ikus) 258 ztra = zu_bbl * ( pt(iid,jj,ikus,jn) - pt(iis,jj,ikus,jn) ) * zbtr 259 pt_rhs(iis,jj,ikus,jn) = pt_rhs(iis,jj,ikus,jn) + ztra 255 260 ! 256 261 DO jk = ikus, ikud-1 ! down-slope upper to down T-point (deep column) 257 zbtr = r1_e1e2t(iid,jj) / e3t_n(iid,jj,jk)258 ztra = zu_bbl * ( pt b(iid,jj,jk+1,jn) - ptb(iid,jj,jk,jn) ) * zbtr259 pt a(iid,jj,jk,jn) = pta(iid,jj,jk,jn) + ztra262 zbtr = r1_e1e2t(iid,jj) / pe3t(iid,jj,jk) 263 ztra = zu_bbl * ( pt(iid,jj,jk+1,jn) - pt(iid,jj,jk,jn) ) * zbtr 264 pt_rhs(iid,jj,jk,jn) = pt_rhs(iid,jj,jk,jn) + ztra 260 265 END DO 261 266 ! 262 zbtr = r1_e1e2t(iid,jj) / e3t_n(iid,jj,ikud)263 ztra = zu_bbl * ( pt b(iis,jj,ikus,jn) - ptb(iid,jj,ikud,jn) ) * zbtr264 pt a(iid,jj,ikud,jn) = pta(iid,jj,ikud,jn) + ztra267 zbtr = r1_e1e2t(iid,jj) / pe3t(iid,jj,ikud) 268 ztra = zu_bbl * ( pt(iis,jj,ikus,jn) - pt(iid,jj,ikud,jn) ) * zbtr 269 pt_rhs(iid,jj,ikud,jn) = pt_rhs(iid,jj,ikud,jn) + ztra 265 270 ENDIF 266 271 ! … … 272 277 ! 273 278 ! up -slope T-point (shelf bottom point) 274 zbtr = r1_e1e2t(ji,ijs) / e3t_n(ji,ijs,ikvs)275 ztra = zv_bbl * ( pt b(ji,ijd,ikvs,jn) - ptb(ji,ijs,ikvs,jn) ) * zbtr276 pt a(ji,ijs,ikvs,jn) = pta(ji,ijs,ikvs,jn) + ztra279 zbtr = r1_e1e2t(ji,ijs) / pe3t(ji,ijs,ikvs) 280 ztra = zv_bbl * ( pt(ji,ijd,ikvs,jn) - pt(ji,ijs,ikvs,jn) ) * zbtr 281 pt_rhs(ji,ijs,ikvs,jn) = pt_rhs(ji,ijs,ikvs,jn) + ztra 277 282 ! 278 283 DO jk = ikvs, ikvd-1 ! down-slope upper to down T-point (deep column) 279 zbtr = r1_e1e2t(ji,ijd) / e3t_n(ji,ijd,jk)280 ztra = zv_bbl * ( pt b(ji,ijd,jk+1,jn) - ptb(ji,ijd,jk,jn) ) * zbtr281 pt a(ji,ijd,jk,jn) = pta(ji,ijd,jk,jn) + ztra284 zbtr = r1_e1e2t(ji,ijd) / pe3t(ji,ijd,jk) 285 ztra = zv_bbl * ( pt(ji,ijd,jk+1,jn) - pt(ji,ijd,jk,jn) ) * zbtr 286 pt_rhs(ji,ijd,jk,jn) = pt_rhs(ji,ijd,jk,jn) + ztra 282 287 END DO 283 288 ! ! down-slope T-point (deep bottom point) 284 zbtr = r1_e1e2t(ji,ijd) / e3t_n(ji,ijd,ikvd)285 ztra = zv_bbl * ( pt b(ji,ijs,ikvs,jn) - ptb(ji,ijd,ikvd,jn) ) * zbtr286 pt a(ji,ijd,ikvd,jn) = pta(ji,ijd,ikvd,jn) + ztra289 zbtr = r1_e1e2t(ji,ijd) / pe3t(ji,ijd,ikvd) 290 ztra = zv_bbl * ( pt(ji,ijs,ikvs,jn) - pt(ji,ijd,ikvd,jn) ) * zbtr 291 pt_rhs(ji,ijd,ikvd,jn) = pt_rhs(ji,ijd,ikvd,jn) + ztra 287 292 ENDIF 288 293 END DO … … 295 300 296 301 297 SUBROUTINE bbl( kt, kit000, cdtype )302 SUBROUTINE bbl( kt, kit000, ktlev1, ktlev2, kt2lev, cdtype ) 298 303 !!---------------------------------------------------------------------- 299 304 !! *** ROUTINE bbl *** … … 323 328 INTEGER , INTENT(in ) :: kt ! ocean time-step index 324 329 INTEGER , INTENT(in ) :: kit000 ! first time step index 330 INTEGER , INTENT(in ) :: ktlev1, ktlev2 ! time level indices for 3-time-levelsource terms 331 INTEGER , INTENT(in ) :: kt2lev ! time level index for 2-time-level source terms 325 332 CHARACTER(len=3), INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 326 333 ! … … 344 351 DO ji = 1, jpi 345 352 ik = mbkt(ji,jj) ! bottom T-level index 346 zts (ji,jj,jp_tem) = ts b(ji,jj,ik,jp_tem) ! bottom before T and S347 zts (ji,jj,jp_sal) = ts b(ji,jj,ik,jp_sal)353 zts (ji,jj,jp_tem) = ts(ji,jj,ik,jp_tem,ktlev1) ! bottom before T and S 354 zts (ji,jj,jp_sal) = ts(ji,jj,ik,jp_sal,ktlev1) 348 355 ! 349 zdep(ji,jj) = gdept _n(ji,jj,ik) ! bottom T-level reference depth350 zub (ji,jj) = u n(ji,jj,mbku(ji,jj)) ! bottom velocity351 zvb (ji,jj) = v n(ji,jj,mbkv(ji,jj))356 zdep(ji,jj) = gdept(ji,jj,ik,kt2lev) ! bottom T-level reference depth 357 zub (ji,jj) = uu(ji,jj,mbku(ji,jj),ktlev2) ! bottom velocity 358 zvb (ji,jj) = vv(ji,jj,mbkv(ji,jj),ktlev2) 352 359 END DO 353 360 END DO -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/tradmp.F90
r10425 r10806 72 72 73 73 74 SUBROUTINE tra_dmp( kt )74 SUBROUTINE tra_dmp( kt, ktlev, kt2lev, pts_rhs ) 75 75 !!---------------------------------------------------------------------- 76 76 !! *** ROUTINE tra_dmp *** … … 90 90 !! ** Action : - tsa: tracer trends updated with the damping trend 91 91 !!---------------------------------------------------------------------- 92 INTEGER, INTENT(in) :: kt ! ocean time-step index 92 INTEGER, INTENT(in) :: kt ! ocean time-step index 93 INTEGER, INTENT(in) :: ktlev ! time level index for 3-time-level source terms 94 INTEGER, INTENT(in) :: kt2lev ! time level index for 2-time-level source terms 95 REAL(wp), INTENT( inout), DIMENSION(jpi,jpj,jpk,jpts) :: pts_rhs ! temperature and salinity trends 93 96 ! 94 97 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 101 104 IF( l_trdtra ) THEN !* Save ta and sa trends 102 105 ALLOCATE( ztrdts(jpi,jpj,jpk,jpts) ) 103 ztrdts(:,:,:,:) = tsa(:,:,:,:)106 ztrdts(:,:,:,:) = pts_rhs(:,:,:,:) 104 107 ENDIF 105 108 ! !== input T-S data at kt ==! … … 113 116 DO jj = 2, jpjm1 114 117 DO ji = fs_2, fs_jpim1 ! vector opt. 115 tsa(ji,jj,jk,jn) = tsa(ji,jj,jk,jn) + resto(ji,jj,jk) * ( zts_dta(ji,jj,jk,jn) - tsb(ji,jj,jk,jn) )118 pts_rhs(ji,jj,jk,jn) = pts_rhs(ji,jj,jk,jn) + resto(ji,jj,jk) * ( zts_dta(ji,jj,jk,jn) - ts(ji,jj,jk,jn,ktlev) ) 116 119 END DO 117 120 END DO … … 124 127 DO ji = fs_2, fs_jpim1 ! vector opt. 125 128 IF( avt(ji,jj,jk) <= avt_c ) THEN 126 tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) &127 & + resto(ji,jj,jk) * ( zts_dta(ji,jj,jk,jp_tem) - ts b(ji,jj,jk,jp_tem) )128 tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) &129 & + resto(ji,jj,jk) * ( zts_dta(ji,jj,jk,jp_sal) - ts b(ji,jj,jk,jp_sal) )129 pts_rhs(ji,jj,jk,jp_tem) = pts_rhs(ji,jj,jk,jp_tem) & 130 & + resto(ji,jj,jk) * ( zts_dta(ji,jj,jk,jp_tem) - ts(ji,jj,jk,jp_tem,ktlev) ) 131 pts_rhs(ji,jj,jk,jp_sal) = pts_rhs(ji,jj,jk,jp_sal) & 132 & + resto(ji,jj,jk) * ( zts_dta(ji,jj,jk,jp_sal) - ts(ji,jj,jk,jp_sal,ktlev) ) 130 133 ENDIF 131 134 END DO … … 137 140 DO jj = 2, jpjm1 138 141 DO ji = fs_2, fs_jpim1 ! vector opt. 139 IF( gdept _n(ji,jj,jk) >= hmlp (ji,jj) ) THEN140 tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) &141 & + resto(ji,jj,jk) * ( zts_dta(ji,jj,jk,jp_tem) - ts b(ji,jj,jk,jp_tem) )142 tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) &143 & + resto(ji,jj,jk) * ( zts_dta(ji,jj,jk,jp_sal) - ts b(ji,jj,jk,jp_sal) )142 IF( gdept(ji,jj,jk,kt2lev) >= hmlp (ji,jj) ) THEN 143 pts_rhs(ji,jj,jk,jp_tem) = pts_rhs(ji,jj,jk,jp_tem) & 144 & + resto(ji,jj,jk) * ( zts_dta(ji,jj,jk,jp_tem) - ts(ji,jj,jk,jp_tem,ktlev) ) 145 pts_rhs(ji,jj,jk,jp_sal) = pts_rhs(ji,jj,jk,jp_sal) & 146 & + resto(ji,jj,jk) * ( zts_dta(ji,jj,jk,jp_sal) - ts(ji,jj,jk,jp_sal,ktlev) ) 144 147 ENDIF 145 148 END DO … … 150 153 ! 151 154 IF( l_trdtra ) THEN ! trend diagnostic 152 ztrdts(:,:,:,:) = tsa(:,:,:,:) - ztrdts(:,:,:,:)155 ztrdts(:,:,:,:) = pts_rhs(:,:,:,:) - ztrdts(:,:,:,:) 153 156 CALL trd_tra( kt, 'TRA', jp_tem, jptra_dmp, ztrdts(:,:,:,jp_tem) ) 154 157 CALL trd_tra( kt, 'TRA', jp_sal, jptra_dmp, ztrdts(:,:,:,jp_sal) ) -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traldf.F90
r10068 r10806 47 47 CONTAINS 48 48 49 SUBROUTINE tra_ldf( kt )49 SUBROUTINE tra_ldf( kt, ktlev1, ktlev2, kt2lev, pts_rhs ) 50 50 !!---------------------------------------------------------------------- 51 51 !! *** ROUTINE tra_ldf *** … … 53 53 !! ** Purpose : compute the lateral ocean tracer physics. 54 54 !!---------------------------------------------------------------------- 55 INTEGER, INTENT( in ) :: kt ! ocean time-step index 55 INTEGER, INTENT( in ) :: kt ! ocean time-step index 56 INTEGER, INTENT( in ) :: ktlev1, ktlev2 ! time level indices for 3-time-level source terms 57 INTEGER, INTENT( in ) :: kt2lev ! time level index for 2-time-level source terms 58 REAL(wp), INTENT( inout), DIMENSION(jpi,jpj,jpk,jpts) :: pts_rhs ! temperature and salinity trends 56 59 !! 57 60 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdt, ztrds … … 62 65 IF( l_trdtra ) THEN !* Save ta and sa trends 63 66 ALLOCATE( ztrdt(jpi,jpj,jpk) , ztrds(jpi,jpj,jpk) ) 64 ztrdt(:,:,:) = tsa(:,:,:,jp_tem)65 ztrds(:,:,:) = tsa(:,:,:,jp_sal)67 ztrdt(:,:,:) = pts_rhs(:,:,:,jp_tem) 68 ztrds(:,:,:) = pts_rhs(:,:,:,jp_sal) 66 69 ENDIF 67 70 ! 68 71 SELECT CASE ( nldf_tra ) !* compute lateral mixing trend and add it to the general trend 69 72 CASE ( np_lap ) ! laplacian: iso-level operator 70 CALL tra_ldf_lap ( kt, nit000, 'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, tsb, tsa, jpts, 1 )73 CALL tra_ldf_lap ( kt, nit000, ktlev2, kt2lev, 'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, ts(:,:,:,:,ktlev1), pts_rhs, jpts, 1 ) 71 74 CASE ( np_lap_i ) ! laplacian: standard iso-neutral operator (Madec) 72 CALL tra_ldf_iso ( kt, nit000, 'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, tsb, tsb, tsa, jpts, 1 )75 CALL tra_ldf_iso ( kt, nit000, ktlev2, kt2lev, 'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, ts(:,:,:,:,ktlev1), ts(:,:,:,:,ktlev1), pts_rhs, jpts, 1 ) 73 76 CASE ( np_lap_it ) ! laplacian: triad iso-neutral operator (griffies) 74 CALL tra_ldf_triad( kt, nit000, 'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, tsb, tsb, tsa, jpts, 1 )77 CALL tra_ldf_triad( kt, nit000, ktlev2, kt2lev, 'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, ts(:,:,:,:,ktlev1), ts(:,:,:,:,ktlev1), pts_rhs, jpts, 1 ) 75 78 CASE ( np_blp , np_blp_i , np_blp_it ) ! bilaplacian: iso-level & iso-neutral operators 76 CALL tra_ldf_blp ( kt, nit000, 'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, tsb , tsa, jpts, nldf_tra )79 CALL tra_ldf_blp ( kt, nit000, ktlev2, kt2lev, 'TRA', ahtu, ahtv, gtsu, gtsv, gtui, gtvi, ts(:,:,:,:,ktlev1) , pts_rhs, jpts, nldf_tra ) 77 80 END SELECT 78 81 ! 79 82 IF( l_trdtra ) THEN !* save the horizontal diffusive trends for further diagnostics 80 ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)81 ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:)83 ztrdt(:,:,:) = pts_rhs(:,:,:,jp_tem) - ztrdt(:,:,:) 84 ztrds(:,:,:) = pts_rhs(:,:,:,jp_sal) - ztrds(:,:,:) 82 85 CALL trd_tra( kt, 'TRA', jp_tem, jptra_ldf, ztrdt ) 83 86 CALL trd_tra( kt, 'TRA', jp_sal, jptra_ldf, ztrds ) -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traldf_iso.F90
r10068 r10806 48 48 CONTAINS 49 49 50 SUBROUTINE tra_ldf_iso( kt, kit000, cdtype, pahu, pahv, pgu , pgv , &50 SUBROUTINE tra_ldf_iso( kt, kit000, ktlev, kt2lev, cdtype, pahu, pahv, pgu , pgv , & 51 51 & pgui, pgvi, & 52 & pt b , ptbb, pta, kjpt, kpass )52 & pt , pt_lev0, pt_rhs , kjpt, kpass ) 53 53 !!---------------------------------------------------------------------- 54 54 !! *** ROUTINE tra_ldf_iso *** … … 87 87 !! difft = 1/(e1e2t*e3t) dk[ zftw ] 88 88 !! Add this trend to the general trend (ta,sa): 89 !! pt a = pta+ difft90 !! 91 !! ** Action : Update pt aarrays with the before rotated diffusion89 !! pt_rhs = pt_rhs + difft 90 !! 91 !! ** Action : Update pt_rhs arrays with the before rotated diffusion 92 92 !!---------------------------------------------------------------------- 93 93 INTEGER , INTENT(in ) :: kt ! ocean time-step index 94 94 INTEGER , INTENT(in ) :: kit000 ! first time step index 95 INTEGER , INTENT(in ) :: ktlev ! time level index for e3t 96 INTEGER , INTENT(in ) :: kt2lev ! time level index for 2-time-level thicknesses 95 97 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 96 98 INTEGER , INTENT(in ) :: kjpt ! number of tracers … … 99 101 REAL(wp), DIMENSION(jpi,jpj ,kjpt), INTENT(in ) :: pgu, pgv ! tracer gradient at pstep levels 100 102 REAL(wp), DIMENSION(jpi,jpj, kjpt), INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 101 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt b! tracer (kpass=1) or laplacian of tracer (kpass=2)102 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt bb! tracer (only used in kpass=2)103 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt a! tracer trend103 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt ! tracer (kpass=1) or laplacian of tracer (kpass=2) 104 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt_lev0 ! tracer (only used in kpass=2) 105 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 104 106 ! 105 107 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 182 184 DO ji = 1, fs_jpim1 183 185 akz(ji,jj,jk) = 16._wp * ah_wslp2(ji,jj,jk) & 184 & * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( e3w _n(ji,jj,jk) * e3w_n(ji,jj,jk) ) )186 & * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( e3w(ji,jj,jk,kt2lev) * e3w(ji,jj,jk,kt2lev) ) ) 185 187 END DO 186 188 END DO … … 190 192 DO jj = 1, jpjm1 191 193 DO ji = 1, fs_jpim1 192 ze3w_2 = e3w _n(ji,jj,jk) * e3w_n(ji,jj,jk)194 ze3w_2 = e3w(ji,jj,jk,kt2lev) * e3w(ji,jj,jk,kt2lev) 193 195 zcoef0 = z2dt * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ze3w_2 ) 194 196 akz(ji,jj,jk) = MAX( zcoef0 - 0.5_wp , 0._wp ) * ze3w_2 * z1_2dt … … 219 221 DO jj = 1, jpjm1 220 222 DO ji = 1, fs_jpim1 ! vector opt. 221 zdit(ji,jj,jk) = ( pt b(ji+1,jj ,jk,jn) - ptb(ji,jj,jk,jn) ) * umask(ji,jj,jk)222 zdjt(ji,jj,jk) = ( pt b(ji ,jj+1,jk,jn) - ptb(ji,jj,jk,jn) ) * vmask(ji,jj,jk)223 zdit(ji,jj,jk) = ( pt(ji+1,jj ,jk,jn) - pt(ji,jj,jk,jn) ) * umask(ji,jj,jk) 224 zdjt(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 223 225 END DO 224 226 END DO … … 248 250 ! 249 251 ! !== Vertical tracer gradient 250 zdk1t(:,:) = ( pt b(:,:,jk,jn) - ptb(:,:,jk+1,jn) ) * wmask(:,:,jk+1) ! level jk+1252 zdk1t(:,:) = ( pt(:,:,jk,jn) - pt(:,:,jk+1,jn) ) * wmask(:,:,jk+1) ! level jk+1 251 253 ! 252 254 IF( jk == 1 ) THEN ; zdkt(:,:) = zdk1t(:,:) ! surface: zdkt(jk=1)=zdkt(jk=2) 253 ELSE ; zdkt(:,:) = ( pt b(:,:,jk-1,jn) - ptb(:,:,jk,jn) ) * wmask(:,:,jk)255 ELSE ; zdkt(:,:) = ( pt(:,:,jk-1,jn) - pt(:,:,jk,jn) ) * wmask(:,:,jk) 254 256 ENDIF 255 257 DO jj = 1 , jpjm1 !== Horizontal fluxes 256 258 DO ji = 1, fs_jpim1 ! vector opt. 257 zabe1 = pahu(ji,jj,jk) * e2_e1u(ji,jj) * e3u _n(ji,jj,jk)258 zabe2 = pahv(ji,jj,jk) * e1_e2v(ji,jj) * e3v _n(ji,jj,jk)259 zabe1 = pahu(ji,jj,jk) * e2_e1u(ji,jj) * e3u(ji,jj,jk,ktlev) 260 zabe2 = pahv(ji,jj,jk) * e1_e2v(ji,jj) * e3v(ji,jj,jk,ktlev) 259 261 ! 260 262 zmsku = 1. / MAX( wmask(ji+1,jj,jk ) + wmask(ji,jj,jk+1) & … … 276 278 END DO 277 279 ! 278 DO jj = 2 , jpjm1 !== horizontal divergence and add to pt a280 DO jj = 2 , jpjm1 !== horizontal divergence and add to pt_rhs 279 281 DO ji = fs_2, fs_jpim1 ! vector opt. 280 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + zsign * ( zftu(ji,jj,jk) - zftu(ji-1,jj,jk) &282 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * ( zftu(ji,jj,jk) - zftu(ji-1,jj,jk) & 281 283 & + zftv(ji,jj,jk) - zftv(ji,jj-1,jk) ) & 282 & * r1_e1e2t(ji,jj) / e3t _n(ji,jj,jk)284 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev) 283 285 END DO 284 286 END DO … … 325 327 DO jj = 1, jpjm1 326 328 DO ji = fs_2, fs_jpim1 327 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w _n(ji,jj,jk) * wmask(ji,jj,jk) &329 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w(ji,jj,jk,kt2lev) * wmask(ji,jj,jk) & 328 330 & * ( ah_wslp2(ji,jj,jk) - akz(ji,jj,jk) ) & 329 & * ( pt b(ji,jj,jk-1,jn) - ptb(ji,jj,jk,jn) )331 & * ( pt(ji,jj,jk-1,jn) - pt(ji,jj,jk,jn) ) 330 332 END DO 331 333 END DO … … 340 342 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) & 341 343 & + ah_wslp2(ji,jj,jk) * e1e2t(ji,jj) & 342 & * ( pt b(ji,jj,jk-1,jn) - ptb(ji,jj,jk,jn) ) / e3w_n(ji,jj,jk) * wmask(ji,jj,jk)344 & * ( pt(ji,jj,jk-1,jn) - pt(ji,jj,jk,jn) ) / e3w(ji,jj,jk,kt2lev) * wmask(ji,jj,jk) 343 345 END DO 344 346 END DO 345 347 END DO 346 CASE( 2 ) ! 2nd pass : eddy flux = ah_wslp2 and akz applied on pt b and ptbbgradients, resp.348 CASE( 2 ) ! 2nd pass : eddy flux = ah_wslp2 and akz applied on pt and pt_lev0 gradients, resp. 347 349 DO jk = 2, jpkm1 348 350 DO jj = 1, jpjm1 349 351 DO ji = fs_2, fs_jpim1 350 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w _n(ji,jj,jk) * wmask(ji,jj,jk) &351 & * ( ah_wslp2(ji,jj,jk) * ( pt b (ji,jj,jk-1,jn) - ptb(ji,jj,jk,jn) ) &352 & + akz (ji,jj,jk) * ( pt bb(ji,jj,jk-1,jn) - ptbb(ji,jj,jk,jn) ) )352 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) + e1e2t(ji,jj) / e3w(ji,jj,jk,kt2lev) * wmask(ji,jj,jk) & 353 & * ( ah_wslp2(ji,jj,jk) * ( pt (ji,jj,jk-1,jn) - pt (ji,jj,jk,jn) ) & 354 & + akz (ji,jj,jk) * ( pt_lev0(ji,jj,jk-1,jn) - pt_lev0(ji,jj,jk,jn) ) ) 353 355 END DO 354 356 END DO … … 357 359 ENDIF 358 360 ! 359 DO jk = 1, jpkm1 !== Divergence of vertical fluxes added to pt a==!361 DO jk = 1, jpkm1 !== Divergence of vertical fluxes added to pt_rhs ==! 360 362 DO jj = 2, jpjm1 361 363 DO ji = fs_2, fs_jpim1 ! vector opt. 362 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + zsign * ( ztfw (ji,jj,jk) - ztfw(ji,jj,jk+1) ) &363 & * r1_e1e2t(ji,jj) / e3t _n(ji,jj,jk)364 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * ( ztfw (ji,jj,jk) - ztfw(ji,jj,jk+1) ) & 365 & * r1_e1e2t(ji,jj) / e3t(ji,jj,jk,ktlev) 364 366 END DO 365 367 END DO -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traldf_lap_blp.F90
r10425 r10806 45 45 CONTAINS 46 46 47 SUBROUTINE tra_ldf_lap( kt, kit000, cdtype, pahu, pahv, pgu , pgv , &47 SUBROUTINE tra_ldf_lap( kt, kit000, ktlev, kt2lev, cdtype, pahu, pahv, pgu , pgv , & 48 48 & pgui, pgvi, & 49 & pt b , pta, kjpt, kpass )49 & pt , pt_rhs , kjpt, kpass ) 50 50 !!---------------------------------------------------------------------- 51 51 !! *** ROUTINE tra_ldf_lap *** … … 59 59 !! difft = 1/(e1e2t*e3t) { di-1[ pahu e2u*e3u/e1u di(tb) ] 60 60 !! + dj-1[ pahv e1v*e3v/e2v dj(tb) ] } 61 !! Add this trend to the general tracer trend pt a:62 !! pt a = pta+ difft63 !! 64 !! ** Action : - Update pt aarrays with the before iso-level61 !! Add this trend to the general tracer trend pt_rhs : 62 !! pt_rhs = pt_rhs + difft 63 !! 64 !! ** Action : - Update pt_rhs arrays with the before iso-level 65 65 !! harmonic mixing trend. 66 66 !!---------------------------------------------------------------------- 67 67 INTEGER , INTENT(in ) :: kt ! ocean time-step index 68 68 INTEGER , INTENT(in ) :: kit000 ! first time step index 69 INTEGER , INTENT(in ) :: ktlev ! time level index for e3t 70 INTEGER , INTENT(in ) :: kt2lev ! time level index for 2-time-level thicknesses 69 71 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 70 72 INTEGER , INTENT(in ) :: kjpt ! number of tracers … … 73 75 REAL(wp), DIMENSION(jpi,jpj ,kjpt), INTENT(in ) :: pgu, pgv ! tracer gradient at pstep levels 74 76 REAL(wp), DIMENSION(jpi,jpj, kjpt), INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 75 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt b! before and now tracer fields76 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt a! tracer trend77 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt ! before and now tracer fields 78 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 77 79 ! 78 80 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 100 102 DO jj = 1, jpjm1 101 103 DO ji = 1, fs_jpim1 ! vector opt. 102 zaheeu(ji,jj,jk) = zsign * pahu(ji,jj,jk) * e2_e1u(ji,jj) * e3u _n(ji,jj,jk) !!gm * umask(ji,jj,jk) pah masked!103 zaheev(ji,jj,jk) = zsign * pahv(ji,jj,jk) * e1_e2v(ji,jj) * e3v _n(ji,jj,jk) !!gm * vmask(ji,jj,jk)104 zaheeu(ji,jj,jk) = zsign * pahu(ji,jj,jk) * e2_e1u(ji,jj) * e3u(ji,jj,jk,ktlev) !!gm * umask(ji,jj,jk) pah masked! 105 zaheev(ji,jj,jk) = zsign * pahv(ji,jj,jk) * e1_e2v(ji,jj) * e3v(ji,jj,jk,ktlev) !!gm * vmask(ji,jj,jk) 104 106 END DO 105 107 END DO … … 113 115 DO jj = 1, jpjm1 114 116 DO ji = 1, fs_jpim1 115 ztu(ji,jj,jk) = zaheeu(ji,jj,jk) * ( pt b(ji+1,jj ,jk,jn) - ptb(ji,jj,jk,jn) )116 ztv(ji,jj,jk) = zaheev(ji,jj,jk) * ( pt b(ji ,jj+1,jk,jn) - ptb(ji,jj,jk,jn) )117 ztu(ji,jj,jk) = zaheeu(ji,jj,jk) * ( pt(ji+1,jj ,jk,jn) - pt(ji,jj,jk,jn) ) 118 ztv(ji,jj,jk) = zaheev(ji,jj,jk) * ( pt(ji ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) 117 119 END DO 118 120 END DO … … 138 140 DO jj = 2, jpjm1 139 141 DO ji = fs_2, fs_jpim1 140 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) &142 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + ( ztu(ji,jj,jk) - ztu(ji-1,jj,jk) & 141 143 & + ztv(ji,jj,jk) - ztv(ji,jj-1,jk) ) & 142 & / ( e1e2t(ji,jj) * e3t _n(ji,jj,jk) )144 & / ( e1e2t(ji,jj) * e3t(ji,jj,jk,ktlev) ) 143 145 END DO 144 146 END DO … … 159 161 160 162 161 SUBROUTINE tra_ldf_blp( kt, kit000, cdtype, pahu, pahv, pgu , pgv , &163 SUBROUTINE tra_ldf_blp( kt, kit000, ktlev, kt2lev, cdtype, pahu, pahv, pgu , pgv , & 162 164 & pgui, pgvi, & 163 & pt b , pta, kjpt, kldf )165 & pt , pt_rhs , kjpt, kldf ) 164 166 !!---------------------------------------------------------------------- 165 167 !! *** ROUTINE tra_ldf_blp *** … … 172 174 !! It is computed by two successive calls to laplacian routine 173 175 !! 174 !! ** Action : pt aupdated with the before rotated bilaplacian diffusion176 !! ** Action : pt_rhs updated with the before rotated bilaplacian diffusion 175 177 !!---------------------------------------------------------------------- 176 178 INTEGER , INTENT(in ) :: kt ! ocean time-step index 177 179 INTEGER , INTENT(in ) :: kit000 ! first time step index 180 INTEGER , INTENT(in ) :: ktlev ! time level index for e3t 181 INTEGER , INTENT(in ) :: kt2lev ! time level index for 2-time-level thicknesses 178 182 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 179 183 INTEGER , INTENT(in ) :: kjpt ! number of tracers … … 182 186 REAL(wp), DIMENSION(jpi,jpj ,kjpt), INTENT(in ) :: pgu, pgv ! tracer gradient at pstep levels 183 187 REAL(wp), DIMENSION(jpi,jpj, kjpt), INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 184 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt b! before and now tracer fields185 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt a! tracer trend188 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt ! before and now tracer fields 189 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 186 190 ! 187 191 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 203 207 zlap(:,:,:,:) = 0._wp 204 208 ! 205 SELECT CASE ( kldf ) !== 1st laplacian applied to pt b(output in zlap) ==!209 SELECT CASE ( kldf ) !== 1st laplacian applied to pt (output in zlap) ==! 206 210 ! 207 211 CASE ( np_blp ) ! iso-level bilaplacian 208 CALL tra_ldf_lap ( kt, kit000, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, ptb, zlap, kjpt, 1 )212 CALL tra_ldf_lap ( kt, kit000, ktlev, kt2lev, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, pt, zlap, kjpt, 1 ) 209 213 CASE ( np_blp_i ) ! rotated bilaplacian : standard operator (Madec) 210 CALL tra_ldf_iso ( kt, kit000, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, ptb, ptb, zlap, kjpt, 1 )214 CALL tra_ldf_iso ( kt, kit000, ktlev, kt2lev, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, pt, pt, zlap, kjpt, 1 ) 211 215 CASE ( np_blp_it ) ! rotated bilaplacian : triad operator (griffies) 212 CALL tra_ldf_triad( kt, kit000, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, ptb, ptb, zlap, kjpt, 1 )216 CALL tra_ldf_triad( kt, kit000, ktlev, kt2lev, cdtype, pahu, pahv, pgu, pgv, pgui, pgvi, pt, pt, zlap, kjpt, 1 ) 213 217 END SELECT 214 218 ! … … 219 223 ENDIF 220 224 ! 221 SELECT CASE ( kldf ) !== 2nd laplacian applied to zlap (output in pt a) ==!225 SELECT CASE ( kldf ) !== 2nd laplacian applied to zlap (output in pt_rhs) ==! 222 226 ! 223 227 CASE ( np_blp ) ! iso-level bilaplacian 224 CALL tra_ldf_lap ( kt, kit000, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, pta, kjpt, 2 )228 CALL tra_ldf_lap ( kt, kit000, ktlev, kt2lev, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, pt_rhs, kjpt, 2 ) 225 229 CASE ( np_blp_i ) ! rotated bilaplacian : standard operator (Madec) 226 CALL tra_ldf_iso ( kt, kit000, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, ptb, pta, kjpt, 2 )230 CALL tra_ldf_iso ( kt, kit000, ktlev, kt2lev, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, pt, pt_rhs, kjpt, 2 ) 227 231 CASE ( np_blp_it ) ! rotated bilaplacian : triad operator (griffies) 228 CALL tra_ldf_triad( kt, kit000, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, ptb, pta, kjpt, 2 )232 CALL tra_ldf_triad( kt, kit000, ktlev, kt2lev, cdtype, pahu, pahv, zglu, zglv, zgui, zgvi, zlap, pt, pt_rhs, kjpt, 2 ) 229 233 END SELECT 230 234 ! -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traldf_triad.F90
r10425 r10806 48 48 CONTAINS 49 49 50 SUBROUTINE tra_ldf_triad( kt, kit000, cdtype, pahu, pahv, pgu , pgv , &50 SUBROUTINE tra_ldf_triad( kt, kit000, ktlev, kt2lev, cdtype, pahu, pahv, pgu , pgv , & 51 51 & pgui, pgvi, & 52 & pt b , ptbb, pta, kjpt, kpass )52 & pt , pt_lev0, pt_rhs , kjpt, kpass ) 53 53 !!---------------------------------------------------------------------- 54 54 !! *** ROUTINE tra_ldf_triad *** … … 66 66 !! see documentation for the desciption 67 67 !! 68 !! ** Action : pt aupdated with the before rotated diffusion68 !! ** Action : pt_rhs updated with the before rotated diffusion 69 69 !! ah_wslp2 .... 70 70 !! akz stabilizing vertical diffusivity coefficient (used in trazdf_imp) … … 72 72 INTEGER , INTENT(in ) :: kt ! ocean time-step index 73 73 INTEGER , INTENT(in ) :: kit000 ! first time step index 74 INTEGER , INTENT(in ) :: ktlev ! time level index for e3t 75 INTEGER , INTENT(in ) :: kt2lev ! time level index for 2-time-level thicknesses 74 76 CHARACTER(len=3) , INTENT(in ) :: cdtype ! =TRA or TRC (tracer indicator) 75 77 INTEGER , INTENT(in ) :: kjpt ! number of tracers … … 78 80 REAL(wp), DIMENSION(jpi,jpj ,kjpt), INTENT(in ) :: pgu , pgv ! tracer gradient at pstep levels 79 81 REAL(wp), DIMENSION(jpi,jpj, kjpt), INTENT(in ) :: pgui, pgvi ! tracer gradient at top levels 80 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt b! tracer (kpass=1) or laplacian of tracer (kpass=2)81 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt bb! tracer (only used in kpass=2)82 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt a! tracer trend82 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt ! tracer (kpass=1) or laplacian of tracer (kpass=2) 83 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(in ) :: pt_lev0 ! tracer (only used in kpass=2) 84 REAL(wp), DIMENSION(jpi,jpj,jpk,kjpt), INTENT(inout) :: pt_rhs ! tracer trend 83 85 ! 84 86 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 142 144 DO jj = 1, jpjm1 143 145 DO ji = 1, fs_jpim1 144 ze3wr = 1._wp / e3w _n(ji+ip,jj,jk+kp)145 zbu = e1e2u(ji,jj) * e3u _n(ji,jj,jk)146 ze3wr = 1._wp / e3w(ji+ip,jj,jk+kp,kt2lev) 147 zbu = e1e2u(ji,jj) * e3u(ji,jj,jk,ktlev) 146 148 zah = 0.25_wp * pahu(ji,jj,jk) 147 149 zslope_skew = triadi_g(ji+ip,jj,jk,1-ip,kp) 148 150 ! Subtract s-coordinate slope at t-points to give slope rel to s-surfaces (do this by *adding* gradient of depth) 149 zslope2 = zslope_skew + ( gdept _n(ji+1,jj,jk) - gdept_n(ji,jj,jk) ) * r1_e1u(ji,jj) * umask(ji,jj,jk+kp)151 zslope2 = zslope_skew + ( gdept(ji+1,jj,jk,kt2lev) - gdept(ji,jj,jk,kt2lev) ) * r1_e1u(ji,jj) * umask(ji,jj,jk+kp) 150 152 zslope2 = zslope2 *zslope2 151 153 ah_wslp2(ji+ip,jj,jk+kp) = ah_wslp2(ji+ip,jj,jk+kp) + zah * zbu * ze3wr * r1_e1e2t(ji+ip,jj) * zslope2 … … 166 168 DO jj = 1, jpjm1 167 169 DO ji = 1, fs_jpim1 168 ze3wr = 1.0_wp / e3w _n(ji,jj+jp,jk+kp)169 zbv = e1e2v(ji,jj) * e3v _n(ji,jj,jk)170 ze3wr = 1.0_wp / e3w(ji,jj+jp,jk+kp,kt2lev) 171 zbv = e1e2v(ji,jj) * e3v(ji,jj,jk,ktlev) 170 172 zah = 0.25_wp * pahv(ji,jj,jk) 171 173 zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 172 174 ! Subtract s-coordinate slope at t-points to give slope rel to s surfaces 173 175 ! (do this by *adding* gradient of depth) 174 zslope2 = zslope_skew + ( gdept _n(ji,jj+1,jk) - gdept_n(ji,jj,jk) ) * r1_e2v(ji,jj) * vmask(ji,jj,jk+kp)176 zslope2 = zslope_skew + ( gdept(ji,jj+1,jk,kt2lev) - gdept(ji,jj,jk,kt2lev) ) * r1_e2v(ji,jj) * vmask(ji,jj,jk+kp) 175 177 zslope2 = zslope2 * zslope2 176 178 ah_wslp2(ji,jj+jp,jk+kp) = ah_wslp2(ji,jj+jp,jk+kp) + zah * zbv * ze3wr * r1_e1e2t(ji,jj+jp) * zslope2 … … 193 195 DO ji = 1, fs_jpim1 194 196 akz(ji,jj,jk) = 16._wp * ah_wslp2(ji,jj,jk) & 195 & * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( e3w _n(ji,jj,jk) * e3w_n(ji,jj,jk) ) )197 & * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ( e3w(ji,jj,jk,kt2lev) * e3w(ji,jj,jk,kt2lev) ) ) 196 198 END DO 197 199 END DO … … 201 203 DO jj = 1, jpjm1 202 204 DO ji = 1, fs_jpim1 203 ze3w_2 = e3w _n(ji,jj,jk) * e3w_n(ji,jj,jk)205 ze3w_2 = e3w(ji,jj,jk,kt2lev) * e3w(ji,jj,jk,kt2lev) 204 206 zcoef0 = z2dt * ( akz(ji,jj,jk) + ah_wslp2(ji,jj,jk) / ze3w_2 ) 205 207 akz(ji,jj,jk) = MAX( zcoef0 - 0.5_wp , 0._wp ) * ze3w_2 * z1_2dt … … 229 231 DO jj = 1, jpjm1 230 232 DO ji = 1, fs_jpim1 ! vector opt. 231 zdit(ji,jj,jk) = ( pt b(ji+1,jj ,jk,jn) - ptb(ji,jj,jk,jn) ) * umask(ji,jj,jk)232 zdjt(ji,jj,jk) = ( pt b(ji ,jj+1,jk,jn) - ptb(ji,jj,jk,jn) ) * vmask(ji,jj,jk)233 zdit(ji,jj,jk) = ( pt(ji+1,jj ,jk,jn) - pt(ji,jj,jk,jn) ) * umask(ji,jj,jk) 234 zdjt(ji,jj,jk) = ( pt(ji ,jj+1,jk,jn) - pt(ji,jj,jk,jn) ) * vmask(ji,jj,jk) 233 235 END DO 234 236 END DO … … 257 259 DO jk = 1, jpkm1 258 260 ! !== Vertical tracer gradient at level jk and jk+1 259 zdkt3d(:,:,1) = ( pt b(:,:,jk,jn) - ptb(:,:,jk+1,jn) ) * tmask(:,:,jk+1)261 zdkt3d(:,:,1) = ( pt(:,:,jk,jn) - pt(:,:,jk+1,jn) ) * tmask(:,:,jk+1) 260 262 ! 261 263 ! ! surface boundary condition: zdkt3d(jk=0)=zdkt3d(jk=1) 262 264 IF( jk == 1 ) THEN ; zdkt3d(:,:,0) = zdkt3d(:,:,1) 263 ELSE ; zdkt3d(:,:,0) = ( pt b(:,:,jk-1,jn) - ptb(:,:,jk,jn) ) * tmask(:,:,jk)265 ELSE ; zdkt3d(:,:,0) = ( pt(:,:,jk-1,jn) - pt(:,:,jk,jn) ) * tmask(:,:,jk) 264 266 ENDIF 265 267 ! … … 273 275 ze1ur = r1_e1u(ji,jj) 274 276 zdxt = zdit(ji,jj,jk) * ze1ur 275 ze3wr = 1._wp / e3w _n(ji+ip,jj,jk+kp)277 ze3wr = 1._wp / e3w(ji+ip,jj,jk+kp,kt2lev) 276 278 zdzt = zdkt3d(ji+ip,jj,kp) * ze3wr 277 279 zslope_skew = triadi_g(ji+ip,jj,jk,1-ip,kp) 278 280 zslope_iso = triadi (ji+ip,jj,jk,1-ip,kp) 279 281 ! 280 zbu = 0.25_wp * e1e2u(ji,jj) * e3u _n(ji,jj,jk)282 zbu = 0.25_wp * e1e2u(ji,jj) * e3u(ji,jj,jk,ktlev) 281 283 ! ln_botmix_triad is .T. don't mask zah for bottom half cells !!gm ????? ahu is masked.... 282 284 zah = pahu(ji,jj,jk) … … 296 298 ze2vr = r1_e2v(ji,jj) 297 299 zdyt = zdjt(ji,jj,jk) * ze2vr 298 ze3wr = 1._wp / e3w _n(ji,jj+jp,jk+kp)300 ze3wr = 1._wp / e3w(ji,jj+jp,jk+kp,kt2lev) 299 301 zdzt = zdkt3d(ji,jj+jp,kp) * ze3wr 300 302 zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 301 303 zslope_iso = triadj(ji,jj+jp,jk,1-jp,kp) 302 zbv = 0.25_wp * e1e2v(ji,jj) * e3v _n(ji,jj,jk)304 zbv = 0.25_wp * e1e2v(ji,jj) * e3v(ji,jj,jk,ktlev) 303 305 ! ln_botmix_triad is .T. don't mask zah for bottom half cells !!gm ????? ahv is masked... 304 306 zah = pahv(ji,jj,jk) … … 320 322 ze1ur = r1_e1u(ji,jj) 321 323 zdxt = zdit(ji,jj,jk) * ze1ur 322 ze3wr = 1._wp / e3w _n(ji+ip,jj,jk+kp)324 ze3wr = 1._wp / e3w(ji+ip,jj,jk+kp,kt2lev) 323 325 zdzt = zdkt3d(ji+ip,jj,kp) * ze3wr 324 326 zslope_skew = triadi_g(ji+ip,jj,jk,1-ip,kp) 325 327 zslope_iso = triadi(ji+ip,jj,jk,1-ip,kp) 326 328 ! 327 zbu = 0.25_wp * e1e2u(ji,jj) * e3u _n(ji,jj,jk)329 zbu = 0.25_wp * e1e2u(ji,jj) * e3u(ji,jj,jk,ktlev) 328 330 ! ln_botmix_triad is .F. mask zah for bottom half cells 329 331 zah = pahu(ji,jj,jk) * umask(ji,jj,jk+kp) ! pahu(ji+ip,jj,jk) ===>> ???? … … 343 345 ze2vr = r1_e2v(ji,jj) 344 346 zdyt = zdjt(ji,jj,jk) * ze2vr 345 ze3wr = 1._wp / e3w _n(ji,jj+jp,jk+kp)347 ze3wr = 1._wp / e3w(ji,jj+jp,jk+kp,kt2lev) 346 348 zdzt = zdkt3d(ji,jj+jp,kp) * ze3wr 347 349 zslope_skew = triadj_g(ji,jj+jp,jk,1-jp,kp) 348 350 zslope_iso = triadj(ji,jj+jp,jk,1-jp,kp) 349 zbv = 0.25_wp * e1e2v(ji,jj) * e3v _n(ji,jj,jk)351 zbv = 0.25_wp * e1e2v(ji,jj) * e3v(ji,jj,jk,ktlev) 350 352 ! ln_botmix_triad is .F. mask zah for bottom half cells 351 353 zah = pahv(ji,jj,jk) * vmask(ji,jj,jk+kp) ! pahv(ji,jj+jp,jk) ???? … … 362 364 DO jj = 2 , jpjm1 363 365 DO ji = fs_2, fs_jpim1 ! vector opt. 364 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + zsign * ( zftu(ji-1,jj,jk) - zftu(ji,jj,jk) &366 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * ( zftu(ji-1,jj,jk) - zftu(ji,jj,jk) & 365 367 & + zftv(ji,jj-1,jk) - zftv(ji,jj,jk) ) & 366 & / ( e1e2t(ji,jj) * e3t _n(ji,jj,jk) )368 & / ( e1e2t(ji,jj) * e3t(ji,jj,jk,ktlev) ) 367 369 END DO 368 370 END DO … … 375 377 DO jj = 1, jpjm1 376 378 DO ji = fs_2, fs_jpim1 377 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) - e1e2t(ji,jj) / e3w _n(ji,jj,jk) * tmask(ji,jj,jk) &379 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) - e1e2t(ji,jj) / e3w(ji,jj,jk,kt2lev) * tmask(ji,jj,jk) & 378 380 & * ( ah_wslp2(ji,jj,jk) - akz(ji,jj,jk) ) & 379 & * ( pt b(ji,jj,jk-1,jn) - ptb(ji,jj,jk,jn) )381 & * ( pt(ji,jj,jk-1,jn) - pt(ji,jj,jk,jn) ) 380 382 END DO 381 383 END DO … … 387 389 DO jj = 1, jpjm1 388 390 DO ji = fs_2, fs_jpim1 389 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) - e1e2t(ji,jj) / e3w _n(ji,jj,jk) * tmask(ji,jj,jk) &390 & * ah_wslp2(ji,jj,jk) * ( pt b(ji,jj,jk-1,jn) - ptb(ji,jj,jk,jn) )391 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) - e1e2t(ji,jj) / e3w(ji,jj,jk,kt2lev) * tmask(ji,jj,jk) & 392 & * ah_wslp2(ji,jj,jk) * ( pt(ji,jj,jk-1,jn) - pt(ji,jj,jk,jn) ) 391 393 END DO 392 394 END DO 393 395 END DO 394 CASE( 2 ) ! 2nd pass : eddy flux = ah_wslp2 and akz applied on pt b and ptbbgradients, resp.396 CASE( 2 ) ! 2nd pass : eddy flux = ah_wslp2 and akz applied on pt and pt_lev0 gradients, resp. 395 397 DO jk = 2, jpkm1 396 398 DO jj = 1, jpjm1 397 399 DO ji = fs_2, fs_jpim1 398 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) - e1e2t(ji,jj) / e3w _n(ji,jj,jk) * tmask(ji,jj,jk) &399 & * ( ah_wslp2(ji,jj,jk) * ( pt b (ji,jj,jk-1,jn) - ptb(ji,jj,jk,jn) ) &400 & + akz (ji,jj,jk) * ( pt bb(ji,jj,jk-1,jn) - ptbb(ji,jj,jk,jn) ) )400 ztfw(ji,jj,jk) = ztfw(ji,jj,jk) - e1e2t(ji,jj) / e3w(ji,jj,jk,kt2lev) * tmask(ji,jj,jk) & 401 & * ( ah_wslp2(ji,jj,jk) * ( pt (ji,jj,jk-1,jn) - pt (ji,jj,jk,jn) ) & 402 & + akz (ji,jj,jk) * ( pt_lev0(ji,jj,jk-1,jn) - pt_lev0(ji,jj,jk,jn) ) ) 401 403 END DO 402 404 END DO … … 405 407 ENDIF 406 408 ! 407 DO jk = 1, jpkm1 !== Divergence of vertical fluxes added to pt a==!409 DO jk = 1, jpkm1 !== Divergence of vertical fluxes added to pt_rhs ==! 408 410 DO jj = 2, jpjm1 409 411 DO ji = fs_2, fs_jpim1 ! vector opt. 410 pt a(ji,jj,jk,jn) = pta(ji,jj,jk,jn) + zsign * ( ztfw(ji,jj,jk+1) - ztfw(ji,jj,jk) ) &411 & / ( e1e2t(ji,jj) * e3t _n(ji,jj,jk) )412 pt_rhs(ji,jj,jk,jn) = pt_rhs(ji,jj,jk,jn) + zsign * ( ztfw(ji,jj,jk+1) - ztfw(ji,jj,jk) ) & 413 & / ( e1e2t(ji,jj) * e3t(ji,jj,jk,ktlev) ) 412 414 END DO 413 415 END DO -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/traqsr.F90
r10425 r10806 75 75 CONTAINS 76 76 77 SUBROUTINE tra_qsr( kt )77 SUBROUTINE tra_qsr( kt, ktlev, kt2lev, pts_rhs ) 78 78 !!---------------------------------------------------------------------- 79 79 !! *** ROUTINE tra_qsr *** … … 102 102 !!---------------------------------------------------------------------- 103 103 INTEGER, INTENT(in) :: kt ! ocean time-step 104 INTEGER, INTENT(in) :: ktlev ! time level index for 3-time-level source terms 105 INTEGER, INTENT(in) :: kt2lev ! time level index for 2-time-level source terms 106 REAL(wp), INTENT( inout), DIMENSION(jpi,jpj,jpk,jpts) :: pts_rhs ! temperature and salinity trends 104 107 ! 105 108 INTEGER :: ji, jj, jk ! dummy loop indices … … 126 129 IF( l_trdtra ) THEN ! trends diagnostic: save the input temperature trend 127 130 ALLOCATE( ztrdt(jpi,jpj,jpk) ) 128 ztrdt(:,:,:) = tsa(:,:,:,jp_tem)131 ztrdt(:,:,:) = pts_rhs(:,:,:,jp_tem) 129 132 ENDIF 130 133 ! … … 172 175 zze = 568.2 * zCtot**(-0.746) 173 176 IF( zze > 102. ) zze = 200.0 * zCtot**(-0.293) 174 zpsi = gdepw _n(ji,jj,jk) / zze177 zpsi = gdepw(ji,jj,jk,kt2lev) / zze 175 178 ! 176 179 zlogc = LOG( zchl ) … … 218 221 DO jj = 2, jpjm1 219 222 DO ji = fs_2, fs_jpim1 220 zc0 = ze0(ji,jj,jk-1) * EXP( - e3t _n(ji,jj,jk-1) * xsi0r )221 zc1 = ze1(ji,jj,jk-1) * EXP( - e3t _n(ji,jj,jk-1) * zekb(ji,jj) )222 zc2 = ze2(ji,jj,jk-1) * EXP( - e3t _n(ji,jj,jk-1) * zekg(ji,jj) )223 zc3 = ze3(ji,jj,jk-1) * EXP( - e3t _n(ji,jj,jk-1) * zekr(ji,jj) )223 zc0 = ze0(ji,jj,jk-1) * EXP( - e3t(ji,jj,jk-1,ktlev) * xsi0r ) 224 zc1 = ze1(ji,jj,jk-1) * EXP( - e3t(ji,jj,jk-1,ktlev) * zekb(ji,jj) ) 225 zc2 = ze2(ji,jj,jk-1) * EXP( - e3t(ji,jj,jk-1,ktlev) * zekg(ji,jj) ) 226 zc3 = ze3(ji,jj,jk-1) * EXP( - e3t(ji,jj,jk-1,ktlev) * zekr(ji,jj) ) 224 227 ze0(ji,jj,jk) = zc0 225 228 ze1(ji,jj,jk) = zc1 … … 248 251 DO jj = 2, jpjm1 249 252 DO ji = fs_2, fs_jpim1 250 zc0 = zz0 * EXP( -gdepw _n(ji,jj,jk )*xsi0r ) + zz1 * EXP( -gdepw_n(ji,jj,jk)*xsi1r )251 zc1 = zz0 * EXP( -gdepw _n(ji,jj,jk+1)*xsi0r ) + zz1 * EXP( -gdepw_n(ji,jj,jk+1)*xsi1r )253 zc0 = zz0 * EXP( -gdepw(ji,jj,jk ,kt2lev)*xsi0r ) + zz1 * EXP( -gdepw(ji,jj,jk ,kt2lev)*xsi1r ) 254 zc1 = zz0 * EXP( -gdepw(ji,jj,jk+1,kt2lev)*xsi0r ) + zz1 * EXP( -gdepw(ji,jj,jk+1,kt2lev)*xsi1r ) 252 255 qsr_hc(ji,jj,jk) = qsr(ji,jj) * ( zc0 * wmask(ji,jj,jk) - zc1 * wmask(ji,jj,jk+1) ) 253 256 END DO … … 261 264 DO jj = 2, jpjm1 !-----------------------------! 262 265 DO ji = fs_2, fs_jpim1 ! vector opt. 263 tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) &264 & + z1_2 * ( qsr_hc_b(ji,jj,jk) + qsr_hc(ji,jj,jk) ) / e3t _n(ji,jj,jk)266 pts_rhs(ji,jj,jk,jp_tem) = pts_rhs(ji,jj,jk,jp_tem) & 267 & + z1_2 * ( qsr_hc_b(ji,jj,jk) + qsr_hc(ji,jj,jk) ) / e3t(ji,jj,jk,ktlev) 265 268 END DO 266 269 END DO … … 295 298 ! 296 299 IF( l_trdtra ) THEN ! qsr tracers trends saved for diagnostics 297 ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)300 ztrdt(:,:,:) = pts_rhs(:,:,:,jp_tem) - ztrdt(:,:,:) 298 301 CALL trd_tra( kt, 'TRA', jp_tem, jptra_qsr, ztrdt ) 299 302 DEALLOCATE( ztrdt ) -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/TRA/trasbc.F90
r10499 r10806 51 51 CONTAINS 52 52 53 SUBROUTINE tra_sbc ( kt )53 SUBROUTINE tra_sbc ( kt, ktlev, pts_rhs ) 54 54 !!---------------------------------------------------------------------- 55 55 !! *** ROUTINE tra_sbc *** … … 63 63 !! (2) Fwe , tracer carried with the water that is exchanged with air+ice. 64 64 !! The input forcing fields (emp, rnf, sfx, isf) contain Fext+Fwe, 65 !! they are simply added to the tracer trend ( tsa).65 !! they are simply added to the tracer trend (pts_rhs). 66 66 !! In linear free surface case (ln_linssh=T), the volume of the 67 67 !! ocean does not change with the water exchanges at the (air+ice)-sea … … 69 69 !! concentration/dilution effect associated with water exchanges. 70 70 !! 71 !! ** Action : - Update tsawith the surface boundary condition trend71 !! ** Action : - Update pts_rhs with the surface boundary condition trend 72 72 !! - send trends to trdtra module for further diagnostics(l_trdtra=T) 73 73 !!---------------------------------------------------------------------- 74 INTEGER, INTENT(in) :: kt ! ocean time-step index 74 INTEGER, INTENT(in) :: kt ! ocean time-step index 75 INTEGER, INTENT(in) :: ktlev ! time level index for source terms 76 REAL(wp), INTENT( inout), DIMENSION(jpi,jpj,jpk,jpts) :: pts_rhs ! temperature and salinity trends 75 77 ! 76 78 INTEGER :: ji, jj, jk, jn ! dummy loop indices … … 90 92 IF( l_trdtra ) THEN !* Save ta and sa trends 91 93 ALLOCATE( ztrdt(jpi,jpj,jpk) , ztrds(jpi,jpj,jpk) ) 92 ztrdt(:,:,:) = tsa(:,:,:,jp_tem)93 ztrds(:,:,:) = tsa(:,:,:,jp_sal)94 ztrdt(:,:,:) = pts_rhs(:,:,:,jp_tem) 95 ztrds(:,:,:) = pts_rhs(:,:,:,jp_sal) 94 96 ENDIF 95 97 ! … … 131 133 DO jj = 2, jpj !==>> add concentration/dilution effect due to constant volume cell 132 134 DO ji = fs_2, fs_jpim1 ! vector opt. 133 sbc_tsc(ji,jj,jp_tem) = sbc_tsc(ji,jj,jp_tem) + r1_rau0 * emp(ji,jj) * ts n(ji,jj,1,jp_tem)134 sbc_tsc(ji,jj,jp_sal) = sbc_tsc(ji,jj,jp_sal) + r1_rau0 * emp(ji,jj) * ts n(ji,jj,1,jp_sal)135 sbc_tsc(ji,jj,jp_tem) = sbc_tsc(ji,jj,jp_tem) + r1_rau0 * emp(ji,jj) * ts(ji,jj,1,jp_tem,ktlev) 136 sbc_tsc(ji,jj,jp_sal) = sbc_tsc(ji,jj,jp_sal) + r1_rau0 * emp(ji,jj) * ts(ji,jj,1,jp_sal,ktlev) 135 137 END DO 136 138 END DO !==>> output c./d. term 137 IF( iom_use('emp_x_sst') ) CALL iom_put( "emp_x_sst", emp (:,:) * ts n(:,:,1,jp_tem) )138 IF( iom_use('emp_x_sss') ) CALL iom_put( "emp_x_sss", emp (:,:) * ts n(:,:,1,jp_sal) )139 IF( iom_use('emp_x_sst') ) CALL iom_put( "emp_x_sst", emp (:,:) * ts(:,:,1,jp_tem,ktlev) ) 140 IF( iom_use('emp_x_sss') ) CALL iom_put( "emp_x_sss", emp (:,:) * ts(:,:,1,jp_sal,ktlev) ) 139 141 ENDIF 140 142 ! … … 142 144 DO jj = 2, jpj 143 145 DO ji = fs_2, fs_jpim1 ! vector opt. 144 tsa(ji,jj,1,jn) = tsa(ji,jj,1,jn) + zfact * ( sbc_tsc_b(ji,jj,jn) + sbc_tsc(ji,jj,jn) ) / e3t_n(ji,jj,1)146 pts_rhs(ji,jj,1,jn) = pts_rhs(ji,jj,1,jn) + zfact * ( sbc_tsc_b(ji,jj,jn) + sbc_tsc(ji,jj,jn) ) / e3t(ji,jj,1,ktlev) 145 147 END DO 146 148 END DO … … 173 175 DO jk = ikt, ikb - 1 174 176 ! compute trend 175 tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) &177 pts_rhs(ji,jj,jk,jp_tem) = pts_rhs(ji,jj,jk,jp_tem) & 176 178 & + zfact * ( risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem) ) & 177 179 & * r1_hisf_tbl(ji,jj) … … 180 182 ! level partially include in ice shelf boundary layer 181 183 ! compute trend 182 tsa(ji,jj,ikb,jp_tem) = tsa(ji,jj,ikb,jp_tem) &184 pts_rhs(ji,jj,ikb,jp_tem) = pts_rhs(ji,jj,ikb,jp_tem) & 183 185 & + zfact * ( risf_tsc_b(ji,jj,jp_tem) + risf_tsc(ji,jj,jp_tem) ) & 184 186 & * r1_hisf_tbl(ji,jj) * ralpha(ji,jj) … … 199 201 zdep = zfact / h_rnf(ji,jj) 200 202 DO jk = 1, nk_rnf(ji,jj) 201 tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) &203 pts_rhs(ji,jj,jk,jp_tem) = pts_rhs(ji,jj,jk,jp_tem) & 202 204 & + ( rnf_tsc_b(ji,jj,jp_tem) + rnf_tsc(ji,jj,jp_tem) ) * zdep 203 IF( ln_rnf_sal ) tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) &205 IF( ln_rnf_sal ) pts_rhs(ji,jj,jk,jp_sal) = pts_rhs(ji,jj,jk,jp_sal) & 204 206 & + ( rnf_tsc_b(ji,jj,jp_sal) + rnf_tsc(ji,jj,jp_sal) ) * zdep 205 207 END DO … … 209 211 ENDIF 210 212 211 IF( iom_use('rnf_x_sst') ) CALL iom_put( "rnf_x_sst", rnf*ts n(:,:,1,jp_tem) ) ! runoff term on sst212 IF( iom_use('rnf_x_sss') ) CALL iom_put( "rnf_x_sss", rnf*ts n(:,:,1,jp_sal) ) ! runoff term on sss213 IF( iom_use('rnf_x_sst') ) CALL iom_put( "rnf_x_sst", rnf*ts(:,:,1,jp_tem,ktlev) ) ! runoff term on sst 214 IF( iom_use('rnf_x_sss') ) CALL iom_put( "rnf_x_sss", rnf*ts(:,:,1,jp_sal,ktlev) ) ! runoff term on sss 213 215 214 216 #if defined key_asminc … … 223 225 DO jj = 2, jpj 224 226 DO ji = fs_2, fs_jpim1 225 ztim = ssh_iau(ji,jj) / e3t _n(ji,jj,1)226 tsa(ji,jj,1,jp_tem) = tsa(ji,jj,1,jp_tem) + tsn(ji,jj,1,jp_tem) * ztim227 tsa(ji,jj,1,jp_sal) = tsa(ji,jj,1,jp_sal) + tsn(ji,jj,1,jp_sal) * ztim227 ztim = ssh_iau(ji,jj) / e3t(ji,jj,1,ktlev) 228 pts_rhs(ji,jj,1,jp_tem) = pts_rhs(ji,jj,1,jp_tem) + ts(ji,jj,1,jp_tem,ktlev) * ztim 229 pts_rhs(ji,jj,1,jp_sal) = pts_rhs(ji,jj,1,jp_sal) + ts(ji,jj,1,jp_sal,ktlev) * ztim 228 230 END DO 229 231 END DO … … 232 234 DO ji = fs_2, fs_jpim1 233 235 ztim = ssh_iau(ji,jj) / ( ht_n(ji,jj) + 1. - ssmask(ji, jj) ) 234 tsa(ji,jj,:,jp_tem) = tsa(ji,jj,:,jp_tem) + tsn(ji,jj,:,jp_tem) * ztim235 tsa(ji,jj,:,jp_sal) = tsa(ji,jj,:,jp_sal) + tsn(ji,jj,:,jp_sal) * ztim236 pts_rhs(ji,jj,:,jp_tem) = pts_rhs(ji,jj,:,jp_tem) + ts(ji,jj,:,jp_tem,ktlev) * ztim 237 pts_rhs(ji,jj,:,jp_sal) = pts_rhs(ji,jj,:,jp_sal) + ts(ji,jj,:,jp_sal,ktlev) * ztim 236 238 END DO 237 239 END DO … … 250 252 DO jj = 2, jpj 251 253 DO ji = fs_2, fs_jpim1 252 zdep = 1._wp / e3t _n(ji,jj,jk)253 tsa(ji,jj,jk,jp_tem) = tsa(ji,jj,jk,jp_tem) - htsc_iscpl(ji,jj,jk,jp_tem) * zdep254 tsa(ji,jj,jk,jp_sal) = tsa(ji,jj,jk,jp_sal) - htsc_iscpl(ji,jj,jk,jp_sal) * zdep254 zdep = 1._wp / e3t(ji,jj,jk,ktlev) 255 pts_rhs(ji,jj,jk,jp_tem) = pts_rhs(ji,jj,jk,jp_tem) - htsc_iscpl(ji,jj,jk,jp_tem) * zdep 256 pts_rhs(ji,jj,jk,jp_sal) = pts_rhs(ji,jj,jk,jp_sal) - htsc_iscpl(ji,jj,jk,jp_sal) * zdep 255 257 END DO 256 258 END DO … … 259 261 260 262 IF( l_trdtra ) THEN ! save the horizontal diffusive trends for further diagnostics 261 ztrdt(:,:,:) = tsa(:,:,:,jp_tem) - ztrdt(:,:,:)262 ztrds(:,:,:) = tsa(:,:,:,jp_sal) - ztrds(:,:,:)263 ztrdt(:,:,:) = pts_rhs(:,:,:,jp_tem) - ztrdt(:,:,:) 264 ztrds(:,:,:) = pts_rhs(:,:,:,jp_sal) - ztrds(:,:,:) 263 265 CALL trd_tra( kt, 'TRA', jp_tem, jptra_nsr, ztrdt ) 264 266 CALL trd_tra( kt, 'TRA', jp_sal, jptra_nsr, ztrds ) -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/oce.F90
r10802 r10806 22 22 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:), TARGET :: ww !: vertical velocity [m/s] 23 23 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: wi !: vertical vel. (adaptive-implicit) [m/s] 24 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hdivn!: horizontal divergence [s-1]24 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:), TARGET :: hdiv !: horizontal divergence [s-1] 25 25 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:,:), TARGET :: ts !: 4D T-S fields [Celsius,psu] 26 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,: ) :: rab_b, rab_n!: thermal/haline expansion coef. [Celsius-1,psu-1]27 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,: ) :: rn2b , rn2!: brunt-vaisala frequency**2 [s-2]26 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:,:), TARGET :: r_ab !: thermal/haline expansion coef. [Celsius-1,psu-1] 27 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:,:),TARGET :: r_n2 !: brunt-vaisala frequency**2 [s-2] 28 28 ! 29 29 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: rhd !: in situ density anomalie rhd=(rho-rau0)/rau0 [no units] … … 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(:,:,:) :: hdivn !: horizontal divergence [s-1] 75 REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:,:,:) :: rab_b, rab_n !: thermal/haline expansion coef. [Celsius-1,psu-1] 76 REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:,:) :: rn2b , rn2 !: brunt-vaisala frequency**2 [s-2] 74 77 REAL(wp), PUBLIC, POINTER, SAVE, DIMENSION(:,:,:,:) :: tsb , tsn , tsa !: 4D T-S fields [Celsius,psu] 75 78 !! TEMPORARY POINTERS FOR DEVELOPMENT ONLY … … 91 94 ierr(:) = 0 92 95 ALLOCATE( uu (jpi,jpj,jpk, jpt) , vv (jpi,jpj,jpk,jpt) , & 93 & ww (jpi,jpj,jpk) , hdiv n(jpi,jpj,jpk), &96 & ww (jpi,jpj,jpk) , hdiv(jpi,jpj,jpk) , & 94 97 & ts (jpi,jpj,jpk,jpts,jpt) , & 95 & rab_b(jpi,jpj,jpk,jpts) , rab_n(jpi,jpj,jpk,jpts) , & 96 & rn2b (jpi,jpj,jpk) , rn2 (jpi,jpj,jpk) , & 97 & rhd (jpi,jpj,jpk) , rhop (jpi,jpj,jpk) , STAT=ierr(1) ) 98 & r_ab (jpi,jpj,jpk,jpts,jpt), r_n2 (jpi,jpj,jpk,jpt), & 99 & rhd (jpi,jpj,jpk) , rhop (jpi,jpj,jpk) , STAT=ierr(1) ) 98 100 ! 99 101 ALLOCATE( sshb(jpi,jpj) , sshn(jpi,jpj) , ssha(jpi,jpj) , & -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/OCE/step.F90
r10802 r10806 178 178 !!jc: fs simplification 179 179 180 u a(:,:,:) = 0._wp ! set dynamics trends to zero181 v a(:,:,:) = 0._wp180 uu(:,:,:,Nrhs) = 0._wp ! set dynamics trends to zero 181 vv(:,:,:,Nrhs) = 0._wp 182 182 183 183 IF( lk_asminc .AND. ln_asmiau .AND. ln_dyninc ) & … … 190 190 CALL dyn_adv ( kstp, Nm1, Nnn, uu(:,:,:,Nrhs), vv(:,:,:,Nrhs) ) ! advection (vector or flux form) 191 191 CALL dyn_vor ( kstp, Nnn, uu(:,:,:,Nrhs), vv(:,:,:,Nrhs) ) ! vorticity term including Coriolis 192 CALL dyn_ldf ( kstp ) ! lateral mixing192 CALL dyn_ldf ( kstp, Nm1, Nnn, uu(:,:,:,Nrhs), vv(:,:,:,Nrhs) ) ! lateral mixing 193 193 IF( ln_zdfosm ) CALL dyn_osm ( kstp ) ! OSMOSIS non-local velocity fluxes 194 194 CALL dyn_hpg ( kstp ) ! horizontal gradient of Hydrostatic pressure … … 233 233 ! Active tracers 234 234 !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 235 ts a(:,:,:,:) = 0._wp ! set tracer trends to zero235 ts(:,:,:,:,Nrhs) = 0._wp ! set tracer trends to zero 236 236 237 237 IF( lk_asminc .AND. ln_asmiau .AND. & 238 238 & ln_trainc ) CALL tra_asm_inc ( kstp ) ! apply tracer assimilation increment 239 CALL tra_sbc ( kstp ) ! surface boundary condition240 IF( ln_traqsr ) CALL tra_qsr ( kstp ) ! penetrative solar radiation qsr241 IF( ln_trabbc ) CALL tra_bbc ( kstp ) ! bottom heat flux242 IF( ln_trabbl ) CALL tra_bbl ( kstp ) ! advective (and/or diffusive) bottom boundary layer scheme243 IF( ln_tradmp ) CALL tra_dmp ( kstp ) ! internal damping trends239 CALL tra_sbc ( kstp, Nnn, ts(:,:,:,:,Nrhs) ) ! surface boundary condition 240 IF( ln_traqsr ) CALL tra_qsr ( kstp, Nnn, Nnn_2lev, ts(:,:,:,:,Nrhs) ) ! penetrative solar radiation qsr 241 IF( ln_trabbc ) CALL tra_bbc ( kstp, Nnn, ts(:,:,:,:,Nrhs) ) ! bottom heat flux 242 IF( ln_trabbl ) CALL tra_bbl ( kstp, Nm1, Nnn, Nnn_2lev, ts(:,:,:,:,Nrhs) ) ! advective (and/or diffusive) bottom boundary layer scheme 243 IF( ln_tradmp ) CALL tra_dmp ( kstp, Nm1, Nnn_2lev, ts(:,:,:,:,Nrhs) ) ! internal damping trends 244 244 IF( ln_bdy ) CALL bdy_tra_dmp ( kstp ) ! bdy damping trends 245 245 #if defined key_agrif … … 247 247 & CALL Agrif_Sponge_tra ! tracers sponge 248 248 #endif 249 CALL tra_adv ( kstp, Nm1, Nnn, Np1, ts(:,:,:,:,Nrhs) ) ! horizontal & vertical advection249 CALL tra_adv ( kstp, Nm1, Nnn, Np1, Nnn_2lev, 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 ) & 252 252 & CALL osm_rst( kstp, 'WRITE' )! write OSMOSIS outputs + wn (so must do here) to restarts 253 CALL tra_ldf ( kstp ) ! lateral mixing253 CALL tra_ldf ( kstp, Nm1, Nnn, Nnn_2lev, ts(:,:,:,:,Nrhs) ) ! lateral mixing 254 254 255 255 !!gm : why CALL to dia_ptr has been moved here??? (use trends info?) … … 343 343 vb => vv(:,:,:,Nm1); vn => vv(:,:,:,Nnn); va => vv(:,:,:,Np1) 344 344 wn => ww(:,:,:) 345 hdivn => hdiv(:,:,:) 346 rab_b => r_ab (:,:,:,:,Nm1_2lev); rab_n => r_ab (:,:,:,:,Nnn_2lev) 347 rn2b => r_n2 (:,:,:,Nm1_2lev) ; rn2 => r_n2 (:,:,:,Nnn_2lev) 345 348 346 349 tsb => ts(:,:,:,:,Nm1); tsn => ts(:,:,:,:,Nnn); tsa => ts(:,:,:,:,Np1) … … 360 363 gde3w_n => gde3w(:,:,:) 361 364 362 ht_n => ht(:,:)363 hu_b => hu(:,:,Nm1); hu_n => hu(:,:,Nnn); hu_a => hu(:,:,Np1)364 hv_b => hv(:,:,Nm1); hv_n => hv(:,:,Nnn); hv_a => hv(:,:,Np1)365 r1_hu_b => r1_hu(:,:,Nm1); r1_hu_n => r1_hu(:,:,Nnn); r1_hu_a => r1_hu(:,:,Np1)366 r1_hv_b => r1_hv(:,:,Nm1); r1_hv_n => r1_hv(:,:,Nnn); r1_hv_a => r1_hv(:,:,Np1)367 368 365 END SUBROUTINE update_pointers 369 366 -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trcadv.F90
r10802 r10806 68 68 CONTAINS 69 69 70 SUBROUTINE trc_adv( kt, ktlev1, ktlev2, ktlev3 )70 SUBROUTINE trc_adv( kt, ktlev1, ktlev2, ktlev3, kt2lev ) 71 71 !!---------------------------------------------------------------------- 72 72 !! *** ROUTINE trc_adv *** … … 78 78 INTEGER, INTENT(in) :: kt ! ocean time-step index 79 79 INTEGER, INTENT(in) :: ktlev1, ktlev2, ktlev3 ! time level indices for source terms 80 INTEGER, INTENT(in) :: kt2lev ! time level index for 2-time-level source terms 80 81 ! 81 82 INTEGER :: jk ! dummy loop index … … 128 129 CALL tra_adv_fct( kt, nittrc000, ktlev1, ktlev2, ktlev3, 'TRC', r2dttrc, zun, zvn, zwn, trb, trn, tra, jptra, nn_fct_h, nn_fct_v ) 129 130 CASE ( np_MUS ) ! MUSCL 130 CALL tra_adv_mus( kt, nittrc000, ktlev2, 'TRC', r2dttrc, zun, zvn, zwn, trb, tra, jptra , ln_mus_ups )131 CALL tra_adv_mus( kt, nittrc000, ktlev2, kt2lev, 'TRC', r2dttrc, zun, zvn, zwn, trb, tra, jptra , ln_mus_ups ) 131 132 CASE ( np_UBS ) ! UBS 132 133 CALL tra_adv_ubs( kt, nittrc000, ktlev2, 'TRC', r2dttrc, zun, zvn, zwn, trb, trn, tra, jptra , nn_ubs_v ) -
NEMO/branches/2019/dev_r10721_KERNEL-02_Storkey_Coward_IMMERSE_first_steps/src/TOP/TRP/trctrp.F90
r10802 r10806 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, Nm1, Nnn, Np1 ) ! horizontal & vertical advection66 CALL trc_adv ( kt, Nm1, Nnn, Np1, Nnn_2lev ) ! 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.