Changeset 11993 for NEMO/trunk/src/OCE/DIA
- Timestamp:
- 2019-11-28T11:20:53+01:00 (4 years ago)
- Location:
- NEMO/trunk/src/OCE/DIA
- Files:
-
- 3 edited
-
diaar5.F90 (modified) (15 diffs)
-
diahth.F90 (modified) (4 diffs)
-
diaptr.F90 (modified) (17 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/OCE/DIA/diaar5.F90
r11989 r11993 71 71 INTEGER, INTENT( in ) :: kt ! ocean time-step index 72 72 ! 73 INTEGER :: ji, jj, jk , iks, ikb! dummy loop arguments74 REAL(wp) :: zvolssh, zvol, zssh_steric, zztmp, zarho, ztemp, zsal, zmass , zsst73 INTEGER :: ji, jj, jk ! dummy loop arguments 74 REAL(wp) :: zvolssh, zvol, zssh_steric, zztmp, zarho, ztemp, zsal, zmass 75 75 REAL(wp) :: zaw, zbw, zrw 76 76 ! 77 77 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zarea_ssh , zbotpres ! 2D workspace 78 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zpe , z2d! 2D workspace79 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zrhd , zrhop , ztpot! 3D workspace78 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zpe ! 2D workspace 79 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zrhd , zrhop ! 3D workspace 80 80 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: ztsn ! 4D workspace 81 81 … … 86 86 87 87 IF( l_ar5 ) THEN 88 ALLOCATE( zarea_ssh(jpi,jpj) , zbotpres(jpi,jpj), z2d(jpi,jpj) )88 ALLOCATE( zarea_ssh(jpi,jpj) , zbotpres(jpi,jpj) ) 89 89 ALLOCATE( zrhd(jpi,jpj,jpk) , zrhop(jpi,jpj,jpk) ) 90 90 ALLOCATE( ztsn(jpi,jpj,jpk,jpts) ) … … 92 92 ENDIF 93 93 ! 94 CALL iom_put( 'e2u' , e2u (:,:) )95 CALL iom_put( 'e1v' , e1v (:,:) )96 CALL iom_put( 'areacello', area(:,:) )97 !98 IF( iom_use( 'volcello' ) .OR. iom_use( 'masscello' ) ) THEN99 zrhd(:,:,jpk) = 0._wp ! ocean volume ; rhd is used as workspace100 DO jk = 1, jpkm1101 zrhd(:,:,jk) = area(:,:) * e3t_n(:,:,jk) * tmask(:,:,jk)102 END DO103 CALL iom_put( 'volcello' , zrhd(:,:,:) ) ! WARNING not consistent with CMIP DR where volcello is at ca. 2000104 CALL iom_put( 'masscello' , rau0 * e3t_n(:,:,:) * tmask(:,:,:) ) ! ocean mass105 ENDIF106 !107 IF( iom_use( 'e3tb' ) ) THEN ! bottom layer thickness108 DO jj = 1, jpj109 DO ji = 1, jpi110 ikb = mbkt(ji,jj)111 z2d(ji,jj) = e3t_n(ji,jj,ikb)112 END DO113 END DO114 CALL iom_put( 'e3tb', z2d )115 ENDIF116 !117 94 IF( iom_use( 'voltot' ) .OR. iom_use( 'sshtot' ) .OR. iom_use( 'sshdyn' ) ) THEN 118 95 ! ! total volume of liquid seawater 119 zvolssh = glob_sum( 'diaar5', zarea_ssh(:,:) ) 120 zvol = vol0 + zvolssh 96 zvolssh = SUM( zarea_ssh(:,:) ) 97 CALL mpp_sum( 'diaar5', zvolssh ) 98 zvol = vol0 + zvolssh 121 99 122 100 CALL iom_put( 'voltot', zvol ) … … 140 118 DO ji = 1, jpi 141 119 DO jj = 1, jpj 142 iks = mikt(ji,jj) 143 zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,iks) + riceload(ji,jj) 120 zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj) 144 121 END DO 145 122 END DO … … 152 129 END IF 153 130 ! 154 zarho = glob_sum( 'diaar5', area(:,:) * zbotpres(:,:) ) 131 zarho = SUM( area(:,:) * zbotpres(:,:) ) 132 CALL mpp_sum( 'diaar5', zarho ) 155 133 zssh_steric = - zarho / area_tot 156 134 CALL iom_put( 'sshthster', zssh_steric ) … … 169 147 DO ji = 1,jpi 170 148 DO jj = 1,jpj 171 iks = mikt(ji,jj) 172 zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,iks) + riceload(ji,jj) 149 zbotpres(ji,jj) = zbotpres(ji,jj) + sshn(ji,jj) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj) 173 150 END DO 174 151 END DO … … 178 155 END IF 179 156 ! 180 zarho = glob_sum( 'diaar5', area(:,:) * zbotpres(:,:) ) 157 zarho = SUM( area(:,:) * zbotpres(:,:) ) 158 CALL mpp_sum( 'diaar5', zarho ) 181 159 zssh_steric = - zarho / area_tot 182 160 CALL iom_put( 'sshsteric', zssh_steric ) 161 183 162 ! ! ocean bottom pressure 184 163 zztmp = rau0 * grav * 1.e-4_wp ! recover pressure from pressure anomaly and cover to dbar = 1.e4 Pa … … 189 168 190 169 IF( iom_use( 'masstot' ) .OR. iom_use( 'temptot' ) .OR. iom_use( 'saltot' ) ) THEN 191 192 ztsn(:,:,:,:) = 0._wp ! ztsn(:,:,1,jp_tem/sal) is used here as 2D Workspace for temperature & salinity193 DO jk = 1, jpkm1194 DO jj = 1, jpj195 DO ji = 1, jpi196 zztmp = area(ji,jj) * e3t_n(ji,jj,jk)197 ztsn(ji,jj,1,jp_tem) = ztsn(ji,jj,1,jp_tem) + zztmp * tsn(ji,jj,jk,jp_tem)198 ztsn(ji,jj,1,jp_sal) = ztsn(ji,jj,1,jp_sal) + zztmp * tsn(ji,jj,jk,jp_sal)199 ENDDO200 ENDDO201 ENDDO202 203 170 ! ! Mean density anomalie, temperature and salinity 171 ztemp = 0._wp 172 zsal = 0._wp 173 DO jk = 1, jpkm1 174 DO jj = 1, jpj 175 DO ji = 1, jpi 176 zztmp = area(ji,jj) * e3t_n(ji,jj,jk) 177 ztemp = ztemp + zztmp * tsn(ji,jj,jk,jp_tem) 178 zsal = zsal + zztmp * tsn(ji,jj,jk,jp_sal) 179 END DO 180 END DO 181 END DO 182 IF( ln_linssh ) THEN 204 183 IF( ln_isfcav ) THEN 205 184 DO ji = 1, jpi 206 185 DO jj = 1, jpj 207 iks = mikt(ji,jj) 208 ztsn(ji,jj,1,jp_tem) = ztsn(ji,jj,1,jp_tem) + zarea_ssh(ji,jj) * tsn(ji,jj,iks,jp_tem) 209 ztsn(ji,jj,1,jp_sal) = ztsn(ji,jj,1,jp_sal) + zarea_ssh(ji,jj) * tsn(ji,jj,iks,jp_sal) 186 ztemp = ztemp + zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_tem) 187 zsal = zsal + zarea_ssh(ji,jj) * tsn(ji,jj,mikt(ji,jj),jp_sal) 210 188 END DO 211 189 END DO 212 190 ELSE 213 zt sn(:,:,1,jp_tem) = ztsn(:,:,1,jp_tem) + zarea_ssh(:,:) * tsn(:,:,1,jp_tem)214 z tsn(:,:,1,jp_sal) = ztsn(:,:,1,jp_sal) + zarea_ssh(:,:) * tsn(:,:,1,jp_sal)191 ztemp = ztemp + SUM( zarea_ssh(:,:) * tsn(:,:,1,jp_tem) ) 192 zsal = zsal + SUM( zarea_ssh(:,:) * tsn(:,:,1,jp_sal) ) 215 193 END IF 216 194 ENDIF 217 ! 218 ztemp = glob_sum( 'diaar5', ztsn(:,:,1,jp_tem) ) 219 zsal = glob_sum( 'diaar5', ztsn(:,:,1,jp_sal) ) 220 zmass = rau0 * ( zarho + zvol ) 195 IF( lk_mpp ) THEN 196 CALL mpp_sum( 'diaar5', ztemp ) 197 CALL mpp_sum( 'diaar5', zsal ) 198 END IF 199 ! 200 zmass = rau0 * ( zarho + zvol ) ! total mass of liquid seawater 201 ztemp = ztemp / zvol ! potential temperature in liquid seawater 202 zsal = zsal / zvol ! Salinity of liquid seawater 221 203 ! 222 204 CALL iom_put( 'masstot', zmass ) 223 CALL iom_put( 'temptot', ztemp / zvol ) 224 CALL iom_put( 'saltot' , zsal / zvol ) 225 ! 226 ENDIF 227 228 IF( ln_teos10 ) THEN ! ! potential temperature (TEOS-10 case) 229 IF( iom_use( 'toce_pot') .OR. iom_use( 'temptot_pot' ) .OR. iom_use( 'sst_pot' ) & 230 .OR. iom_use( 'ssttot' ) .OR. iom_use( 'tosmint_pot' ) ) THEN 231 ! 232 ALLOCATE( ztpot(jpi,jpj,jpk) ) 233 ztpot(:,:,jpk) = 0._wp 234 ztpot(:,:,:) = eos_pt_from_ct( tsn(:,:,:,jp_tem), tsn(:,:,:,jp_sal) ) 235 ! 236 CALL iom_put( 'toce_pot', ztpot(:,:,:) ) ! potential temperature (TEOS-10 case) 237 CALL iom_put( 'sst_pot' , ztpot(:,:,1) ) ! surface temperature 238 ! 239 IF( iom_use( 'temptot_pot' ) ) THEN ! Output potential temperature in case we use TEOS-10 240 z2d(:,:) = 0._wp 241 DO jk = 1, jpkm1 242 z2d(:,:) = z2d(:,:) + area(:,:) * e3t_n(:,:,jk) * ztpot(:,:,jk) 243 END DO 244 ztemp = glob_sum( 'diaar5', z2d(:,:) ) 245 CALL iom_put( 'temptot_pot', ztemp / zvol ) 246 ENDIF 247 ! 248 IF( iom_use( 'ssttot' ) ) THEN ! Output potential temperature in case we use TEOS-10 249 zsst = glob_sum( 'diaar5', area(:,:) * ztpot(:,:,1) ) 250 CALL iom_put( 'ssttot', zsst / area_tot ) 251 ENDIF 252 ! Vertical integral of temperature 253 IF( iom_use( 'tosmint_pot') ) THEN 254 z2d(:,:) = 0._wp 255 DO jk = 1, jpkm1 256 DO jj = 1, jpj 257 DO ji = 1, jpi ! vector opt. 258 z2d(ji,jj) = z2d(ji,jj) + rau0 * e3t_n(ji,jj,jk) * ztpot(ji,jj,jk) 259 END DO 260 END DO 261 END DO 262 CALL iom_put( 'tosmint_pot', z2d ) 263 ENDIF 264 DEALLOCATE( ztpot ) 265 ENDIF 266 ELSE 267 IF( iom_use('ssttot') ) THEN ! Output sst in case we use EOS-80 268 zsst = glob_sum( 'diaar5', area(:,:) * tsn(:,:,1,jp_tem) ) 269 CALL iom_put('ssttot', zsst / area_tot ) 270 ENDIF 205 CALL iom_put( 'temptot', ztemp ) 206 CALL iom_put( 'saltot' , zsal ) 207 ! 271 208 ENDIF 272 209 273 210 IF( iom_use( 'tnpeo' )) THEN 274 275 276 211 ! Work done against stratification by vertical mixing 212 ! Exclude points where rn2 is negative as convection kicks in here and 213 ! work is not being done against stratification 277 214 ALLOCATE( zpe(jpi,jpj) ) 278 215 zpe(:,:) = 0._wp … … 282 219 DO ji = 1, jpi 283 220 IF( rn2(ji,jj,jk) > 0._wp ) THEN 284 zrw = ( gdept_n(ji,jj,jk) - gdepw_n(ji,jj,jk) ) / e3w_n(ji,jj,jk) 221 zrw = ( gdepw_n(ji,jj,jk ) - gdept_n(ji,jj,jk) ) & 222 & / ( gdept_n(ji,jj,jk-1) - gdept_n(ji,jj,jk) ) 223 !!gm this can be reduced to : (depw-dept) / e3w (NB idem dans bn2 !) 224 ! zrw = ( gdept_n(ji,jj,jk) - gdepw_n(ji,jj,jk) ) / e3w_n(ji,jj,jk) 225 !!gm end 285 226 ! 286 227 zaw = rab_n(ji,jj,jk,jp_tem) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_tem)* zrw 287 228 zbw = rab_n(ji,jj,jk,jp_sal) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_sal)* zrw 288 229 ! 289 zpe(ji, jj) = zpe(ji, jj)&230 zpe(ji, jj) = zpe(ji, jj) & 290 231 & - grav * ( avt(ji,jj,jk) * zaw * (tsn(ji,jj,jk-1,jp_tem) - tsn(ji,jj,jk,jp_tem) ) & 291 232 & - avs(ji,jj,jk) * zbw * (tsn(ji,jj,jk-1,jp_sal) - tsn(ji,jj,jk,jp_sal) ) ) … … 298 239 DO ji = 1, jpi 299 240 DO jj = 1, jpj 300 zpe(ji,jj) = zpe(ji,jj) + avt(ji, jj,jk) * MIN(0._wp,rn2(ji,jj,jk)) * rau0 * e3w_n(ji,jj,jk)241 zpe(ji,jj) = zpe(ji,jj) + avt(ji, jj, jk) * MIN(0._wp,rn2(ji, jj, jk)) * rau0 * e3w_n(ji, jj, jk) 301 242 END DO 302 243 END DO 303 244 END DO 304 245 ENDIF 246 !!gm useless lbc_lnk since the computation above is performed over 1:jpi & 1:jpj 247 !!gm CALL lbc_lnk( 'diaar5', zpe, 'T', 1._wp) 305 248 CALL iom_put( 'tnpeo', zpe ) 306 249 DEALLOCATE( zpe ) … … 308 251 309 252 IF( l_ar5 ) THEN 310 DEALLOCATE( zarea_ssh , zbotpres , z2d)253 DEALLOCATE( zarea_ssh , zbotpres ) 311 254 DEALLOCATE( zrhd , zrhop ) 312 255 DEALLOCATE( ztsn ) … … 344 287 CALL lbc_lnk( 'diaar5', z2d, 'U', -1. ) 345 288 IF( cptr == 'adv' ) THEN 346 IF( ktra == jp_tem ) CALL iom_put( 'uadv_heattr', rau0_rcp * z2d ) ! advective heat transport in i-direction347 IF( ktra == jp_sal ) CALL iom_put( 'uadv_salttr', rau0 * z2d ) ! advective salt transport in i-direction289 IF( ktra == jp_tem ) CALL iom_put( "uadv_heattr" , rau0_rcp * z2d ) ! advective heat transport in i-direction 290 IF( ktra == jp_sal ) CALL iom_put( "uadv_salttr" , rau0 * z2d ) ! advective salt transport in i-direction 348 291 ENDIF 349 292 IF( cptr == 'ldf' ) THEN 350 IF( ktra == jp_tem ) CALL iom_put( 'udiff_heattr', rau0_rcp * z2d ) ! diffusive heat transport in i-direction351 IF( ktra == jp_sal ) CALL iom_put( 'udiff_salttr', rau0 * z2d ) ! diffusive salt transport in i-direction293 IF( ktra == jp_tem ) CALL iom_put( "udiff_heattr" , rau0_rcp * z2d ) ! diffusive heat transport in i-direction 294 IF( ktra == jp_sal ) CALL iom_put( "udiff_salttr" , rau0 * z2d ) ! diffusive salt transport in i-direction 352 295 ENDIF 353 296 ! … … 362 305 CALL lbc_lnk( 'diaar5', z2d, 'V', -1. ) 363 306 IF( cptr == 'adv' ) THEN 364 IF( ktra == jp_tem ) CALL iom_put( 'vadv_heattr', rau0_rcp * z2d ) ! advective heat transport in j-direction365 IF( ktra == jp_sal ) CALL iom_put( 'vadv_salttr', rau0 * z2d ) ! advective salt transport in j-direction307 IF( ktra == jp_tem ) CALL iom_put( "vadv_heattr" , rau0_rcp * z2d ) ! advective heat transport in j-direction 308 IF( ktra == jp_sal ) CALL iom_put( "vadv_salttr" , rau0 * z2d ) ! advective salt transport in j-direction 366 309 ENDIF 367 310 IF( cptr == 'ldf' ) THEN 368 IF( ktra == jp_tem ) CALL iom_put( 'vdiff_heattr', rau0_rcp * z2d ) ! diffusive heat transport in j-direction369 IF( ktra == jp_sal ) CALL iom_put( 'vdiff_salttr', rau0 * z2d ) ! diffusive salt transport in j-direction311 IF( ktra == jp_tem ) CALL iom_put( "vdiff_heattr" , rau0_rcp * z2d ) ! diffusive heat transport in j-direction 312 IF( ktra == jp_sal ) CALL iom_put( "vdiff_salttr" , rau0 * z2d ) ! diffusive salt transport in j-direction 370 313 ENDIF 371 314 … … 380 323 !!---------------------------------------------------------------------- 381 324 INTEGER :: inum 382 INTEGER :: ik , idep325 INTEGER :: ik 383 326 INTEGER :: ji, jj, jk ! dummy loop indices 384 327 REAL(wp) :: zztmp 385 328 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: zsaldta ! Jan/Dec levitus salinity 386 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zvol0387 329 ! 388 330 !!---------------------------------------------------------------------- … … 398 340 IF( dia_ar5_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'dia_ar5_init : unable to allocate arrays' ) 399 341 400 area(:,:) = e1e2t(:,:) 401 area_tot = glob_sum( 'diaar5', area(:,:) ) 402 403 ALLOCATE( zvol0(jpi,jpj) ) 404 zvol0 (:,:)= 0._wp342 area(:,:) = e1e2t(:,:) * tmask_i(:,:) 343 344 area_tot = SUM( area(:,:) ) ; CALL mpp_sum( 'diaar5', area_tot ) 345 346 vol0 = 0._wp 405 347 thick0(:,:) = 0._wp 406 348 DO jk = 1, jpkm1 407 DO jj = 1, jpj ! interpolation of salinity at the last ocean level (i.e. the partial step) 408 DO ji = 1, jpi 409 idep = tmask(ji,jj,jk) * e3t_0(ji,jj,jk) 410 zvol0 (ji,jj) = zvol0 (ji,jj) + idep * area(ji,jj) 411 thick0(ji,jj) = thick0(ji,jj) + idep 412 END DO 413 END DO 414 END DO 415 vol0 = glob_sum( 'diaar5', zvol0 ) 416 DEALLOCATE( zvol0 ) 349 vol0 = vol0 + SUM( area (:,:) * tmask(:,:,jk) * e3t_0(:,:,jk) ) 350 thick0(:,:) = thick0(:,:) + tmask_i(:,:) * tmask(:,:,jk) * e3t_0(:,:,jk) 351 END DO 352 CALL mpp_sum( 'diaar5', vol0 ) 417 353 418 354 IF( iom_use( 'sshthster' ) ) THEN 419 ALLOCATE( zsaldta(jpi,jpj,jp k,jpts) )355 ALLOCATE( zsaldta(jpi,jpj,jpj,jpts) ) 420 356 CALL iom_open ( 'sali_ref_clim_monthly', inum ) 421 357 CALL iom_get ( inum, jpdom_data, 'vosaline' , zsaldta(:,:,:,1), 1 ) -
NEMO/trunk/src/OCE/DIA/diahth.F90
r11989 r11993 11 11 !! 3.2 ! 2009-07 (S. Masson) hc300 bugfix + cleaning + add new diag 12 12 !!---------------------------------------------------------------------- 13 #if defined key_diahth 14 !!---------------------------------------------------------------------- 15 !! 'key_diahth' : thermocline depth diag. 16 !!---------------------------------------------------------------------- 13 17 !! dia_hth : Compute varius diagnostics associated with the mixed layer 14 18 !!---------------------------------------------------------------------- … … 28 32 PUBLIC dia_hth_alloc ! routine called by nemogcm.F90 29 33 30 LOGICAL , SAVE :: l_hth!: thermocline-20d depths flag34 LOGICAL , PUBLIC, PARAMETER :: lk_diahth = .TRUE. !: thermocline-20d depths flag 31 35 32 36 ! note: following variables should move to local variables once iom_put is always used 33 37 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hth !: depth of the max vertical temperature gradient [m] 34 38 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hd20 !: depth of 20 C isotherm [m] 35 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hd26 !: depth of 26 C isotherm [m]36 39 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hd28 !: depth of 28 C isotherm [m] 37 40 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: htc3 !: heat content of first 300 m [W] 38 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: htc7 !: heat content of first 700 m [W]39 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: htc20 !: heat content of first 2000 m [W]40 41 41 42 42 !!---------------------------------------------------------------------- … … 52 52 !!--------------------------------------------------------------------- 53 53 ! 54 ALLOCATE( hth(jpi,jpj), hd20(jpi,jpj), hd26(jpi,jpj), hd28(jpi,jpj), & 55 & htc3(jpi,jpj), htc7(jpi,jpj), htc20(jpi,jpj), STAT=dia_hth_alloc ) 54 ALLOCATE( hth(jpi,jpj), hd20(jpi,jpj), hd28(jpi,jpj), htc3(jpi,jpj), STAT=dia_hth_alloc ) 56 55 ! 57 56 CALL mpp_sum ( 'diahth', dia_hth_alloc ) … … 83 82 INTEGER, INTENT( in ) :: kt ! ocean time-step index 84 83 !! 85 INTEGER :: ji, jj, jk ! dummy loop arguments 86 REAL(wp) :: zrho3 = 0.03_wp ! density criterion for mixed layer depth 87 REAL(wp) :: zrho1 = 0.01_wp ! density criterion for mixed layer depth 88 REAL(wp) :: ztem2 = 0.2_wp ! temperature criterion for mixed layer depth 89 REAL(wp) :: zztmp, zzdep ! temporary scalars inside do loop 90 REAL(wp) :: zu, zv, zw, zut, zvt ! temporary workspace 91 REAL(wp), DIMENSION(jpi,jpj) :: zabs2 ! MLD: abs( tn - tn(10m) ) = ztem2 92 REAL(wp), DIMENSION(jpi,jpj) :: ztm2 ! Top of thermocline: tn = tn(10m) - ztem2 93 REAL(wp), DIMENSION(jpi,jpj) :: zrho10_3 ! MLD: rho = rho10m + zrho3 94 REAL(wp), DIMENSION(jpi,jpj) :: zpycn ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC) 95 REAL(wp), DIMENSION(jpi,jpj) :: ztinv ! max of temperature inversion 96 REAL(wp), DIMENSION(jpi,jpj) :: zdepinv ! depth of temperature inversion 97 REAL(wp), DIMENSION(jpi,jpj) :: zrho0_3 ! MLD rho = rho(surf) = 0.03 98 REAL(wp), DIMENSION(jpi,jpj) :: zrho0_1 ! MLD rho = rho(surf) = 0.01 99 REAL(wp), DIMENSION(jpi,jpj) :: zmaxdzT ! max of dT/dz 100 REAL(wp), DIMENSION(jpi,jpj) :: zdelr ! delta rho equivalent to deltaT = 0.2 84 INTEGER :: ji, jj, jk ! dummy loop arguments 85 INTEGER :: iid, ilevel ! temporary integers 86 INTEGER, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ik20, ik28 ! levels 87 REAL(wp) :: zavt5 = 5.e-4_wp ! Kz criterion for the turbocline depth 88 REAL(wp) :: zrho3 = 0.03_wp ! density criterion for mixed layer depth 89 REAL(wp) :: zrho1 = 0.01_wp ! density criterion for mixed layer depth 90 REAL(wp) :: ztem2 = 0.2_wp ! temperature criterion for mixed layer depth 91 REAL(wp) :: zthick_0, zcoef ! temporary scalars 92 REAL(wp) :: zztmp, zzdep ! temporary scalars inside do loop 93 REAL(wp) :: zu, zv, zw, zut, zvt ! temporary workspace 94 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zabs2 ! MLD: abs( tn - tn(10m) ) = ztem2 95 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ztm2 ! Top of thermocline: tn = tn(10m) - ztem2 96 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zrho10_3 ! MLD: rho = rho10m + zrho3 97 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zpycn ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC) 98 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: ztinv ! max of temperature inversion 99 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zdepinv ! depth of temperature inversion 100 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zrho0_3 ! MLD rho = rho(surf) = 0.03 101 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zrho0_1 ! MLD rho = rho(surf) = 0.01 102 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zmaxdzT ! max of dT/dz 103 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zthick ! vertical integration thickness 104 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zdelr ! delta rho equivalent to deltaT = 0.2 101 105 !!---------------------------------------------------------------------- 102 106 IF( ln_timing ) CALL timing_start('dia_hth') 103 107 104 108 IF( kt == nit000 ) THEN 105 l_hth = .FALSE.106 IF( iom_use( 'mlddzt' ) .OR. iom_use( 'mldr0_3' ) .OR. iom_use( 'mldr0_1' ) .OR. &107 & iom_use( 'mld_dt02' ) .OR. iom_use( 'topthdep' ) .OR. iom_use( 'mldr10_3' ) .OR. &108 & iom_use( '20d' ) .OR. iom_use( '26d' ) .OR. iom_use( '28d' ) .OR. &109 & iom_use( 'hc300' ) .OR. iom_use( 'hc700' ) .OR. iom_use( 'hc2000' ) .OR. &110 & iom_use( 'pycndep' ) .OR. iom_use( 'tinv' ) .OR. iom_use( 'depti' ) ) l_hth = .TRUE.111 109 ! ! allocate dia_hth array 112 IF( l_hth ) THEN 113 IF( dia_hth_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'dia_hth : unable to allocate standard arrays' ) 114 IF(lwp) WRITE(numout,*) 115 IF(lwp) WRITE(numout,*) 'dia_hth : diagnostics of the thermocline depth' 116 IF(lwp) WRITE(numout,*) '~~~~~~~ ' 117 IF(lwp) WRITE(numout,*) 118 ENDIF 110 IF( dia_hth_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'dia_hth : unable to allocate standard arrays' ) 111 112 IF(.NOT. ALLOCATED(ik20) ) THEN 113 ALLOCATE(ik20(jpi,jpj), ik28(jpi,jpj), & 114 & zabs2(jpi,jpj), & 115 & ztm2(jpi,jpj), & 116 & zrho10_3(jpi,jpj),& 117 & zpycn(jpi,jpj), & 118 & ztinv(jpi,jpj), & 119 & zdepinv(jpi,jpj), & 120 & zrho0_3(jpi,jpj), & 121 & zrho0_1(jpi,jpj), & 122 & zmaxdzT(jpi,jpj), & 123 & zthick(jpi,jpj), & 124 & zdelr(jpi,jpj), STAT=ji) 125 CALL mpp_sum('diahth', ji) 126 IF( ji /= 0 ) CALL ctl_stop( 'STOP', 'dia_hth : unable to allocate standard ocean arrays' ) 127 END IF 128 129 IF(lwp) WRITE(numout,*) 130 IF(lwp) WRITE(numout,*) 'dia_hth : diagnostics of the thermocline depth' 131 IF(lwp) WRITE(numout,*) '~~~~~~~ ' 132 IF(lwp) WRITE(numout,*) 119 133 ENDIF 120 134 121 IF( l_hth ) THEN 122 ! 123 IF( iom_use( 'mlddzt' ) .OR. iom_use( 'mldr0_3' ) .OR. iom_use( 'mldr0_1' ) ) THEN 124 ! initialization 125 ztinv (:,:) = 0._wp 126 zdepinv(:,:) = 0._wp 127 zmaxdzT(:,:) = 0._wp 128 DO jj = 1, jpj 129 DO ji = 1, jpi 130 zztmp = gdepw_n(ji,jj,mbkt(ji,jj)+1) 131 hth (ji,jj) = zztmp 132 zabs2 (ji,jj) = zztmp 133 ztm2 (ji,jj) = zztmp 134 zrho10_3(ji,jj) = zztmp 135 zpycn (ji,jj) = zztmp 136 END DO 135 ! initialization 136 ztinv (:,:) = 0._wp 137 zdepinv(:,:) = 0._wp 138 zmaxdzT(:,:) = 0._wp 139 DO jj = 1, jpj 140 DO ji = 1, jpi 141 zztmp = gdepw_n(ji,jj,mbkt(ji,jj)+1) 142 hth (ji,jj) = zztmp 143 zabs2 (ji,jj) = zztmp 144 ztm2 (ji,jj) = zztmp 145 zrho10_3(ji,jj) = zztmp 146 zpycn (ji,jj) = zztmp 147 END DO 148 END DO 149 IF( nla10 > 1 ) THEN 150 DO jj = 1, jpj 151 DO ji = 1, jpi 152 zztmp = gdepw_n(ji,jj,mbkt(ji,jj)+1) 153 zrho0_3(ji,jj) = zztmp 154 zrho0_1(ji,jj) = zztmp 137 155 END DO 138 IF( nla10 > 1 ) THEN 139 DO jj = 1, jpj 140 DO ji = 1, jpi 141 zztmp = gdepw_n(ji,jj,mbkt(ji,jj)+1) 142 zrho0_3(ji,jj) = zztmp 143 zrho0_1(ji,jj) = zztmp 144 END DO 145 END DO 156 END DO 157 ENDIF 158 159 ! Preliminary computation 160 ! computation of zdelr = (dr/dT)(T,S,10m)*(-0.2 degC) 161 DO jj = 1, jpj 162 DO ji = 1, jpi 163 IF( tmask(ji,jj,nla10) == 1. ) THEN 164 zu = 1779.50 + 11.250 * tsn(ji,jj,nla10,jp_tem) - 3.80 * tsn(ji,jj,nla10,jp_sal) & 165 & - 0.0745 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_tem) & 166 & - 0.0100 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_sal) 167 zv = 5891.00 + 38.000 * tsn(ji,jj,nla10,jp_tem) + 3.00 * tsn(ji,jj,nla10,jp_sal) & 168 & - 0.3750 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_tem) 169 zut = 11.25 - 0.149 * tsn(ji,jj,nla10,jp_tem) - 0.01 * tsn(ji,jj,nla10,jp_sal) 170 zvt = 38.00 - 0.750 * tsn(ji,jj,nla10,jp_tem) 171 zw = (zu + 0.698*zv) * (zu + 0.698*zv) 172 zdelr(ji,jj) = ztem2 * (1000.*(zut*zv - zvt*zu)/zw) 173 ELSE 174 zdelr(ji,jj) = 0._wp 146 175 ENDIF 176 END DO 177 END DO 178 179 ! ------------------------------------------------------------- ! 180 ! thermocline depth: strongest vertical gradient of temperature ! 181 ! turbocline depth (mixing layer depth): avt = zavt5 ! 182 ! MLD: rho = rho(1) + zrho3 ! 183 ! MLD: rho = rho(1) + zrho1 ! 184 ! ------------------------------------------------------------- ! 185 DO jk = jpkm1, 2, -1 ! loop from bottom to 2 186 DO jj = 1, jpj 187 DO ji = 1, jpi 188 ! 189 zzdep = gdepw_n(ji,jj,jk) 190 zztmp = ( tsn(ji,jj,jk-1,jp_tem) - tsn(ji,jj,jk,jp_tem) ) / zzdep * tmask(ji,jj,jk) ! vertical gradient of temperature (dT/dz) 191 zzdep = zzdep * tmask(ji,jj,1) 192 193 IF( zztmp > zmaxdzT(ji,jj) ) THEN 194 zmaxdzT(ji,jj) = zztmp ; hth (ji,jj) = zzdep ! max and depth of dT/dz 195 ENDIF 196 197 IF( nla10 > 1 ) THEN 198 zztmp = rhop(ji,jj,jk) - rhop(ji,jj,1) ! delta rho(1) 199 IF( zztmp > zrho3 ) zrho0_3(ji,jj) = zzdep ! > 0.03 200 IF( zztmp > zrho1 ) zrho0_1(ji,jj) = zzdep ! > 0.01 201 ENDIF 202 203 END DO 204 END DO 205 END DO 147 206 148 ! Preliminary computation 149 ! computation of zdelr = (dr/dT)(T,S,10m)*(-0.2 degC) 150 DO jj = 1, jpj 151 DO ji = 1, jpi 152 IF( tmask(ji,jj,nla10) == 1. ) THEN 153 zu = 1779.50 + 11.250 * tsn(ji,jj,nla10,jp_tem) - 3.80 * tsn(ji,jj,nla10,jp_sal) & 154 & - 0.0745 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_tem) & 155 & - 0.0100 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_sal) 156 zv = 5891.00 + 38.000 * tsn(ji,jj,nla10,jp_tem) + 3.00 * tsn(ji,jj,nla10,jp_sal) & 157 & - 0.3750 * tsn(ji,jj,nla10,jp_tem) * tsn(ji,jj,nla10,jp_tem) 158 zut = 11.25 - 0.149 * tsn(ji,jj,nla10,jp_tem) - 0.01 * tsn(ji,jj,nla10,jp_sal) 159 zvt = 38.00 - 0.750 * tsn(ji,jj,nla10,jp_tem) 160 zw = (zu + 0.698*zv) * (zu + 0.698*zv) 161 zdelr(ji,jj) = ztem2 * (1000.*(zut*zv - zvt*zu)/zw) 162 ELSE 163 zdelr(ji,jj) = 0._wp 164 ENDIF 165 END DO 207 CALL iom_put( "mlddzt", hth ) ! depth of the thermocline 208 IF( nla10 > 1 ) THEN 209 CALL iom_put( "mldr0_3", zrho0_3 ) ! MLD delta rho(surf) = 0.03 210 CALL iom_put( "mldr0_1", zrho0_1 ) ! MLD delta rho(surf) = 0.01 211 ENDIF 212 213 ! ------------------------------------------------------------- ! 214 ! MLD: abs( tn - tn(10m) ) = ztem2 ! 215 ! Top of thermocline: tn = tn(10m) - ztem2 ! 216 ! MLD: rho = rho10m + zrho3 ! 217 ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC) ! 218 ! temperature inversion: max( 0, max of tn - tn(10m) ) ! 219 ! depth of temperature inversion ! 220 ! ------------------------------------------------------------- ! 221 DO jk = jpkm1, nlb10, -1 ! loop from bottom to nlb10 222 DO jj = 1, jpj 223 DO ji = 1, jpi 224 ! 225 zzdep = gdepw_n(ji,jj,jk) * tmask(ji,jj,1) 226 ! 227 zztmp = tsn(ji,jj,nla10,jp_tem) - tsn(ji,jj,jk,jp_tem) ! - delta T(10m) 228 IF( ABS(zztmp) > ztem2 ) zabs2 (ji,jj) = zzdep ! abs > 0.2 229 IF( zztmp > ztem2 ) ztm2 (ji,jj) = zzdep ! > 0.2 230 zztmp = -zztmp ! delta T(10m) 231 IF( zztmp > ztinv(ji,jj) ) THEN ! temperature inversion 232 ztinv(ji,jj) = zztmp ; zdepinv (ji,jj) = zzdep ! max value and depth 233 ENDIF 234 235 zztmp = rhop(ji,jj,jk) - rhop(ji,jj,nla10) ! delta rho(10m) 236 IF( zztmp > zrho3 ) zrho10_3(ji,jj) = zzdep ! > 0.03 237 IF( zztmp > zdelr(ji,jj) ) zpycn (ji,jj) = zzdep ! > equi. delta T(10m) - 0.2 238 ! 166 239 END DO 167 168 ! ------------------------------------------------------------- ! 169 ! thermocline depth: strongest vertical gradient of temperature ! 170 ! turbocline depth (mixing layer depth): avt = zavt5 ! 171 ! MLD: rho = rho(1) + zrho3 ! 172 ! MLD: rho = rho(1) + zrho1 ! 173 ! ------------------------------------------------------------- ! 174 DO jk = jpkm1, 2, -1 ! loop from bottom to 2 175 DO jj = 1, jpj 176 DO ji = 1, jpi 177 ! 178 zzdep = gdepw_n(ji,jj,jk) 179 zztmp = ( tsn(ji,jj,jk-1,jp_tem) - tsn(ji,jj,jk,jp_tem) ) & 180 & / zzdep * tmask(ji,jj,jk) ! vertical gradient of temperature (dT/dz) 181 zzdep = zzdep * tmask(ji,jj,1) 182 183 IF( zztmp > zmaxdzT(ji,jj) ) THEN 184 zmaxdzT(ji,jj) = zztmp 185 hth (ji,jj) = zzdep ! max and depth of dT/dz 186 ENDIF 187 188 IF( nla10 > 1 ) THEN 189 zztmp = rhop(ji,jj,jk) - rhop(ji,jj,1) ! delta rho(1) 190 IF( zztmp > zrho3 ) zrho0_3(ji,jj) = zzdep ! > 0.03 191 IF( zztmp > zrho1 ) zrho0_1(ji,jj) = zzdep ! > 0.01 192 ENDIF 193 END DO 194 END DO 195 END DO 196 197 CALL iom_put( 'mlddzt', hth ) ! depth of the thermocline 198 IF( nla10 > 1 ) THEN 199 CALL iom_put( 'mldr0_3', zrho0_3 ) ! MLD delta rho(surf) = 0.03 200 CALL iom_put( 'mldr0_1', zrho0_1 ) ! MLD delta rho(surf) = 0.01 201 ENDIF 202 ! 203 ENDIF 204 ! 205 IF( iom_use( 'mld_dt02' ) .OR. iom_use( 'topthdep' ) .OR. iom_use( 'mldr10_3' ) .OR. & 206 & iom_use( 'pycndep' ) .OR. iom_use( 'tinv' ) .OR. iom_use( 'depti' ) ) THEN 207 ! ------------------------------------------------------------- ! 208 ! MLD: abs( tn - tn(10m) ) = ztem2 ! 209 ! Top of thermocline: tn = tn(10m) - ztem2 ! 210 ! MLD: rho = rho10m + zrho3 ! 211 ! pycnocline: rho = rho10m + (dr/dT)(T,S,10m)*(-0.2 degC) ! 212 ! temperature inversion: max( 0, max of tn - tn(10m) ) ! 213 ! depth of temperature inversion ! 214 ! ------------------------------------------------------------- ! 215 DO jk = jpkm1, nlb10, -1 ! loop from bottom to nlb10 216 DO jj = 1, jpj 217 DO ji = 1, jpi 218 ! 219 zzdep = gdepw_n(ji,jj,jk) * tmask(ji,jj,1) 220 ! 221 zztmp = tsn(ji,jj,nla10,jp_tem) - tsn(ji,jj,jk,jp_tem) ! - delta T(10m) 222 IF( ABS(zztmp) > ztem2 ) zabs2 (ji,jj) = zzdep ! abs > 0.2 223 IF( zztmp > ztem2 ) ztm2 (ji,jj) = zzdep ! > 0.2 224 zztmp = -zztmp ! delta T(10m) 225 IF( zztmp > ztinv(ji,jj) ) THEN ! temperature inversion 226 ztinv(ji,jj) = zztmp 227 zdepinv (ji,jj) = zzdep ! max value and depth 228 ENDIF 229 230 zztmp = rhop(ji,jj,jk) - rhop(ji,jj,nla10) ! delta rho(10m) 231 IF( zztmp > zrho3 ) zrho10_3(ji,jj) = zzdep ! > 0.03 232 IF( zztmp > zdelr(ji,jj) ) zpycn (ji,jj) = zzdep ! > equi. delta T(10m) - 0.2 233 ! 234 END DO 235 END DO 236 END DO 237 238 CALL iom_put( 'mld_dt02', zabs2 ) ! MLD abs(delta t) - 0.2 239 CALL iom_put( 'topthdep', ztm2 ) ! T(10) - 0.2 240 CALL iom_put( 'mldr10_3', zrho10_3 ) ! MLD delta rho(10m) = 0.03 241 CALL iom_put( 'pycndep' , zpycn ) ! MLD delta rho equi. delta T(10m) = 0.2 242 CALL iom_put( 'tinv' , ztinv ) ! max. temp. inv. (t10 ref) 243 CALL iom_put( 'depti' , zdepinv ) ! depth of max. temp. inv. (t10 ref) 244 ! 245 ENDIF 246 247 ! ------------------------------- ! 248 ! Depth of 20C/26C/28C isotherm ! 249 ! ------------------------------- ! 250 IF( iom_use ('20d') ) THEN ! depth of the 20 isotherm 251 ztem2 = 20. 252 CALL dia_hth_dep( ztem2, hd20 ) 253 CALL iom_put( '20d', hd20 ) 254 ENDIF 255 ! 256 IF( iom_use ('26d') ) THEN ! depth of the 26 isotherm 257 ztem2 = 26. 258 CALL dia_hth_dep( ztem2, hd26 ) 259 CALL iom_put( '26d', hd26 ) 260 ENDIF 261 ! 262 IF( iom_use ('28d') ) THEN ! depth of the 28 isotherm 263 ztem2 = 28. 264 CALL dia_hth_dep( ztem2, hd28 ) 265 CALL iom_put( '28d', hd28 ) 266 ENDIF 267 268 ! ----------------------------- ! 269 ! Heat content of first 300 m ! 270 ! ----------------------------- ! 271 IF( iom_use ('hc300') ) THEN 272 zzdep = 300. 273 CALL dia_hth_htc( zzdep, tsn(:,:,:,jp_tem), htc3 ) 274 CALL iom_put( 'hc300', rau0_rcp * htc3 ) ! vertically integrated heat content (J/m2) 275 ENDIF 276 ! 277 ! ----------------------------- ! 278 ! Heat content of first 700 m ! 279 ! ----------------------------- ! 280 IF( iom_use ('hc700') ) THEN 281 zzdep = 700. 282 CALL dia_hth_htc( zzdep, tsn(:,:,:,jp_tem), htc7 ) 283 CALL iom_put( 'hc700', rau0_rcp * htc7 ) ! vertically integrated heat content (J/m2) 284 285 ENDIF 286 ! 287 ! ----------------------------- ! 288 ! Heat content of first 2000 m ! 289 ! ----------------------------- ! 290 IF( iom_use ('hc2000') ) THEN 291 zzdep = 2000. 292 CALL dia_hth_htc( zzdep, tsn(:,:,:,jp_tem), htc20 ) 293 CALL iom_put( 'hc2000', rau0_rcp * htc20 ) ! vertically integrated heat content (J/m2) 294 ENDIF 295 ! 296 ENDIF 297 298 ! 299 IF( ln_timing ) CALL timing_stop('dia_hth') 300 ! 301 END SUBROUTINE dia_hth 302 303 SUBROUTINE dia_hth_dep( ptem, pdept ) 304 ! 305 REAL(wp), INTENT(in) :: ptem 306 REAL(wp), DIMENSION(jpi,jpj), INTENT(out) :: pdept 307 ! 308 INTEGER :: ji, jj, jk, iid 309 REAL(wp) :: zztmp, zzdep 310 INTEGER, DIMENSION(jpi,jpj) :: iktem 311 312 ! --------------------------------------- ! 313 ! search deepest level above ptem ! 314 ! --------------------------------------- ! 315 iktem(:,:) = 1 240 END DO 241 END DO 242 243 CALL iom_put( "mld_dt02", zabs2 ) ! MLD abs(delta t) - 0.2 244 CALL iom_put( "topthdep", ztm2 ) ! T(10) - 0.2 245 CALL iom_put( "mldr10_3", zrho10_3 ) ! MLD delta rho(10m) = 0.03 246 CALL iom_put( "pycndep" , zpycn ) ! MLD delta rho equi. delta T(10m) = 0.2 247 CALL iom_put( "tinv" , ztinv ) ! max. temp. inv. (t10 ref) 248 CALL iom_put( "depti" , zdepinv ) ! depth of max. temp. inv. (t10 ref) 249 250 251 ! ----------------------------------- ! 252 ! search deepest level above 20C/28C ! 253 ! ----------------------------------- ! 254 ik20(:,:) = 1 255 ik28(:,:) = 1 316 256 DO jk = 1, jpkm1 ! beware temperature is not always decreasing with depth => loop from top to bottom 317 257 DO jj = 1, jpj 318 258 DO ji = 1, jpi 319 259 zztmp = tsn(ji,jj,jk,jp_tem) 320 IF( zztmp >= ptem ) iktem(ji,jj) = jk 260 IF( zztmp >= 20. ) ik20(ji,jj) = jk 261 IF( zztmp >= 28. ) ik28(ji,jj) = jk 321 262 END DO 322 263 END DO 323 264 END DO 324 265 325 ! --------------------------- ----!326 ! Depth of ptem isotherm!327 ! --------------------------- ----!266 ! --------------------------- ! 267 ! Depth of 20C/28C isotherm ! 268 ! --------------------------- ! 328 269 DO jj = 1, jpj 329 270 DO ji = 1, jpi 330 271 ! 331 zzdep = gdepw_n(ji,jj,mbkt(ji,jj)+1) ! depth of the o cean bottom272 zzdep = gdepw_n(ji,jj,mbkt(ji,jj)+1) ! depth of the oean bottom 332 273 ! 333 iid = ik tem(ji,jj)274 iid = ik20(ji,jj) 334 275 IF( iid /= 1 ) THEN 335 zztmp =gdept_n(ji,jj,iid ) & ! linear interpolation276 zztmp = gdept_n(ji,jj,iid ) & ! linear interpolation 336 277 & + ( gdept_n(ji,jj,iid+1) - gdept_n(ji,jj,iid) ) & 337 278 & * ( 20.*tmask(ji,jj,iid+1) - tsn(ji,jj,iid,jp_tem) ) & 338 279 & / ( tsn(ji,jj,iid+1,jp_tem) - tsn(ji,jj,iid,jp_tem) + (1.-tmask(ji,jj,1)) ) 339 pdept(ji,jj) = MIN( zztmp , zzdep) * tmask(ji,jj,1) ! bound by the ocean depth280 hd20(ji,jj) = MIN( zztmp , zzdep) * tmask(ji,jj,1) ! bound by the ocean depth 340 281 ELSE 341 pdept(ji,jj) = 0._wp282 hd20(ji,jj) = 0._wp 342 283 ENDIF 343 END DO 344 END DO 345 ! 346 END SUBROUTINE dia_hth_dep 347 348 349 SUBROUTINE dia_hth_htc( pdep, ptn, phtc ) 350 ! 351 REAL(wp), INTENT(in) :: pdep ! depth over the heat content 352 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in) :: ptn 353 REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: phtc 354 ! 355 INTEGER :: ji, jj, jk, ik 356 REAL(wp), DIMENSION(jpi,jpj) :: zthick 357 INTEGER , DIMENSION(jpi,jpj) :: ilevel 358 359 284 ! 285 iid = ik28(ji,jj) 286 IF( iid /= 1 ) THEN 287 zztmp = gdept_n(ji,jj,iid ) & ! linear interpolation 288 & + ( gdept_n(ji,jj,iid+1) - gdept_n(ji,jj,iid) ) & 289 & * ( 28.*tmask(ji,jj,iid+1) - tsn(ji,jj,iid,jp_tem) ) & 290 & / ( tsn(ji,jj,iid+1,jp_tem) - tsn(ji,jj,iid,jp_tem) + (1.-tmask(ji,jj,1)) ) 291 hd28(ji,jj) = MIN( zztmp , zzdep ) * tmask(ji,jj,1) ! bound by the ocean depth 292 ELSE 293 hd28(ji,jj) = 0._wp 294 ENDIF 295 296 END DO 297 END DO 298 CALL iom_put( "20d", hd20 ) ! depth of the 20 isotherm 299 CALL iom_put( "28d", hd28 ) ! depth of the 28 isotherm 300 301 ! ----------------------------- ! 302 ! Heat content of first 300 m ! 303 ! ----------------------------- ! 304 305 ! find ilevel with (ilevel+1) the deepest W-level above 300m (we assume we can use e3t_1d to do this search...) 306 ilevel = 0 307 zthick_0 = 0._wp 308 DO jk = 1, jpkm1 309 zthick_0 = zthick_0 + e3t_1d(jk) 310 IF( zthick_0 < 300. ) ilevel = jk 311 END DO 360 312 ! surface boundary condition 361 362 IF( .NOT. ln_linssh ) THEN ; zthick(:,:) = 0._wp ; phtc(:,:) = 0._wp 363 ELSE ; zthick(:,:) = sshn(:,:) ; phtc(:,:) = ptn(:,:,1) * sshn(:,:) * tmask(:,:,1) 313 IF( ln_linssh ) THEN ; zthick(:,:) = sshn(:,:) ; htc3(:,:) = tsn(:,:,1,jp_tem) * sshn(:,:) * tmask(:,:,1) 314 ELSE ; zthick(:,:) = 0._wp ; htc3(:,:) = 0._wp 364 315 ENDIF 365 ! 366 ilevel(:,:) = 1 367 DO jk = 2, jpkm1 368 DO jj = 1, jpj 369 DO ji = 1, jpi 370 IF( ( gdept_n(ji,jj,jk) < pdep ) .AND. ( tmask(ji,jj,jk) == 1 ) ) THEN 371 ilevel(ji,jj) = jk 372 zthick(ji,jj) = zthick(ji,jj) + e3t_n(ji,jj,jk) 373 phtc (ji,jj) = phtc (ji,jj) + e3t_n(ji,jj,jk) * ptn(ji,jj,jk) 374 ENDIF 375 ENDDO 376 ENDDO 377 ENDDO 378 ! 316 ! integration down to ilevel 317 DO jk = 1, ilevel 318 zthick(:,:) = zthick(:,:) + e3t_n(:,:,jk) 319 htc3 (:,:) = htc3 (:,:) + e3t_n(:,:,jk) * tsn(:,:,jk,jp_tem) * tmask(:,:,jk) 320 END DO 321 ! deepest layer 322 zthick(:,:) = 300. - zthick(:,:) ! remaining thickness to reach 300m 379 323 DO jj = 1, jpj 380 324 DO ji = 1, jpi 381 ik = ilevel(ji,jj) 382 zthick(ji,jj) = pdep - zthick(ji,jj) ! remaining thickness to reach depht pdep 383 phtc(ji,jj) = phtc(ji,jj) + ptn(ji,jj,ik+1) * MIN( e3t_n(ji,jj,ik+1), zthick(ji,jj) ) & 384 * tmask(ji,jj,ik+1) 385 END DO 386 ENDDO 387 ! 388 ! 389 END SUBROUTINE dia_hth_htc 325 htc3(ji,jj) = htc3(ji,jj) + tsn(ji,jj,ilevel+1,jp_tem) & 326 & * MIN( e3t_n(ji,jj,ilevel+1), zthick(ji,jj) ) * tmask(ji,jj,ilevel+1) 327 END DO 328 END DO 329 ! from temperature to heat contain 330 zcoef = rau0 * rcp 331 htc3(:,:) = zcoef * htc3(:,:) 332 CALL iom_put( "hc300", htc3 ) ! first 300m heat content 333 ! 334 IF( ln_timing ) CALL timing_stop('dia_hth') 335 ! 336 END SUBROUTINE dia_hth 337 338 #else 339 !!---------------------------------------------------------------------- 340 !! Default option : Empty module 341 !!---------------------------------------------------------------------- 342 LOGICAL , PUBLIC, PARAMETER :: lk_diahth = .FALSE. !: thermocline-20d depths flag 343 CONTAINS 344 SUBROUTINE dia_hth( kt ) ! Empty routine 345 IMPLICIT NONE 346 INTEGER, INTENT( in ) :: kt 347 WRITE(*,*) 'dia_hth: You should not have seen this print! error?', kt 348 END SUBROUTINE dia_hth 349 #endif 390 350 391 351 !!====================================================================== -
NEMO/trunk/src/OCE/DIA/diaptr.F90
r11989 r11993 10 10 !! 3.6 ! 2014-12 (C. Ethe) use of IOM 11 11 !! 3.6 ! 2016-06 (T. Graham) Addition of diagnostics for CMIP6 12 !! 4.0 ! 2010-08 ( C. Ethe, J. Deshayes ) Improvment13 12 !!---------------------------------------------------------------------- 14 13 … … 43 42 44 43 ! !!** namelist namptr ** 45 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hstr_adv, hstr_ldf, hstr_eiv !: Heat/Salt TRansports(adv, diff, Bolus.) 46 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: hstr_ove, hstr_btr, hstr_vtr !: heat Salt TRansports(overturn, baro, merional) 47 48 LOGICAL , PUBLIC :: l_diaptr !: tracers trend flag (set from namelist in trdini) 49 INTEGER, PARAMETER, PUBLIC :: nptr = 5 ! (glo, atl, pac, ind, ipc) 44 REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) :: htr_adv, htr_ldf, htr_eiv !: Heat TRansports (adv, diff, Bolus.) 45 REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) :: str_adv, str_ldf, str_eiv !: Salt TRansports (adv, diff, Bolus.) 46 REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) :: htr_ove, str_ove !: heat Salt TRansports ( overturn.) 47 REAL(wp), ALLOCATABLE, SAVE, PUBLIC, DIMENSION(:,:) :: htr_btr, str_btr !: heat Salt TRansports ( barotropic ) 48 49 LOGICAL, PUBLIC :: ln_diaptr ! Poleward transport flag (T) or not (F) 50 LOGICAL, PUBLIC :: ln_subbas ! Atlantic/Pacific/Indian basins calculation 51 INTEGER, PUBLIC :: nptr ! = 1 (l_subbas=F) or = 5 (glo, atl, pac, ind, ipc) (l_subbas=T) 50 52 51 53 REAL(wp) :: rc_sv = 1.e-6_wp ! conversion from m3/s to Sverdrup 52 54 REAL(wp) :: rc_pwatt = 1.e-15_wp ! conversion from W to PW (further x rau0 x Cp) 53 REAL(wp) :: rc_ggram = 1.e- 9_wp ! conversion from g to Gg (further x rau0)54 55 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: btmsk ! T-point basin interior masks56 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: btmsk34 ! mask out Southern Ocean (=0 south of 34°S)57 58 REAL(wp), TARGET, ALLOCATABLE, SAVE, DIMENSION(:) :: p_fval1d 59 REAL(wp), TARGET, ALLOCATABLE, SAVE, DIMENSION(: ,:) :: p_fval2d60 61 LOGICAL :: ll_init = .TRUE. !: tracers trend flag (set from namelist in trdini) 55 REAL(wp) :: rc_ggram = 1.e-6_wp ! conversion from g to Pg 56 57 CHARACTER(len=3), ALLOCATABLE, SAVE, DIMENSION(:) :: clsubb 58 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: btmsk ! T-point basin interior masks 59 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: btm30 ! mask out Southern Ocean (=0 south of 30°S) 60 61 REAL(wp), TARGET, ALLOCATABLE, SAVE, DIMENSION(:) :: p_fval1d 62 REAL(wp), TARGET, ALLOCATABLE, SAVE, DIMENSION(:,:) :: p_fval2d 63 62 64 !! * Substitutions 63 65 # include "vectopt_loop_substitute.h90" … … 69 71 CONTAINS 70 72 71 SUBROUTINE dia_ptr( kt,pvtr )73 SUBROUTINE dia_ptr( pvtr ) 72 74 !!---------------------------------------------------------------------- 73 75 !! *** ROUTINE dia_ptr *** 74 76 !!---------------------------------------------------------------------- 75 INTEGER, INTENT( in ) :: kt ! ocean time-step index76 77 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in), OPTIONAL :: pvtr ! j-effective transport 77 78 ! … … 79 80 REAL(wp) :: zsfc,zvfc ! local scalar 80 81 REAL(wp), DIMENSION(jpi,jpj) :: z2d ! 2D workspace 82 REAL(wp), DIMENSION(jpi,jpj,jpk) :: z3d ! 3D workspace 81 83 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zmask ! 3D workspace 82 REAL(wp), DIMENSION(jpi,jpj,jpk) :: z3d ! 3D workspace83 84 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts) :: zts ! 3D workspace 84 REAL(wp), DIMENSION(jpj) :: zvsum, ztsum, zssum ! 1D workspace 85 REAL(wp), DIMENSION(jpj) :: vsum ! 1D workspace 86 REAL(wp), DIMENSION(jpj,jpts) :: tssum ! 1D workspace 87 85 88 ! 86 89 !overturning calculation 87 REAL(wp), DIMENSION(jpj,jpk,nptr) :: sjk, r1_sjk, v_msf ! i-mean i-k-surface and its inverse 88 REAL(wp), DIMENSION(jpj,jpk,nptr) :: zt_jk, zs_jk ! i-mean T and S, j-Stream-Function 89 90 REAL(wp), DIMENSION(jpi,jpj,jpk,nptr) :: z4d1, z4d2 91 REAL(wp), DIMENSION(jpi,jpj,nptr) :: z3dtr ! i-mean T and S, j-Stream-Function 90 REAL(wp), DIMENSION(jpj,jpk,nptr) :: sjk , r1_sjk ! i-mean i-k-surface and its inverse 91 REAL(wp), DIMENSION(jpj,jpk,nptr) :: v_msf, sn_jk , tn_jk ! i-mean T and S, j-Stream-Function 92 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zvn ! 3D workspace 93 94 95 CHARACTER( len = 12 ) :: cl1 92 96 !!---------------------------------------------------------------------- 93 97 ! 94 98 IF( ln_timing ) CALL timing_start('dia_ptr') 95 99 96 IF( kt == nit000 .AND. ll_init ) CALL dia_ptr_init 97 ! 98 IF( .NOT. l_diaptr ) RETURN 99 100 ! 100 101 IF( PRESENT( pvtr ) ) THEN 101 IF( iom_use( 'zomsf' ) ) THEN ! effective MSF 102 DO jn = 1, nptr ! by sub-basins 103 z4d1(1,:,:,jn) = ptr_sjk( pvtr(:,:,:), btmsk34(:,:,jn) ) ! zonal cumulative effective transport excluding closed seas 104 DO jk = jpkm1, 1, -1 105 z4d1(1,:,jk,jn) = z4d1(1,:,jk+1,jn) - z4d1(1,:,jk,jn) ! effective j-Stream-Function (MSF) 102 IF( iom_use("zomsfglo") ) THEN ! effective MSF 103 z3d(1,:,:) = ptr_sjk( pvtr(:,:,:) ) ! zonal cumulative effective transport 104 DO jk = 2, jpkm1 105 z3d(1,:,jk) = z3d(1,:,jk-1) + z3d(1,:,jk) ! effective j-Stream-Function (MSF) 106 END DO 107 DO ji = 1, jpi 108 z3d(ji,:,:) = z3d(1,:,:) 109 ENDDO 110 cl1 = TRIM('zomsf'//clsubb(1) ) 111 CALL iom_put( cl1, z3d * rc_sv ) 112 DO jn = 2, nptr ! by sub-basins 113 z3d(1,:,:) = ptr_sjk( pvtr(:,:,:), btmsk(:,:,jn)*btm30(:,:) ) 114 DO jk = 2, jpkm1 115 z3d(1,:,jk) = z3d(1,:,jk-1) + z3d(1,:,jk) ! effective j-Stream-Function (MSF) 106 116 END DO 107 117 DO ji = 1, jpi 108 z 4d1(ji,:,:,jn) = z4d1(1,:,:,jn)109 ENDDO 110 END DO111 CALL iom_put( 'zomsf', z4d1* rc_sv )112 ENDIF113 IF( iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) .OR. &114 & iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr') ) THEN118 z3d(ji,:,:) = z3d(1,:,:) 119 ENDDO 120 cl1 = TRIM('zomsf'//clsubb(jn) ) 121 CALL iom_put( cl1, z3d * rc_sv ) 122 END DO 123 ENDIF 124 IF( iom_use("sopstove") .OR. iom_use("sophtove") .OR. iom_use("sopstbtr") .OR. iom_use("sophtbtr") ) THEN 115 125 ! define fields multiplied by scalar 116 126 zmask(:,:,:) = 0._wp 117 127 zts(:,:,:,:) = 0._wp 128 zvn(:,:,:) = 0._wp 118 129 DO jk = 1, jpkm1 119 130 DO jj = 1, jpjm1 … … 123 134 zts(ji,jj,jk,jp_tem) = (tsn(ji,jj,jk,jp_tem)+tsn(ji,jj+1,jk,jp_tem)) * 0.5 * zvfc !Tracers averaged onto V grid 124 135 zts(ji,jj,jk,jp_sal) = (tsn(ji,jj,jk,jp_sal)+tsn(ji,jj+1,jk,jp_sal)) * 0.5 * zvfc 136 zvn(ji,jj,jk) = vn(ji,jj,jk) * zvfc 125 137 ENDDO 126 138 ENDDO 127 139 ENDDO 128 140 ENDIF 129 IF( iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) ) THEN 130 DO jn = 1, nptr 131 sjk(:,:,jn) = ptr_sjk( zmask(:,:,:), btmsk(:,:,jn) ) 132 r1_sjk(:,:,jn) = 0._wp 133 WHERE( sjk(:,:,jn) /= 0._wp ) r1_sjk(:,:,jn) = 1._wp / sjk(:,:,jn) 134 ! i-mean T and S, j-Stream-Function, basin 135 zt_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) * r1_sjk(:,:,jn) 136 zs_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) * r1_sjk(:,:,jn) 137 v_msf(:,:,jn) = ptr_sjk( pvtr(:,:,:), btmsk34(:,:,jn) ) 138 hstr_ove(:,jp_tem,jn) = SUM( v_msf(:,:,jn)*zt_jk(:,:,jn), 2 ) 139 hstr_ove(:,jp_sal,jn) = SUM( v_msf(:,:,jn)*zs_jk(:,:,jn), 2 ) 140 ! 141 ENDDO 142 DO jn = 1, nptr 143 z3dtr(1,:,jn) = hstr_ove(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) 144 DO ji = 1, jpi 145 z3dtr(ji,:,jn) = z3dtr(1,:,jn) 146 ENDDO 147 ENDDO 148 CALL iom_put( 'sophtove', z3dtr ) 149 DO jn = 1, nptr 150 z3dtr(1,:,jn) = hstr_ove(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) 151 DO ji = 1, jpi 152 z3dtr(ji,:,jn) = z3dtr(1,:,jn) 153 ENDDO 154 ENDDO 155 CALL iom_put( 'sopstove', z3dtr ) 156 ENDIF 157 158 IF( iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) ) THEN 159 ! Calculate barotropic heat and salt transport here 160 DO jn = 1, nptr 161 sjk(:,1,jn) = ptr_sj( zmask(:,:,:), btmsk(:,:,jn) ) 162 r1_sjk(:,1,jn) = 0._wp 163 WHERE( sjk(:,1,jn) /= 0._wp ) r1_sjk(:,1,jn) = 1._wp / sjk(:,1,jn) 164 ! 165 zvsum(:) = ptr_sj( pvtr(:,:,:), btmsk34(:,:,jn) ) 166 ztsum(:) = ptr_sj( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) 167 zssum(:) = ptr_sj( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) 168 hstr_btr(:,jp_tem,jn) = zvsum(:) * ztsum(:) * r1_sjk(:,1,jn) 169 hstr_btr(:,jp_sal,jn) = zvsum(:) * zssum(:) * r1_sjk(:,1,jn) 170 ! 171 ENDDO 172 DO jn = 1, nptr 173 z3dtr(1,:,jn) = hstr_btr(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) 174 DO ji = 1, jpi 175 z3dtr(ji,:,jn) = z3dtr(1,:,jn) 176 ENDDO 177 ENDDO 178 CALL iom_put( 'sophtbtr', z3dtr ) 179 DO jn = 1, nptr 180 z3dtr(1,:,jn) = hstr_btr(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) 181 DO ji = 1, jpi 182 z3dtr(ji,:,jn) = z3dtr(1,:,jn) 183 ENDDO 184 ENDDO 185 CALL iom_put( 'sopstbtr', z3dtr ) 186 ENDIF 141 IF( iom_use("sopstove") .OR. iom_use("sophtove") ) THEN 142 sjk(:,:,1) = ptr_sjk( zmask(:,:,:), btmsk(:,:,1) ) 143 r1_sjk(:,:,1) = 0._wp 144 WHERE( sjk(:,:,1) /= 0._wp ) r1_sjk(:,:,1) = 1._wp / sjk(:,:,1) 145 146 ! i-mean T and S, j-Stream-Function, global 147 tn_jk(:,:,1) = ptr_sjk( zts(:,:,:,jp_tem) ) * r1_sjk(:,:,1) 148 sn_jk(:,:,1) = ptr_sjk( zts(:,:,:,jp_sal) ) * r1_sjk(:,:,1) 149 v_msf(:,:,1) = ptr_sjk( zvn(:,:,:) ) 150 151 htr_ove(:,1) = SUM( v_msf(:,:,1)*tn_jk(:,:,1) ,2 ) 152 str_ove(:,1) = SUM( v_msf(:,:,1)*sn_jk(:,:,1) ,2 ) 153 154 z2d(1,:) = htr_ove(:,1) * rc_pwatt ! (conversion in PW) 155 DO ji = 1, jpi 156 z2d(ji,:) = z2d(1,:) 157 ENDDO 158 cl1 = 'sophtove' 159 CALL iom_put( TRIM(cl1), z2d ) 160 z2d(1,:) = str_ove(:,1) * rc_ggram ! (conversion in Gg) 161 DO ji = 1, jpi 162 z2d(ji,:) = z2d(1,:) 163 ENDDO 164 cl1 = 'sopstove' 165 CALL iom_put( TRIM(cl1), z2d ) 166 IF( ln_subbas ) THEN 167 DO jn = 2, nptr 168 sjk(:,:,jn) = ptr_sjk( zmask(:,:,:), btmsk(:,:,jn) ) 169 r1_sjk(:,:,jn) = 0._wp 170 WHERE( sjk(:,:,jn) /= 0._wp ) r1_sjk(:,:,jn) = 1._wp / sjk(:,:,jn) 171 172 ! i-mean T and S, j-Stream-Function, basin 173 tn_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) * r1_sjk(:,:,jn) 174 sn_jk(:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) * r1_sjk(:,:,jn) 175 v_msf(:,:,jn) = ptr_sjk( zvn(:,:,:), btmsk(:,:,jn) ) 176 htr_ove(:,jn) = SUM( v_msf(:,:,jn)*tn_jk(:,:,jn) ,2 ) 177 str_ove(:,jn) = SUM( v_msf(:,:,jn)*sn_jk(:,:,jn) ,2 ) 178 179 z2d(1,:) = htr_ove(:,jn) * rc_pwatt ! (conversion in PW) 180 DO ji = 1, jpi 181 z2d(ji,:) = z2d(1,:) 182 ENDDO 183 cl1 = TRIM('sophtove_'//clsubb(jn)) 184 CALL iom_put( cl1, z2d ) 185 z2d(1,:) = str_ove(:,jn) * rc_ggram ! (conversion in Gg) 186 DO ji = 1, jpi 187 z2d(ji,:) = z2d(1,:) 188 ENDDO 189 cl1 = TRIM('sopstove_'//clsubb(jn)) 190 CALL iom_put( cl1, z2d ) 191 END DO 192 ENDIF 193 ENDIF 194 IF( iom_use("sopstbtr") .OR. iom_use("sophtbtr") ) THEN 195 ! Calculate barotropic heat and salt transport here 196 sjk(:,1,1) = ptr_sj( zmask(:,:,:), btmsk(:,:,1) ) 197 r1_sjk(:,1,1) = 0._wp 198 WHERE( sjk(:,1,1) /= 0._wp ) r1_sjk(:,1,1) = 1._wp / sjk(:,1,1) 199 200 vsum = ptr_sj( zvn(:,:,:), btmsk(:,:,1)) 201 tssum(:,jp_tem) = ptr_sj( zts(:,:,:,jp_tem), btmsk(:,:,1) ) 202 tssum(:,jp_sal) = ptr_sj( zts(:,:,:,jp_sal), btmsk(:,:,1) ) 203 htr_btr(:,1) = vsum * tssum(:,jp_tem) * r1_sjk(:,1,1) 204 str_btr(:,1) = vsum * tssum(:,jp_sal) * r1_sjk(:,1,1) 205 z2d(1,:) = htr_btr(:,1) * rc_pwatt ! (conversion in PW) 206 DO ji = 2, jpi 207 z2d(ji,:) = z2d(1,:) 208 ENDDO 209 cl1 = 'sophtbtr' 210 CALL iom_put( TRIM(cl1), z2d ) 211 z2d(1,:) = str_btr(:,1) * rc_ggram ! (conversion in Gg) 212 DO ji = 2, jpi 213 z2d(ji,:) = z2d(1,:) 214 ENDDO 215 cl1 = 'sopstbtr' 216 CALL iom_put( TRIM(cl1), z2d ) 217 IF( ln_subbas ) THEN 218 DO jn = 2, nptr 219 sjk(:,1,jn) = ptr_sj( zmask(:,:,:), btmsk(:,:,jn) ) 220 r1_sjk(:,1,jn) = 0._wp 221 WHERE( sjk(:,1,jn) /= 0._wp ) r1_sjk(:,1,jn) = 1._wp / sjk(:,1,jn) 222 vsum = ptr_sj( zvn(:,:,:), btmsk(:,:,jn)) 223 tssum(:,jp_tem) = ptr_sj( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) 224 tssum(:,jp_sal) = ptr_sj( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) 225 htr_btr(:,jn) = vsum * tssum(:,jp_tem) * r1_sjk(:,1,jn) 226 str_btr(:,jn) = vsum * tssum(:,jp_sal) * r1_sjk(:,1,jn) 227 z2d(1,:) = htr_btr(:,jn) * rc_pwatt ! (conversion in PW) 228 DO ji = 1, jpi 229 z2d(ji,:) = z2d(1,:) 230 ENDDO 231 cl1 = TRIM('sophtbtr_'//clsubb(jn)) 232 CALL iom_put( cl1, z2d ) 233 z2d(1,:) = str_btr(:,jn) * rc_ggram ! (conversion in Gg) 234 DO ji = 1, jpi 235 z2d(ji,:) = z2d(1,:) 236 ENDDO 237 cl1 = TRIM('sopstbtr_'//clsubb(jn)) 238 CALL iom_put( cl1, z2d ) 239 ENDDO 240 ENDIF !ln_subbas 241 ENDIF !iom_use("sopstbtr....) 187 242 ! 188 243 ELSE 189 244 ! 190 zmask(:,:,:) = 0._wp 191 zts(:,:,:,:) = 0._wp 192 IF( iom_use( 'zotem' ) .OR. iom_use( 'zosal' ) .OR. iom_use( 'zosrf' ) ) THEN ! i-mean i-k-surface 245 IF( iom_use("zotemglo") ) THEN ! i-mean i-k-surface 193 246 DO jk = 1, jpkm1 194 247 DO jj = 1, jpj … … 201 254 END DO 202 255 END DO 203 !204 256 DO jn = 1, nptr 205 257 zmask(1,:,:) = ptr_sjk( zmask(:,:,:), btmsk(:,:,jn) ) 206 z4d1(:,:,:,jn) = zmask(:,:,:) 207 ENDDO 208 CALL iom_put( 'zosrf', z4d1 ) 209 ! 210 DO jn = 1, nptr 211 z4d2(1,:,:,jn) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) & 212 & / MAX( z4d1(1,:,:,jn), 10.e-15 ) 258 cl1 = TRIM('zosrf'//clsubb(jn) ) 259 CALL iom_put( cl1, zmask ) 260 ! 261 z3d(1,:,:) = ptr_sjk( zts(:,:,:,jp_tem), btmsk(:,:,jn) ) & 262 & / MAX( zmask(1,:,:), 10.e-15 ) 213 263 DO ji = 1, jpi 214 z4d2(ji,:,:,jn) = z4d2(1,:,:,jn) 215 ENDDO 216 ENDDO 217 CALL iom_put( 'zotem', z4d2 ) 218 ! 219 DO jn = 1, nptr 220 z4d2(1,:,:,jn) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) & 221 & / MAX( z4d1(1,:,:,jn), 10.e-15 ) 264 z3d(ji,:,:) = z3d(1,:,:) 265 ENDDO 266 cl1 = TRIM('zotem'//clsubb(jn) ) 267 CALL iom_put( cl1, z3d ) 268 ! 269 z3d(1,:,:) = ptr_sjk( zts(:,:,:,jp_sal), btmsk(:,:,jn) ) & 270 & / MAX( zmask(1,:,:), 10.e-15 ) 222 271 DO ji = 1, jpi 223 z 4d2(ji,:,:,jn) = z4d2(1,:,:,jn)224 ENDDO 225 ENDDO226 CALL iom_put( 'zosal', z4d2)227 !272 z3d(ji,:,:) = z3d(1,:,:) 273 ENDDO 274 cl1 = TRIM('zosal'//clsubb(jn) ) 275 CALL iom_put( cl1, z3d ) 276 END DO 228 277 ENDIF 229 278 ! 230 279 ! ! Advective and diffusive heat and salt transport 231 IF( iom_use( 'sophtadv' ) .OR. iom_use( 'sopstadv' ) ) THEN 232 ! 233 DO jn = 1, nptr 234 z3dtr(1,:,jn) = hstr_adv(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) 280 IF( iom_use("sophtadv") .OR. iom_use("sopstadv") ) THEN 281 z2d(1,:) = htr_adv(:,1) * rc_pwatt ! (conversion in PW) 282 DO ji = 1, jpi 283 z2d(ji,:) = z2d(1,:) 284 ENDDO 285 cl1 = 'sophtadv' 286 CALL iom_put( TRIM(cl1), z2d ) 287 z2d(1,:) = str_adv(:,1) * rc_ggram ! (conversion in Gg) 288 DO ji = 1, jpi 289 z2d(ji,:) = z2d(1,:) 290 ENDDO 291 cl1 = 'sopstadv' 292 CALL iom_put( TRIM(cl1), z2d ) 293 IF( ln_subbas ) THEN 294 DO jn=2,nptr 295 z2d(1,:) = htr_adv(:,jn) * rc_pwatt ! (conversion in PW) 235 296 DO ji = 1, jpi 236 z3dtr(ji,:,jn) = z3dtr(1,:,jn) 237 ENDDO 238 ENDDO 239 CALL iom_put( 'sophtadv', z3dtr ) 240 DO jn = 1, nptr 241 z3dtr(1,:,jn) = hstr_adv(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) 297 z2d(ji,:) = z2d(1,:) 298 ENDDO 299 cl1 = TRIM('sophtadv_'//clsubb(jn)) 300 CALL iom_put( cl1, z2d ) 301 z2d(1,:) = str_adv(:,jn) * rc_ggram ! (conversion in Gg) 242 302 DO ji = 1, jpi 243 z3dtr(ji,:,jn) = z3dtr(1,:,jn) 244 ENDDO 245 ENDDO 246 CALL iom_put( 'sopstadv', z3dtr ) 247 ENDIF 248 ! 249 IF( iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) ) THEN 250 ! 251 DO jn = 1, nptr 252 z3dtr(1,:,jn) = hstr_ldf(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) 303 z2d(ji,:) = z2d(1,:) 304 ENDDO 305 cl1 = TRIM('sopstadv_'//clsubb(jn)) 306 CALL iom_put( cl1, z2d ) 307 ENDDO 308 ENDIF 309 ENDIF 310 ! 311 IF( iom_use("sophtldf") .OR. iom_use("sopstldf") ) THEN 312 z2d(1,:) = htr_ldf(:,1) * rc_pwatt ! (conversion in PW) 313 DO ji = 1, jpi 314 z2d(ji,:) = z2d(1,:) 315 ENDDO 316 cl1 = 'sophtldf' 317 CALL iom_put( TRIM(cl1), z2d ) 318 z2d(1,:) = str_ldf(:,1) * rc_ggram ! (conversion in Gg) 319 DO ji = 1, jpi 320 z2d(ji,:) = z2d(1,:) 321 ENDDO 322 cl1 = 'sopstldf' 323 CALL iom_put( TRIM(cl1), z2d ) 324 IF( ln_subbas ) THEN 325 DO jn=2,nptr 326 z2d(1,:) = htr_ldf(:,jn) * rc_pwatt ! (conversion in PW) 253 327 DO ji = 1, jpi 254 z3dtr(ji,:,jn) = z3dtr(1,:,jn) 255 ENDDO 256 ENDDO 257 CALL iom_put( 'sophtldf', z3dtr ) 258 DO jn = 1, nptr 259 z3dtr(1,:,jn) = hstr_ldf(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) 328 z2d(ji,:) = z2d(1,:) 329 ENDDO 330 cl1 = TRIM('sophtldf_'//clsubb(jn)) 331 CALL iom_put( cl1, z2d ) 332 z2d(1,:) = str_ldf(:,jn) * rc_ggram ! (conversion in Gg) 260 333 DO ji = 1, jpi 261 z3dtr(ji,:,jn) = z3dtr(1,:,jn) 262 ENDDO 263 ENDDO 264 CALL iom_put( 'sopstldf', z3dtr ) 265 ENDIF 266 ! 267 IF( iom_use( 'sophteiv' ) .OR. iom_use( 'sopsteiv' ) ) THEN 268 ! 269 DO jn = 1, nptr 270 z3dtr(1,:,jn) = hstr_eiv(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) 271 DO ji = 1, jpi 272 z3dtr(ji,:,jn) = z3dtr(1,:,jn) 273 ENDDO 274 ENDDO 275 CALL iom_put( 'sophteiv', z3dtr ) 276 DO jn = 1, nptr 277 z3dtr(1,:,jn) = hstr_eiv(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) 278 DO ji = 1, jpi 279 z3dtr(ji,:,jn) = z3dtr(1,:,jn) 280 ENDDO 281 ENDDO 282 CALL iom_put( 'sopsteiv', z3dtr ) 283 ENDIF 284 ! 285 IF( iom_use( 'sopstvtr' ) .OR. iom_use( 'sophtvtr' ) ) THEN 286 zts(:,:,:,:) = 0._wp 287 DO jk = 1, jpkm1 288 DO jj = 1, jpjm1 334 z2d(ji,:) = z2d(1,:) 335 ENDDO 336 cl1 = TRIM('sopstldf_'//clsubb(jn)) 337 CALL iom_put( cl1, z2d ) 338 ENDDO 339 ENDIF 340 ENDIF 341 342 IF( iom_use("sophteiv") .OR. iom_use("sopsteiv") ) THEN 343 z2d(1,:) = htr_eiv(:,1) * rc_pwatt ! (conversion in PW) 344 DO ji = 1, jpi 345 z2d(ji,:) = z2d(1,:) 346 ENDDO 347 cl1 = 'sophteiv' 348 CALL iom_put( TRIM(cl1), z2d ) 349 z2d(1,:) = str_eiv(:,1) * rc_ggram ! (conversion in Gg) 350 DO ji = 1, jpi 351 z2d(ji,:) = z2d(1,:) 352 ENDDO 353 cl1 = 'sopsteiv' 354 CALL iom_put( TRIM(cl1), z2d ) 355 IF( ln_subbas ) THEN 356 DO jn=2,nptr 357 z2d(1,:) = htr_eiv(:,jn) * rc_pwatt ! (conversion in PW) 289 358 DO ji = 1, jpi 290 zvfc = e1v(ji,jj) * e3v_n(ji,jj,jk) 291 zts(ji,jj,jk,jp_tem) = (tsn(ji,jj,jk,jp_tem)+tsn(ji,jj+1,jk,jp_tem)) * 0.5 * zvfc !Tracers averaged onto V grid 292 zts(ji,jj,jk,jp_sal) = (tsn(ji,jj,jk,jp_sal)+tsn(ji,jj+1,jk,jp_sal)) * 0.5 * zvfc 359 z2d(ji,:) = z2d(1,:) 293 360 ENDDO 294 ENDDO 295 ENDDO 296 CALL dia_ptr_hst( jp_tem, 'vtr', zts(:,:,:,jp_tem) ) 297 CALL dia_ptr_hst( jp_sal, 'vtr', zts(:,:,:,jp_sal) ) 298 DO jn = 1, nptr 299 z3dtr(1,:,jn) = hstr_vtr(:,jp_tem,jn) * rc_pwatt ! (conversion in PW) 300 DO ji = 1, jpi 301 z3dtr(ji,:,jn) = z3dtr(1,:,jn) 302 ENDDO 303 ENDDO 304 CALL iom_put( 'sophtvtr', z3dtr ) 305 DO jn = 1, nptr 306 z3dtr(1,:,jn) = hstr_vtr(:,jp_sal,jn) * rc_ggram ! (conversion in Gg) 307 DO ji = 1, jpi 308 z3dtr(ji,:,jn) = z3dtr(1,:,jn) 309 ENDDO 310 ENDDO 311 CALL iom_put( 'sopstvtr', z3dtr ) 312 ENDIF 313 ! 314 IF( iom_use( 'uocetr_vsum_cumul' ) ) THEN 315 CALL iom_get_var( 'uocetr_vsum_op', z2d ) ! get uocetr_vsum_op from xml 316 z2d(:,:) = ptr_ci_2d( z2d(:,:) ) 317 CALL iom_put( 'uocetr_vsum_cumul', z2d ) 361 cl1 = TRIM('sophteiv_'//clsubb(jn)) 362 CALL iom_put( cl1, z2d ) 363 z2d(1,:) = str_eiv(:,jn) * rc_ggram ! (conversion in Gg) 364 DO ji = 1, jpi 365 z2d(ji,:) = z2d(1,:) 366 ENDDO 367 cl1 = TRIM('sopsteiv_'//clsubb(jn)) 368 CALL iom_put( cl1, z2d ) 369 ENDDO 370 ENDIF 318 371 ENDIF 319 372 ! … … 331 384 !! ** Purpose : Initialization, namelist read 332 385 !!---------------------------------------------------------------------- 333 INTEGER :: inum, jn ! local integers 334 !! 335 REAL(wp), DIMENSION(jpi,jpj) :: zmsk 336 !!---------------------------------------------------------------------- 337 338 l_diaptr = .FALSE. 339 IF( iom_use( 'zomsf' ) .OR. iom_use( 'zotem' ) .OR. iom_use( 'zosal' ) .OR. & 340 & iom_use( 'zosrf' ) .OR. iom_use( 'sopstove' ) .OR. iom_use( 'sophtove' ) .OR. & 341 & iom_use( 'sopstbtr' ) .OR. iom_use( 'sophtbtr' ) .OR. iom_use( 'sophtadv' ) .OR. & 342 & iom_use( 'sopstadv' ) .OR. iom_use( 'sophtldf' ) .OR. iom_use( 'sopstldf' ) .OR. & 343 & iom_use( 'sophteiv' ) .OR. iom_use( 'sopsteiv' ) .OR. iom_use( 'sopstvtr' ) .OR. & 344 & iom_use( 'sophtvtr' ) .OR. iom_use( 'uocetr_vsum_cumul' ) ) l_diaptr = .TRUE. 345 346 386 INTEGER :: jn ! local integers 387 INTEGER :: inum, ierr ! local integers 388 INTEGER :: ios ! Local integer output status for namelist read 389 !! 390 NAMELIST/namptr/ ln_diaptr, ln_subbas 391 !!---------------------------------------------------------------------- 392 393 REWIND( numnam_ref ) ! Namelist namptr in reference namelist : Poleward transport 394 READ ( numnam_ref, namptr, IOSTAT = ios, ERR = 901) 395 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namptr in reference namelist' ) 396 397 REWIND( numnam_cfg ) ! Namelist namptr in configuration namelist : Poleward transport 398 READ ( numnam_cfg, namptr, IOSTAT = ios, ERR = 902 ) 399 902 IF( ios > 0 ) CALL ctl_nam ( ios , 'namptr in configuration namelist' ) 400 IF(lwm) WRITE ( numond, namptr ) 401 347 402 IF(lwp) THEN ! Control print 348 403 WRITE(numout,*) … … 350 405 WRITE(numout,*) '~~~~~~~~~~~~' 351 406 WRITE(numout,*) ' Namelist namptr : set ptr parameters' 352 WRITE(numout,*) ' Poleward heat & salt transport (T) or not (F) l_diaptr = ', l_diaptr 407 WRITE(numout,*) ' Poleward heat & salt transport (T) or not (F) ln_diaptr = ', ln_diaptr 408 WRITE(numout,*) ' Global (F) or glo/Atl/Pac/Ind/Indo-Pac basins ln_subbas = ', ln_subbas 353 409 ENDIF 354 410 355 IF( l_diaptr ) THEN 356 ! 411 IF( ln_diaptr ) THEN 412 ! 413 IF( ln_subbas ) THEN 414 nptr = 5 ! Global, Atlantic, Pacific, Indian, Indo-Pacific 415 ALLOCATE( clsubb(nptr) ) 416 clsubb(1) = 'glo' ; clsubb(2) = 'atl' ; clsubb(3) = 'pac' ; clsubb(4) = 'ind' ; clsubb(5) = 'ipc' 417 ELSE 418 nptr = 1 ! Global only 419 ALLOCATE( clsubb(nptr) ) 420 clsubb(1) = 'glo' 421 ENDIF 422 423 ! ! allocate dia_ptr arrays 357 424 IF( dia_ptr_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'dia_ptr_init : unable to allocate arrays' ) 358 425 359 426 rc_pwatt = rc_pwatt * rau0_rcp ! conversion from K.s-1 to PetaWatt 360 rc_ggram = rc_ggram * rau0 ! conversion from m3/s to Gg/s361 427 362 428 IF( lk_mpp ) CALL mpp_ini_znl( numout ) ! Define MPI communicator for zonal sum 363 429 364 btmsk(:,:,1) = tmask_i(:,:) 365 CALL iom_open( 'subbasins', inum, ldstop = .FALSE. ) 366 CALL iom_get( inum, jpdom_data, 'atlmsk', btmsk(:,:,2) ) ! Atlantic basin 367 CALL iom_get( inum, jpdom_data, 'pacmsk', btmsk(:,:,3) ) ! Pacific basin 368 CALL iom_get( inum, jpdom_data, 'indmsk', btmsk(:,:,4) ) ! Indian basin 369 CALL iom_close( inum ) 370 btmsk(:,:,5) = MAX ( btmsk(:,:,3), btmsk(:,:,4) ) ! Indo-Pacific basin 371 DO jn = 2, nptr 430 IF( ln_subbas ) THEN ! load sub-basin mask 431 CALL iom_open( 'subbasins', inum, ldstop = .FALSE. ) 432 CALL iom_get( inum, jpdom_data, 'atlmsk', btmsk(:,:,2) ) ! Atlantic basin 433 CALL iom_get( inum, jpdom_data, 'pacmsk', btmsk(:,:,3) ) ! Pacific basin 434 CALL iom_get( inum, jpdom_data, 'indmsk', btmsk(:,:,4) ) ! Indian basin 435 CALL iom_close( inum ) 436 btmsk(:,:,5) = MAX ( btmsk(:,:,3), btmsk(:,:,4) ) ! Indo-Pacific basin 437 WHERE( gphit(:,:) < -30._wp) ; btm30(:,:) = 0._wp ! mask out Southern Ocean 438 ELSE WHERE ; btm30(:,:) = ssmask(:,:) 439 END WHERE 440 ENDIF 441 442 btmsk(:,:,1) = tmask_i(:,:) ! global ocean 443 444 DO jn = 1, nptr 372 445 btmsk(:,:,jn) = btmsk(:,:,jn) * tmask_i(:,:) ! interior domain only 373 446 END DO 374 ! JD : modification so that overturning streamfunction is available in Atlantic at 34S to compare with observations375 WHERE( gphit(:,:)*tmask_i(:,:) < -34._wp)376 zmsk(:,:) = 0._wp ! mask out Southern Ocean377 ELSE WHERE378 zmsk(:,:) = ssmask(:,:)379 END WHERE380 btmsk34(:,:,1) = btmsk(:,:,1)381 DO jn = 2, nptr382 btmsk34(:,:,jn) = btmsk(:,:,jn) * zmsk(:,:) ! interior domain only383 ENDDO384 447 385 448 ! Initialise arrays to zero because diatpr is called before they are first calculated 386 449 ! Note that this means diagnostics will not be exactly correct when model run is restarted. 387 hstr_adv(:,:,:) = 0._wp 388 hstr_ldf(:,:,:) = 0._wp 389 hstr_eiv(:,:,:) = 0._wp 390 hstr_ove(:,:,:) = 0._wp 391 hstr_btr(:,:,:) = 0._wp ! 392 hstr_vtr(:,:,:) = 0._wp ! 393 ! 394 ll_init = .FALSE. 450 htr_adv(:,:) = 0._wp ; str_adv(:,:) = 0._wp 451 htr_ldf(:,:) = 0._wp ; str_ldf(:,:) = 0._wp 452 htr_eiv(:,:) = 0._wp ; str_eiv(:,:) = 0._wp 453 htr_ove(:,:) = 0._wp ; str_ove(:,:) = 0._wp 454 htr_btr(:,:) = 0._wp ; str_btr(:,:) = 0._wp 395 455 ! 396 456 ENDIF … … 411 471 INTEGER :: jn ! 412 472 413 !414 473 IF( cptr == 'adv' ) THEN 415 IF( ktra == jp_tem ) THEN 416 DO jn = 1, nptr 417 hstr_adv(:,jp_tem,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 418 ENDDO 419 ENDIF 420 IF( ktra == jp_sal ) THEN 421 DO jn = 1, nptr 422 hstr_adv(:,jp_sal,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 423 ENDDO 424 ENDIF 474 IF( ktra == jp_tem ) htr_adv(:,1) = ptr_sj( pva(:,:,:) ) 475 IF( ktra == jp_sal ) str_adv(:,1) = ptr_sj( pva(:,:,:) ) 425 476 ENDIF 426 !427 477 IF( cptr == 'ldf' ) THEN 428 IF( ktra == jp_tem ) THEN 429 DO jn = 1, nptr 430 hstr_ldf(:,jp_tem,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 431 ENDDO 432 ENDIF 433 IF( ktra == jp_sal ) THEN 434 DO jn = 1, nptr 435 hstr_ldf(:,jp_sal,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 436 ENDDO 437 ENDIF 478 IF( ktra == jp_tem ) htr_ldf(:,1) = ptr_sj( pva(:,:,:) ) 479 IF( ktra == jp_sal ) str_ldf(:,1) = ptr_sj( pva(:,:,:) ) 438 480 ENDIF 439 !440 481 IF( cptr == 'eiv' ) THEN 441 IF( ktra == jp_tem ) THEN 442 DO jn = 1, nptr 443 hstr_eiv(:,jp_tem,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 444 ENDDO 445 ENDIF 446 IF( ktra == jp_sal ) THEN 447 DO jn = 1, nptr 448 hstr_eiv(:,jp_sal,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 449 ENDDO 450 ENDIF 482 IF( ktra == jp_tem ) htr_eiv(:,1) = ptr_sj( pva(:,:,:) ) 483 IF( ktra == jp_sal ) str_eiv(:,1) = ptr_sj( pva(:,:,:) ) 451 484 ENDIF 452 485 ! 453 IF( cptr == 'vtr' ) THEN 454 IF( ktra == jp_tem ) THEN 455 DO jn = 1, nptr 456 hstr_vtr(:,jp_tem,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 457 ENDDO 458 ENDIF 459 IF( ktra == jp_sal ) THEN 460 DO jn = 1, nptr 461 hstr_vtr(:,jp_sal,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 462 ENDDO 463 ENDIF 486 IF( ln_subbas ) THEN 487 ! 488 IF( cptr == 'adv' ) THEN 489 IF( ktra == jp_tem ) THEN 490 DO jn = 2, nptr 491 htr_adv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 492 END DO 493 ENDIF 494 IF( ktra == jp_sal ) THEN 495 DO jn = 2, nptr 496 str_adv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 497 END DO 498 ENDIF 499 ENDIF 500 IF( cptr == 'ldf' ) THEN 501 IF( ktra == jp_tem ) THEN 502 DO jn = 2, nptr 503 htr_ldf(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 504 END DO 505 ENDIF 506 IF( ktra == jp_sal ) THEN 507 DO jn = 2, nptr 508 str_ldf(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 509 END DO 510 ENDIF 511 ENDIF 512 IF( cptr == 'eiv' ) THEN 513 IF( ktra == jp_tem ) THEN 514 DO jn = 2, nptr 515 htr_eiv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 516 END DO 517 ENDIF 518 IF( ktra == jp_sal ) THEN 519 DO jn = 2, nptr 520 str_eiv(:,jn) = ptr_sj( pva(:,:,:), btmsk(:,:,jn) ) 521 END DO 522 ENDIF 523 ENDIF 524 ! 464 525 ENDIF 465 !466 526 END SUBROUTINE dia_ptr_hst 467 527 … … 476 536 ierr(:) = 0 477 537 ! 478 IF( .NOT. ALLOCATED( btmsk ) ) THEN 479 ALLOCATE( btmsk(jpi,jpj,nptr) , btmsk34(jpi,jpj,nptr), & 480 & hstr_adv(jpj,jpts,nptr), hstr_eiv(jpj,jpts,nptr), & 481 & hstr_ove(jpj,jpts,nptr), hstr_btr(jpj,jpts,nptr), & 482 & hstr_ldf(jpj,jpts,nptr), hstr_vtr(jpj,jpts,nptr), STAT=ierr(1) ) 483 ! 484 ALLOCATE( p_fval1d(jpj), p_fval2d(jpj,jpk), Stat=ierr(2)) 485 ! 486 dia_ptr_alloc = MAXVAL( ierr ) 487 CALL mpp_sum( 'diaptr', dia_ptr_alloc ) 488 ENDIF 538 ALLOCATE( btmsk(jpi,jpj,nptr) , & 539 & htr_adv(jpj,nptr) , str_adv(jpj,nptr) , & 540 & htr_eiv(jpj,nptr) , str_eiv(jpj,nptr) , & 541 & htr_ove(jpj,nptr) , str_ove(jpj,nptr) , & 542 & htr_btr(jpj,nptr) , str_btr(jpj,nptr) , & 543 & htr_ldf(jpj,nptr) , str_ldf(jpj,nptr) , STAT=ierr(1) ) 544 ! 545 ALLOCATE( p_fval1d(jpj), p_fval2d(jpj,jpk), Stat=ierr(2)) 546 ! 547 ALLOCATE( btm30(jpi,jpj), STAT=ierr(3) ) 548 549 ! 550 dia_ptr_alloc = MAXVAL( ierr ) 551 CALL mpp_sum( 'diaptr', dia_ptr_alloc ) 489 552 ! 490 553 END FUNCTION dia_ptr_alloc … … 502 565 !! ** Action : - p_fval: i-k-mean poleward flux of pva 503 566 !!---------------------------------------------------------------------- 504 REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) :: pva ! mask flux array at V-point505 REAL(wp), INTENT(in), DIMENSION(jpi,jpj) 567 REAL(wp), INTENT(in), DIMENSION(jpi,jpj,jpk) :: pva ! mask flux array at V-point 568 REAL(wp), INTENT(in), DIMENSION(jpi,jpj), OPTIONAL :: pmsk ! Optional 2D basin mask 506 569 ! 507 570 INTEGER :: ji, jj, jk ! dummy loop arguments … … 514 577 ijpj = jpj 515 578 p_fval(:) = 0._wp 516 DO jk = 1, jpkm1 517 DO jj = 2, jpjm1 518 DO ji = fs_2, fs_jpim1 ! Vector opt. 519 p_fval(jj) = p_fval(jj) + pva(ji,jj,jk) * pmsk(ji,jj) * tmask_i(ji,jj) 579 IF( PRESENT( pmsk ) ) THEN 580 DO jk = 1, jpkm1 581 DO jj = 2, jpjm1 582 DO ji = fs_2, fs_jpim1 ! Vector opt. 583 p_fval(jj) = p_fval(jj) + pva(ji,jj,jk) * tmask_i(ji,jj) * pmsk(ji,jj) 584 END DO 520 585 END DO 521 586 END DO 522 END DO 587 ELSE 588 DO jk = 1, jpkm1 589 DO jj = 2, jpjm1 590 DO ji = fs_2, fs_jpim1 ! Vector opt. 591 p_fval(jj) = p_fval(jj) + pva(ji,jj,jk) * tmask_i(ji,jj) 592 END DO 593 END DO 594 END DO 595 ENDIF 523 596 #if defined key_mpp_mpi 524 597 CALL mpp_sum( 'diaptr', p_fval, ijpj, ncomm_znl) … … 539 612 !! ** Action : - p_fval: i-k-mean poleward flux of pva 540 613 !!---------------------------------------------------------------------- 541 REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) :: pva ! mask flux array at V-point542 REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) :: pmsk ! Optional 2D basin mask614 REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) :: pva ! mask flux array at V-point 615 REAL(wp) , INTENT(in), DIMENSION(jpi,jpj), OPTIONAL :: pmsk ! Optional 2D basin mask 543 616 ! 544 617 INTEGER :: ji,jj ! dummy loop arguments … … 551 624 ijpj = jpj 552 625 p_fval(:) = 0._wp 553 DO jj = 2, jpjm1 554 DO ji = fs_2, fs_jpim1 ! Vector opt. 555 p_fval(jj) = p_fval(jj) + pva(ji,jj) * pmsk(ji,jj) * tmask_i(ji,jj) 626 IF( PRESENT( pmsk ) ) THEN 627 DO jj = 2, jpjm1 628 DO ji = nldi, nlei ! No vector optimisation here. Better use a mask ? 629 p_fval(jj) = p_fval(jj) + pva(ji,jj) * tmask_i(ji,jj) * pmsk(ji,jj) 630 END DO 556 631 END DO 557 END DO 632 ELSE 633 DO jj = 2, jpjm1 634 DO ji = nldi, nlei ! No vector optimisation here. Better use a mask ? 635 p_fval(jj) = p_fval(jj) + pva(ji,jj) * tmask_i(ji,jj) 636 END DO 637 END DO 638 ENDIF 558 639 #if defined key_mpp_mpi 559 640 CALL mpp_sum( 'diaptr', p_fval, ijpj, ncomm_znl ) … … 562 643 END FUNCTION ptr_sj_2d 563 644 564 FUNCTION ptr_ci_2d( pva ) RESULT ( p_fval )565 !!----------------------------------------------------------------------566 !! *** ROUTINE ptr_ci_2d ***567 !!568 !! ** Purpose : "meridional" cumulated sum computation of a j-flux array569 !!570 !! ** Method : - j cumulated sum of pva using the interior 2D vmask (umask_i).571 !!572 !! ** Action : - p_fval: j-cumulated sum of pva573 !!----------------------------------------------------------------------574 REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) :: pva ! mask flux array at V-point575 !576 INTEGER :: ji,jj,jc ! dummy loop arguments577 INTEGER :: ijpj ! ???578 REAL(wp), DIMENSION(jpi,jpj) :: p_fval ! function value579 !!--------------------------------------------------------------------580 !581 ijpj = jpj ! ???582 p_fval(:,:) = 0._wp583 DO jc = 1, jpnj ! looping over all processors in j axis584 DO jj = 2, jpjm1585 DO ji = fs_2, fs_jpim1 ! Vector opt.586 p_fval(ji,jj) = p_fval(ji,jj-1) + pva(ji,jj) * tmask_i(ji,jj)587 END DO588 END DO589 CALL lbc_lnk( 'diaptr', p_fval, 'U', -1. )590 END DO591 !592 END FUNCTION ptr_ci_2d593 594 595 645 596 646 FUNCTION ptr_sjk( pta, pmsk ) RESULT ( p_fval ) … … 606 656 !! 607 657 IMPLICIT none 608 REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) :: pta ! mask flux array at V-point609 REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) :: pmsk ! Optional 2D basin mask658 REAL(wp) , INTENT(in), DIMENSION(jpi,jpj,jpk) :: pta ! mask flux array at V-point 659 REAL(wp) , INTENT(in), DIMENSION(jpi,jpj) , OPTIONAL :: pmsk ! Optional 2D basin mask 610 660 !! 611 661 INTEGER :: ji, jj, jk ! dummy loop arguments … … 623 673 p_fval(:,:) = 0._wp 624 674 ! 625 DO jk = 1, jpkm1 626 DO jj = 2, jpjm1 627 DO ji = fs_2, fs_jpim1 ! Vector opt. 628 p_fval(jj,jk) = p_fval(jj,jk) + pta(ji,jj,jk) * pmsk(ji,jj) * tmask_i(ji,jj) 675 IF( PRESENT( pmsk ) ) THEN 676 DO jk = 1, jpkm1 677 DO jj = 2, jpjm1 678 !!gm here, use of tmask_i ==> no need of loop over nldi, nlei.... 679 DO ji = nldi, nlei ! No vector optimisation here. Better use a mask ? 680 p_fval(jj,jk) = p_fval(jj,jk) + pta(ji,jj,jk) * pmsk(ji,jj) 681 END DO 629 682 END DO 630 683 END DO 631 END DO 684 ELSE 685 DO jk = 1, jpkm1 686 DO jj = 2, jpjm1 687 DO ji = nldi, nlei ! No vector optimisation here. Better use a mask ? 688 p_fval(jj,jk) = p_fval(jj,jk) + pta(ji,jj,jk) * tmask_i(ji,jj) 689 END DO 690 END DO 691 END DO 692 END IF 632 693 ! 633 694 #if defined key_mpp_mpi
Note: See TracChangeset
for help on using the changeset viewer.
