- Timestamp:
- 2019-12-11T17:15:54+01:00 (4 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r11943_MERGE_2019/src/OCE/DIA/diaar5.F90
r11949 r12193 72 72 INTEGER, INTENT( in ) :: Kmm ! ocean time level index 73 73 ! 74 INTEGER :: ji, jj, jk ! dummy loop arguments75 REAL(wp) :: zvolssh, zvol, zssh_steric, zztmp, zarho, ztemp, zsal, zmass 74 INTEGER :: ji, jj, jk, iks, ikb ! dummy loop arguments 75 REAL(wp) :: zvolssh, zvol, zssh_steric, zztmp, zarho, ztemp, zsal, zmass, zsst 76 76 REAL(wp) :: zaw, zbw, zrw 77 77 ! 78 78 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zarea_ssh , zbotpres ! 2D workspace 79 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zpe 80 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zrhd , zrhop 79 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zpe, z2d ! 2D workspace 80 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zrhd , zrhop, ztpot ! 3D workspace 81 81 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: ztsn ! 4D workspace 82 82 … … 87 87 88 88 IF( l_ar5 ) THEN 89 ALLOCATE( zarea_ssh(jpi,jpj) , zbotpres(jpi,jpj) )89 ALLOCATE( zarea_ssh(jpi,jpj), zbotpres(jpi,jpj), z2d(jpi,jpj) ) 90 90 ALLOCATE( zrhd(jpi,jpj,jpk) , zrhop(jpi,jpj,jpk) ) 91 91 ALLOCATE( ztsn(jpi,jpj,jpk,jpts) ) … … 93 93 ENDIF 94 94 ! 95 CALL iom_put( 'e2u' , e2u (:,:) ) 96 CALL iom_put( 'e1v' , e1v (:,:) ) 97 CALL iom_put( 'areacello', area(:,:) ) 98 ! 99 IF( iom_use( 'volcello' ) .OR. iom_use( 'masscello' ) ) THEN 100 zrhd(:,:,jpk) = 0._wp ! ocean volume ; rhd is used as workspace 101 DO jk = 1, jpkm1 102 zrhd(:,:,jk) = area(:,:) * e3t(:,:,jk,Kmm) * tmask(:,:,jk) 103 END DO 104 CALL iom_put( 'volcello' , zrhd(:,:,:) ) ! WARNING not consistent with CMIP DR where volcello is at ca. 2000 105 CALL iom_put( 'masscello' , rau0 * e3t(:,:,:,Kmm) * tmask(:,:,:) ) ! ocean mass 106 ENDIF 107 ! 108 IF( iom_use( 'e3tb' ) ) THEN ! bottom layer thickness 109 DO jj = 1, jpj 110 DO ji = 1, jpi 111 ikb = mbkt(ji,jj) 112 z2d(ji,jj) = e3t(ji,jj,ikb,Kmm) 113 END DO 114 END DO 115 CALL iom_put( 'e3tb', z2d ) 116 ENDIF 117 ! 95 118 IF( iom_use( 'voltot' ) .OR. iom_use( 'sshtot' ) .OR. iom_use( 'sshdyn' ) ) THEN 96 119 ! ! total volume of liquid seawater 97 zvolssh = SUM( zarea_ssh(:,:) ) 98 CALL mpp_sum( 'diaar5', zvolssh ) 99 zvol = vol0 + zvolssh 120 zvolssh = glob_sum( 'diaar5', zarea_ssh(:,:) ) 121 zvol = vol0 + zvolssh 100 122 101 123 CALL iom_put( 'voltot', zvol ) … … 119 141 DO ji = 1, jpi 120 142 DO jj = 1, jpj 121 zbotpres(ji,jj) = zbotpres(ji,jj) + ssh(ji,jj,Kmm) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj) 143 iks = mikt(ji,jj) 144 zbotpres(ji,jj) = zbotpres(ji,jj) + ssh(ji,jj,Kmm) * zrhd(ji,jj,iks) + riceload(ji,jj) 122 145 END DO 123 146 END DO … … 130 153 END IF 131 154 ! 132 zarho = SUM( area(:,:) * zbotpres(:,:) ) 133 CALL mpp_sum( 'diaar5', zarho ) 155 zarho = glob_sum( 'diaar5', area(:,:) * zbotpres(:,:) ) 134 156 zssh_steric = - zarho / area_tot 135 157 CALL iom_put( 'sshthster', zssh_steric ) … … 148 170 DO ji = 1,jpi 149 171 DO jj = 1,jpj 150 zbotpres(ji,jj) = zbotpres(ji,jj) + ssh(ji,jj,Kmm) * zrhd(ji,jj,mikt(ji,jj)) + riceload(ji,jj) 172 iks = mikt(ji,jj) 173 zbotpres(ji,jj) = zbotpres(ji,jj) + ssh(ji,jj,Kmm) * zrhd(ji,jj,iks) + riceload(ji,jj) 151 174 END DO 152 175 END DO … … 156 179 END IF 157 180 ! 158 zarho = SUM( area(:,:) * zbotpres(:,:) ) 159 CALL mpp_sum( 'diaar5', zarho ) 181 zarho = glob_sum( 'diaar5', area(:,:) * zbotpres(:,:) ) 160 182 zssh_steric = - zarho / area_tot 161 183 CALL iom_put( 'sshsteric', zssh_steric ) 162 163 184 ! ! ocean bottom pressure 164 185 zztmp = rau0 * grav * 1.e-4_wp ! recover pressure from pressure anomaly and cover to dbar = 1.e4 Pa … … 169 190 170 191 IF( iom_use( 'masstot' ) .OR. iom_use( 'temptot' ) .OR. iom_use( 'saltot' ) ) THEN 171 ! ! Mean density anomalie, temperature and salinity172 ztemp = 0._wp173 zsal = 0._wp174 DO jk = 1, jpkm1175 DO jj = 1, jpj176 DO ji = 1, jpi177 zztmp = area(ji,jj) * e3t(ji,jj,jk,Kmm)178 ztemp = ztemp + zztmp * ts(ji,jj,jk,jp_tem,Kmm)179 zsal = zsal + zztmp * ts(ji,jj,jk,jp_sal,Kmm)180 ENDDO181 ENDDO182 END DO 183 IF( ln_linssh ) THEN192 ! ! Mean density anomalie, temperature and salinity 193 ztsn(:,:,:,:) = 0._wp ! ztsn(:,:,1,jp_tem/sal) is used here as 2D Workspace for temperature & salinity 194 DO jk = 1, jpkm1 195 DO jj = 1, jpj 196 DO ji = 1, jpi 197 zztmp = area(ji,jj) * e3t(ji,jj,jk,Kmm) 198 ztsn(ji,jj,1,jp_tem) = ztsn(ji,jj,1,jp_tem) + zztmp * ts(ji,jj,jk,jp_tem,Kmm) 199 ztsn(ji,jj,1,jp_sal) = ztsn(ji,jj,1,jp_sal) + zztmp * ts(ji,jj,jk,jp_sal,Kmm) 200 ENDDO 201 ENDDO 202 ENDDO 203 204 IF( ln_linssh ) THEN 184 205 IF( ln_isfcav ) THEN 185 206 DO ji = 1, jpi 186 207 DO jj = 1, jpj 187 ztemp = ztemp + zarea_ssh(ji,jj) * ts(ji,jj,mikt(ji,jj),jp_tem,Kmm) 188 zsal = zsal + zarea_ssh(ji,jj) * ts(ji,jj,mikt(ji,jj),jp_sal,Kmm) 208 iks = mikt(ji,jj) 209 ztsn(ji,jj,1,jp_tem) = ztsn(ji,jj,1,jp_tem) + zarea_ssh(ji,jj) * ts(ji,jj,iks,jp_tem,Kmm) 210 ztsn(ji,jj,1,jp_sal) = ztsn(ji,jj,1,jp_sal) + zarea_ssh(ji,jj) * ts(ji,jj,iks,jp_sal,Kmm) 189 211 END DO 190 212 END DO 191 213 ELSE 192 zt emp = ztemp + SUM( zarea_ssh(:,:) * ts(:,:,1,jp_tem,Kmm) )193 z sal = zsal + SUM( zarea_ssh(:,:) * ts(:,:,1,jp_sal,Kmm) )214 ztsn(:,:,1,jp_tem) = ztsn(:,:,1,jp_tem) + zarea_ssh(:,:) * ts(:,:,1,jp_tem,Kmm) 215 ztsn(:,:,1,jp_sal) = ztsn(:,:,1,jp_sal) + zarea_ssh(:,:) * ts(:,:,1,jp_sal,Kmm) 194 216 END IF 195 217 ENDIF 196 IF( lk_mpp ) THEN 197 CALL mpp_sum( 'diaar5', ztemp ) 198 CALL mpp_sum( 'diaar5', zsal ) 199 END IF 200 ! 201 zmass = rau0 * ( zarho + zvol ) ! total mass of liquid seawater 202 ztemp = ztemp / zvol ! potential temperature in liquid seawater 203 zsal = zsal / zvol ! Salinity of liquid seawater 218 ! 219 ztemp = glob_sum( 'diaar5', ztsn(:,:,1,jp_tem) ) 220 zsal = glob_sum( 'diaar5', ztsn(:,:,1,jp_sal) ) 221 zmass = rau0 * ( zarho + zvol ) 204 222 ! 205 223 CALL iom_put( 'masstot', zmass ) 206 CALL iom_put( 'temptot', ztemp ) 207 CALL iom_put( 'saltot' , zsal ) 208 ! 224 CALL iom_put( 'temptot', ztemp / zvol ) 225 CALL iom_put( 'saltot' , zsal / zvol ) 226 ! 227 ENDIF 228 229 IF( ln_teos10 ) THEN ! ! potential temperature (TEOS-10 case) 230 IF( iom_use( 'toce_pot') .OR. iom_use( 'temptot_pot' ) .OR. iom_use( 'sst_pot' ) & 231 .OR. iom_use( 'ssttot' ) .OR. iom_use( 'tosmint_pot' ) ) THEN 232 ! 233 ALLOCATE( ztpot(jpi,jpj,jpk) ) 234 ztpot(:,:,jpk) = 0._wp 235 DO jk = 1, jpkm1 236 ztpot(:,:,jk) = eos_pt_from_ct( ts(:,:,jk,jp_tem,Kmm), ts(:,:,jk,jp_sal,Kmm) ) 237 END DO 238 ! 239 CALL iom_put( 'toce_pot', ztpot(:,:,:) ) ! potential temperature (TEOS-10 case) 240 CALL iom_put( 'sst_pot' , ztpot(:,:,1) ) ! surface temperature 241 ! 242 IF( iom_use( 'temptot_pot' ) ) THEN ! Output potential temperature in case we use TEOS-10 243 z2d(:,:) = 0._wp 244 DO jk = 1, jpkm1 245 z2d(:,:) = z2d(:,:) + area(:,:) * e3t(:,:,jk,Kmm) * ztpot(:,:,jk) 246 END DO 247 ztemp = glob_sum( 'diaar5', z2d(:,:) ) 248 CALL iom_put( 'temptot_pot', ztemp / zvol ) 249 ENDIF 250 ! 251 IF( iom_use( 'ssttot' ) ) THEN ! Output potential temperature in case we use TEOS-10 252 zsst = glob_sum( 'diaar5', area(:,:) * ztpot(:,:,1) ) 253 CALL iom_put( 'ssttot', zsst / area_tot ) 254 ENDIF 255 ! Vertical integral of temperature 256 IF( iom_use( 'tosmint_pot') ) THEN 257 z2d(:,:) = 0._wp 258 DO jk = 1, jpkm1 259 DO jj = 1, jpj 260 DO ji = 1, jpi ! vector opt. 261 z2d(ji,jj) = z2d(ji,jj) + rau0 * e3t(ji,jj,jk,Kmm) * ztpot(ji,jj,jk) 262 END DO 263 END DO 264 END DO 265 CALL iom_put( 'tosmint_pot', z2d ) 266 ENDIF 267 DEALLOCATE( ztpot ) 268 ENDIF 269 ELSE 270 IF( iom_use('ssttot') ) THEN ! Output sst in case we use EOS-80 271 zsst = glob_sum( 'diaar5', area(:,:) * ts(:,:,1,jp_tem,Kmm) ) 272 CALL iom_put('ssttot', zsst / area_tot ) 273 ENDIF 209 274 ENDIF 210 275 211 276 IF( iom_use( 'tnpeo' )) THEN 212 ! Work done against stratification by vertical mixing213 ! Exclude points where rn2 is negative as convection kicks in here and214 ! work is not being done against stratification277 ! Work done against stratification by vertical mixing 278 ! Exclude points where rn2 is negative as convection kicks in here and 279 ! work is not being done against stratification 215 280 ALLOCATE( zpe(jpi,jpj) ) 216 281 zpe(:,:) = 0._wp … … 220 285 DO ji = 1, jpi 221 286 IF( rn2(ji,jj,jk) > 0._wp ) THEN 222 zrw = ( gdepw(ji,jj,jk ,Kmm) - gdept(ji,jj,jk,Kmm) ) & 223 & / ( gdept(ji,jj,jk-1,Kmm) - gdept(ji,jj,jk,Kmm) ) 224 !!gm this can be reduced to : (depw-dept) / e3w (NB idem dans bn2 !) 225 ! zrw = ( gdept(ji,jj,jk,Kmm) - gdepw(ji,jj,jk,Kmm) ) / e3w(ji,jj,jk,Kmm) 226 !!gm end 287 zrw = ( gdept(ji,jj,jk,Kmm) - gdepw(ji,jj,jk,Kmm) ) / e3w(ji,jj,jk,Kmm) 227 288 ! 228 289 zaw = rab_n(ji,jj,jk,jp_tem) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_tem)* zrw 229 290 zbw = rab_n(ji,jj,jk,jp_sal) * (1. - zrw) + rab_n(ji,jj,jk-1,jp_sal)* zrw 230 291 ! 231 zpe(ji, jj) = zpe(ji, jj)&292 zpe(ji, jj) = zpe(ji,jj) & 232 293 & - grav * ( avt(ji,jj,jk) * zaw * (ts(ji,jj,jk-1,jp_tem,Kmm) - ts(ji,jj,jk,jp_tem,Kmm) ) & 233 294 & - avs(ji,jj,jk) * zbw * (ts(ji,jj,jk-1,jp_sal,Kmm) - ts(ji,jj,jk,jp_sal,Kmm) ) ) … … 240 301 DO ji = 1, jpi 241 302 DO jj = 1, jpj 242 zpe(ji,jj) = zpe(ji,jj) + avt(ji, jj, jk) * MIN(0._wp,rn2(ji, jj, jk)) * rau0 * e3w(ji, jj,jk,Kmm)303 zpe(ji,jj) = zpe(ji,jj) + avt(ji,jj,jk) * MIN(0._wp,rn2(ji,jj,jk)) * rau0 * e3w(ji,jj,jk,Kmm) 243 304 END DO 244 305 END DO 245 306 END DO 246 307 ENDIF 247 !!gm useless lbc_lnk since the computation above is performed over 1:jpi & 1:jpj248 !!gm CALL lbc_lnk( 'diaar5', zpe, 'T', 1._wp)249 308 CALL iom_put( 'tnpeo', zpe ) 250 309 DEALLOCATE( zpe ) … … 252 311 253 312 IF( l_ar5 ) THEN 254 DEALLOCATE( zarea_ssh , zbotpres )313 DEALLOCATE( zarea_ssh , zbotpres, z2d ) 255 314 DEALLOCATE( zrhd , zrhop ) 256 315 DEALLOCATE( ztsn ) … … 288 347 CALL lbc_lnk( 'diaar5', z2d, 'U', -1. ) 289 348 IF( cptr == 'adv' ) THEN 290 IF( ktra == jp_tem ) CALL iom_put( "uadv_heattr", rau0_rcp * z2d ) ! advective heat transport in i-direction291 IF( ktra == jp_sal ) CALL iom_put( "uadv_salttr", rau0 * z2d ) ! advective salt transport in i-direction349 IF( ktra == jp_tem ) CALL iom_put( 'uadv_heattr' , rau0_rcp * z2d ) ! advective heat transport in i-direction 350 IF( ktra == jp_sal ) CALL iom_put( 'uadv_salttr' , rau0 * z2d ) ! advective salt transport in i-direction 292 351 ENDIF 293 352 IF( cptr == 'ldf' ) THEN 294 IF( ktra == jp_tem ) CALL iom_put( "udiff_heattr", rau0_rcp * z2d ) ! diffusive heat transport in i-direction295 IF( ktra == jp_sal ) CALL iom_put( "udiff_salttr", rau0 * z2d ) ! diffusive salt transport in i-direction353 IF( ktra == jp_tem ) CALL iom_put( 'udiff_heattr' , rau0_rcp * z2d ) ! diffusive heat transport in i-direction 354 IF( ktra == jp_sal ) CALL iom_put( 'udiff_salttr' , rau0 * z2d ) ! diffusive salt transport in i-direction 296 355 ENDIF 297 356 ! … … 306 365 CALL lbc_lnk( 'diaar5', z2d, 'V', -1. ) 307 366 IF( cptr == 'adv' ) THEN 308 IF( ktra == jp_tem ) CALL iom_put( "vadv_heattr", rau0_rcp * z2d ) ! advective heat transport in j-direction309 IF( ktra == jp_sal ) CALL iom_put( "vadv_salttr", rau0 * z2d ) ! advective salt transport in j-direction367 IF( ktra == jp_tem ) CALL iom_put( 'vadv_heattr' , rau0_rcp * z2d ) ! advective heat transport in j-direction 368 IF( ktra == jp_sal ) CALL iom_put( 'vadv_salttr' , rau0 * z2d ) ! advective salt transport in j-direction 310 369 ENDIF 311 370 IF( cptr == 'ldf' ) THEN 312 IF( ktra == jp_tem ) CALL iom_put( "vdiff_heattr", rau0_rcp * z2d ) ! diffusive heat transport in j-direction313 IF( ktra == jp_sal ) CALL iom_put( "vdiff_salttr", rau0 * z2d ) ! diffusive salt transport in j-direction371 IF( ktra == jp_tem ) CALL iom_put( 'vdiff_heattr' , rau0_rcp * z2d ) ! diffusive heat transport in j-direction 372 IF( ktra == jp_sal ) CALL iom_put( 'vdiff_salttr' , rau0 * z2d ) ! diffusive salt transport in j-direction 314 373 ENDIF 315 374 … … 324 383 !!---------------------------------------------------------------------- 325 384 INTEGER :: inum 326 INTEGER :: ik 385 INTEGER :: ik, idep 327 386 INTEGER :: ji, jj, jk ! dummy loop indices 328 387 REAL(wp) :: zztmp 329 388 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: zsaldta ! Jan/Dec levitus salinity 389 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zvol0 330 390 ! 331 391 !!---------------------------------------------------------------------- … … 341 401 IF( dia_ar5_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'dia_ar5_init : unable to allocate arrays' ) 342 402 343 area(:,:) = e1e2t(:,:) * tmask_i(:,:)344 345 area_tot = SUM( area(:,:) ) ; CALL mpp_sum( 'diaar5', area_tot ) 346 347 vol0= 0._wp403 area(:,:) = e1e2t(:,:) 404 area_tot = glob_sum( 'diaar5', area(:,:) ) 405 406 ALLOCATE( zvol0(jpi,jpj) ) 407 zvol0 (:,:) = 0._wp 348 408 thick0(:,:) = 0._wp 349 409 DO jk = 1, jpkm1 350 vol0 = vol0 + SUM( area (:,:) * tmask(:,:,jk) * e3t_0(:,:,jk) ) 351 thick0(:,:) = thick0(:,:) + tmask_i(:,:) * tmask(:,:,jk) * e3t_0(:,:,jk) 352 END DO 353 CALL mpp_sum( 'diaar5', vol0 ) 410 DO jj = 1, jpj ! interpolation of salinity at the last ocean level (i.e. the partial step) 411 DO ji = 1, jpi 412 idep = tmask(ji,jj,jk) * e3t_0(ji,jj,jk) 413 zvol0 (ji,jj) = zvol0 (ji,jj) + idep * area(ji,jj) 414 thick0(ji,jj) = thick0(ji,jj) + idep 415 END DO 416 END DO 417 END DO 418 vol0 = glob_sum( 'diaar5', zvol0 ) 419 DEALLOCATE( zvol0 ) 354 420 355 421 IF( iom_use( 'sshthster' ) ) THEN 356 ALLOCATE( zsaldta(jpi,jpj,jp j,jpts) )422 ALLOCATE( zsaldta(jpi,jpj,jpk,jpts) ) 357 423 CALL iom_open ( 'sali_ref_clim_monthly', inum ) 358 424 CALL iom_get ( inum, jpdom_data, 'vosaline' , zsaldta(:,:,:,1), 1 )
Note: See TracChangeset
for help on using the changeset viewer.