Changeset 5034 for branches/2014/dev_r4650_UKMO14.12_STAND_ALONE_OBSOPER/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90
- Timestamp:
- 2015-01-15T14:48:42+01:00 (9 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r4650_UKMO14.12_STAND_ALONE_OBSOPER/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90
r4624 r5034 9 9 10 10 !!---------------------------------------------------------------------- 11 !! dia_hsb : Diagnose the conservation of ocean heat and salt contents, and volume 12 !! dia_hsb_rst : Read or write DIA file in restart file 13 !! dia_hsb_init : Initialization of the conservation diagnostic 14 !!---------------------------------------------------------------------- 11 15 USE oce ! ocean dynamics and tracers 12 16 USE dom_oce ! ocean space and time domain 13 17 USE phycst ! physical constants 14 18 USE sbc_oce ! surface thermohaline fluxes 15 USE in_out_manager ! I/O manager 19 USE sbcrnf ! river runoff 20 USE sbcisf ! ice shelves 16 21 USE domvvl ! vertical scale factors 17 22 USE traqsr ! penetrative solar radiation 18 23 USE trabbc ! bottom boundary condition 19 USE lib_mpp ! distributed memory computing library20 24 USE trabbc ! bottom boundary condition 21 25 USE bdy_par ! (for lk_bdy) 26 USE restart ! ocean restart 27 ! 28 USE iom ! I/O manager 29 USE in_out_manager ! I/O manager 30 USE lib_fortran ! glob_sum 31 USE lib_mpp ! distributed memory computing library 22 32 USE timing ! preformance summary 23 USE iom ! I/O manager 24 USE lib_fortran ! glob_sum 25 USE restart ! ocean restart 26 USE wrk_nemo ! work arrays 27 USE sbcrnf ! river runoffd 33 USE wrk_nemo ! work arrays 28 34 29 35 IMPLICIT NONE … … 36 42 LOGICAL, PUBLIC :: ln_diahsb !: check the heat and salt budgets 37 43 38 REAL(dp) :: surf_tot ! 39 REAL(dp) :: frc_t , frc_s , frc_v ! global forcing trends 40 REAL(dp) :: frc_wn_t , frc_wn_s ! global forcing trends 41 REAL(dp), DIMENSION(:,:) , ALLOCATABLE :: surf , ssh_ini ! 42 REAL(dp), DIMENSION(:,:,:), ALLOCATABLE :: hc_loc_ini, sc_loc_ini, e3t_ini ! 43 REAL(dp), DIMENSION(:,:) , ALLOCATABLE :: ssh_hc_loc_ini, ssh_sc_loc_ini 44 REAL(wp) :: surf_tot ! ocean surface 45 REAL(wp) :: frc_t, frc_s, frc_v ! global forcing trends 46 REAL(wp) :: frc_wn_t, frc_wn_s ! global forcing trends 47 ! 48 REAL(wp), DIMENSION(:,:) , ALLOCATABLE :: surf , ssh_ini ! 49 REAL(wp), DIMENSION(:,:) , ALLOCATABLE :: ssh_hc_loc_ini, ssh_sc_loc_ini ! 50 REAL(wp), DIMENSION(:,:,:), ALLOCATABLE :: hc_loc_ini, sc_loc_ini, e3t_ini ! 44 51 45 52 !! * Substitutions 46 53 # include "domzgr_substitute.h90" 47 54 # include "vectopt_loop_substitute.h90" 48 49 55 !!---------------------------------------------------------------------- 50 56 !! NEMO/OPA 3.3 , NEMO Consortium (2010) … … 52 58 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 53 59 !!---------------------------------------------------------------------- 54 55 60 CONTAINS 56 61 … … 67 72 !!--------------------------------------------------------------------------- 68 73 INTEGER, INTENT(in) :: kt ! ocean time-step index 69 !! 70 INTEGER :: jk ! dummy loop indice 71 REAL(dp) :: zdiff_hc , zdiff_sc ! heat and salt content variations 72 REAL(dp) :: zdiff_hc1 , zdiff_sc1 ! - - - - - - - - 73 REAL(dp) :: zdiff_v1 , zdiff_v2 ! volume variation 74 REAL(dp) :: zerr_hc1 , zerr_sc1 ! heat and salt content misfit 75 REAL(dp) :: zvol_tot ! volume 76 REAL(dp) :: z_frc_trd_t , z_frc_trd_s ! - - 77 REAL(dp) :: z_frc_trd_v ! - - 78 REAL(dp) :: z_wn_trd_t , z_wn_trd_s ! - - 79 REAL(dp) :: z_ssh_hc , z_ssh_sc ! - - 74 ! 75 INTEGER :: ji, jj, jk ! dummy loop indice 76 REAL(wp) :: zdiff_hc , zdiff_sc ! heat and salt content variations 77 REAL(wp) :: zdiff_hc1 , zdiff_sc1 ! - - - - 78 REAL(wp) :: zdiff_v1 , zdiff_v2 ! volume variation 79 REAL(wp) :: zerr_hc1 , zerr_sc1 ! heat and salt content misfit 80 REAL(wp) :: zvol_tot ! volume 81 REAL(wp) :: z_frc_trd_t , z_frc_trd_s ! - - 82 REAL(wp) :: z_frc_trd_v ! - - 83 REAL(wp) :: z_wn_trd_t , z_wn_trd_s ! - - 84 REAL(wp) :: z_ssh_hc , z_ssh_sc ! - - 85 REAL(wp), DIMENSION(:,:), POINTER :: z2d0, z2d1 80 86 !!--------------------------------------------------------------------------- 81 87 IF( nn_timing == 1 ) CALL timing_start('dia_hsb') 82 88 CALL wrk_alloc( jpi,jpj, z2d0, z2d1 ) 89 ! 90 tsn(:,:,:,1) = tsn(:,:,:,1) * tmask(:,:,:) ; tsb(:,:,:,1) = tsb(:,:,:,1) * tmask(:,:,:) ; 91 tsn(:,:,:,2) = tsn(:,:,:,2) * tmask(:,:,:) ; tsb(:,:,:,2) = tsb(:,:,:,2) * tmask(:,:,:) ; 83 92 ! ------------------------- ! 84 93 ! 1 - Trends due to forcing ! 85 94 ! ------------------------- ! 86 z_frc_trd_v = r1_rau0 * glob_sum( - ( emp(:,:) - rnf(:,:) ) * surf(:,:) ) ! volume fluxes87 z_frc_trd_t = glob_sum( sbc_tsc(:,:,jp_tem) * surf(:,:) ) ! heat fluxes88 z_frc_trd_s = glob_sum( sbc_tsc(:,:,jp_sal) * surf(:,:) ) ! salt fluxes95 z_frc_trd_v = r1_rau0 * glob_sum( - ( emp(:,:) - rnf(:,:) + rdivisf * fwfisf(:,:) ) * surf(:,:) ) ! volume fluxes 96 z_frc_trd_t = glob_sum( sbc_tsc(:,:,jp_tem) * surf(:,:) ) ! heat fluxes 97 z_frc_trd_s = glob_sum( sbc_tsc(:,:,jp_sal) * surf(:,:) ) ! salt fluxes 89 98 ! Add runoff heat & salt input 90 99 IF( ln_rnf ) z_frc_trd_t = z_frc_trd_t + glob_sum( rnf_tsc(:,:,jp_tem) * surf(:,:) ) 91 100 IF( ln_rnf_sal) z_frc_trd_s = z_frc_trd_s + glob_sum( rnf_tsc(:,:,jp_sal) * surf(:,:) ) 101 ! Add geothermal ice shelf 102 IF( nn_isf .GE. 1 ) THEN 103 z_frc_trd_t = z_frc_trd_t & 104 & + glob_sum( ( risf_tsc(:,:,jp_tem) - rdivisf * fwfisf(:,:) * (-1.9) * r1_rau0 ) * surf(:,:) ) 105 z_frc_trd_s = z_frc_trd_s + (1.0_wp - rdivisf) * glob_sum( risf_tsc(:,:,jp_sal) * surf(:,:) ) 106 ENDIF 92 107 93 108 ! Add penetrative solar radiation … … 97 112 ! 98 113 IF( .NOT. lk_vvl ) THEN 99 z_wn_trd_t = - glob_sum( surf(:,:) * wn(:,:,1) * tsb(:,:,1,jp_tem) ) 100 z_wn_trd_s = - glob_sum( surf(:,:) * wn(:,:,1) * tsb(:,:,1,jp_sal) ) 114 z2d0=0.0_wp ; z2d1=0.0_wp 115 DO ji=1,jpi 116 DO jj=1,jpj 117 z2d0(ji,jj) = surf(ji,jj) * wn(ji,jj,mikt(ji,jj)) * tsb(ji,jj,mikt(ji,jj),jp_tem) 118 z2d1(ji,jj) = surf(ji,jj) * wn(ji,jj,mikt(ji,jj)) * tsb(ji,jj,mikt(ji,jj),jp_sal) 119 ENDDO 120 ENDDO 121 z_wn_trd_t = - glob_sum( z2d0 ) 122 z_wn_trd_s = - glob_sum( z2d1 ) 101 123 ENDIF 102 124 … … 113 135 ! 2 - Content variations ! 114 136 ! ------------------------ ! 115 zdiff_v2 = 0. d0116 zdiff_hc = 0. d0117 zdiff_sc = 0. d0137 zdiff_v2 = 0._wp 138 zdiff_hc = 0._wp 139 zdiff_sc = 0._wp 118 140 119 141 ! volume variation (calculated with ssh) … … 122 144 ! heat & salt content variation (associated with ssh) 123 145 IF( .NOT. lk_vvl ) THEN 124 z_ssh_hc = glob_sum( surf(:,:) * ( tsn(:,:,1,jp_tem) * sshn(:,:) - ssh_hc_loc_ini(:,:) ) ) 125 z_ssh_sc = glob_sum( surf(:,:) * ( tsn(:,:,1,jp_sal) * sshn(:,:) - ssh_sc_loc_ini(:,:) ) ) 146 z2d0 = 0._wp ; z2d1 = 0._wp 147 DO ji = 1, jpi 148 DO jj = 1, jpj 149 z2d0(ji,jj) = surf(ji,jj) * ( tsn(ji,jj,mikt(ji,jj),jp_tem) * sshn(ji,jj) - ssh_hc_loc_ini(ji,jj) ) 150 z2d1(ji,jj) = surf(ji,jj) * ( tsn(ji,jj,mikt(ji,jj),jp_sal) * sshn(ji,jj) - ssh_sc_loc_ini(ji,jj) ) 151 END DO 152 END DO 153 z_ssh_hc = glob_sum( z2d0 ) 154 z_ssh_sc = glob_sum( z2d1 ) 126 155 ENDIF 127 156 … … 153 182 ! 3 - Diagnostics writing ! 154 183 ! ----------------------- ! 155 zvol_tot = 0.d0 ! total ocean volume184 zvol_tot = 0._wp ! total ocean volume (calculated with scale factors) 156 185 DO jk = 1, jpkm1 157 186 zvol_tot = zvol_tot + glob_sum( surf(:,:) * tmask(:,:,jk) * fse3t_n(:,:,jk) ) 158 187 END DO 188 189 !!gm to be added ? 190 ! IF( .NOT. lk_vvl ) THEN ! fixed volume, add the ssh contribution 191 ! zvol_tot = zvol_tot + glob_sum( surf(:,:) * sshn(:,:) ) 192 ! ENDIF 193 !!gm end 194 159 195 160 196 IF( lk_vvl ) THEN … … 183 219 IF( lrst_oce ) CALL dia_hsb_rst( kt, 'WRITE' ) 184 220 221 CALL wrk_dealloc( jpi,jpj, z2d0, z2d1 ) 222 185 223 IF( nn_timing == 1 ) CALL timing_stop('dia_hsb') 186 !224 ! 187 225 END SUBROUTINE dia_hsb 188 226 189 190 SUBROUTINE dia_hsb_init191 !!---------------------------------------------------------------------------192 !! *** ROUTINE dia_hsb ***193 !!194 !! ** Purpose: Initialization for the heat salt volume budgets195 !!196 !! ** Method : Compute initial heat content, salt content and volume197 !!198 !! ** Action : - Compute initial heat content, salt content and volume199 !! - Initialize forcing trends200 !! - Compute coefficients for conversion201 !!---------------------------------------------------------------------------202 INTEGER :: jk ! dummy loop indice203 INTEGER :: ierror ! local integer204 !!205 NAMELIST/namhsb/ ln_diahsb206 !207 INTEGER :: ios208 !!----------------------------------------------------------------------209 210 IF(lwp) THEN211 WRITE(numout,*)212 WRITE(numout,*) 'dia_hsb_init : check the heat and salt budgets'213 WRITE(numout,*) '~~~~~~~~ '214 ENDIF215 216 REWIND( numnam_ref ) ! Namelist namhsb in reference namelist217 READ ( numnam_ref, namhsb, IOSTAT = ios, ERR = 901)218 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namhsb in reference namelist', lwp )219 220 REWIND( numnam_cfg ) ! Namelist namhsb in configuration namelist221 READ ( numnam_cfg, namhsb, IOSTAT = ios, ERR = 902 )222 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namhsb in configuration namelist', lwp )223 IF(lwm) WRITE ( numond, namhsb )224 225 !226 IF(lwp) THEN ! Control print227 WRITE(numout,*)228 WRITE(numout,*) 'dia_hsb_init : check the heat and salt budgets'229 WRITE(numout,*) '~~~~~~~~~~~~'230 WRITE(numout,*) ' Namelist namhsb : set hsb parameters'231 WRITE(numout,*) ' Switch for hsb diagnostic (T) or not (F) ln_diahsb = ', ln_diahsb232 WRITE(numout,*)233 ENDIF234 235 IF( .NOT. ln_diahsb ) RETURN236 ! IF( .NOT. lk_mpp_rep ) &237 ! CALL ctl_stop (' Your global mpp_sum if performed in single precision - 64 bits -', &238 ! & ' whereas the global sum to be precise must be done in double precision ',&239 ! & ' please add key_mpp_rep')240 241 ! ------------------- !242 ! 1 - Allocate memory !243 ! ------------------- !244 ALLOCATE( hc_loc_ini(jpi,jpj,jpk), sc_loc_ini(jpi,jpj,jpk), &245 & e3t_ini(jpi,jpj,jpk), surf(jpi,jpj), ssh_ini(jpi,jpj), STAT=ierror )246 IF( ierror > 0 ) THEN247 CALL ctl_stop( 'dia_hsb: unable to allocate hc_loc_ini' ) ; RETURN248 ENDIF249 250 IF(.NOT. lk_vvl ) ALLOCATE( ssh_hc_loc_ini(jpi,jpj), ssh_sc_loc_ini(jpi,jpj),STAT=ierror )251 IF( ierror > 0 ) THEN252 CALL ctl_stop( 'dia_hsb: unable to allocate hc_loc_ini' ) ; RETURN253 ENDIF254 255 ! ----------------------------------------------- !256 ! 2 - Time independant variables and file opening !257 ! ----------------------------------------------- !258 IF(lwp) WRITE(numout,*) "dia_hsb: heat salt volume budgets activated"259 IF(lwp) WRITE(numout,*) '~~~~~~~'260 surf(:,:) = e1t(:,:) * e2t(:,:) * tmask(:,:,1) * tmask_i(:,:) ! masked surface grid cell area261 surf_tot = glob_sum( surf(:,:) ) ! total ocean surface area262 263 IF( lk_bdy ) CALL ctl_warn( 'dia_hsb does not take open boundary fluxes into account' )264 !265 ! ---------------------------------- !266 ! 4 - initial conservation variables !267 ! ---------------------------------- !268 CALL dia_hsb_rst( nit000, 'READ' ) !* read or initialize all required files269 !270 END SUBROUTINE dia_hsb_init271 227 272 228 SUBROUTINE dia_hsb_rst( kt, cdrw ) … … 281 237 CHARACTER(len=*), INTENT(in) :: cdrw ! "READ"/"WRITE" flag 282 238 ! 283 INTEGER :: j k !284 INTEGER :: id1 ! local integers239 INTEGER :: ji, jj, jk ! dummy loop indices 240 INTEGER :: id1 ! local integers 285 241 !!---------------------------------------------------------------------- 286 242 ! … … 317 273 sc_loc_ini(:,:,jk) = tsn(:,:,jk,jp_sal) * fse3t_n(:,:,jk) ! initial salt content 318 274 END DO 319 frc_v = 0. d0! volume trend due to forcing320 frc_t = 0. d0! heat content - - - -321 frc_s = 0. d0! salt content - - - -275 frc_v = 0._wp ! volume trend due to forcing 276 frc_t = 0._wp ! heat content - - - - 277 frc_s = 0._wp ! salt content - - - - 322 278 IF( .NOT. lk_vvl ) THEN 323 ssh_hc_loc_ini(:,:) = tsn(:,:,1,jp_tem) * sshn(:,:) ! initial heat content in ssh 324 ssh_sc_loc_ini(:,:) = tsn(:,:,1,jp_sal) * sshn(:,:) ! initial salt content in ssh 325 frc_wn_t = 0.d0 ! initial heat content misfit due to free surface 326 frc_wn_s = 0.d0 ! initial salt content misfit due to free surface 279 DO ji=1,jpi 280 DO jj=1,jpj 281 ssh_hc_loc_ini(ji,jj) = tsn(ji,jj,mikt(ji,jj),jp_tem) * sshn(ji,jj) ! initial heat content in ssh 282 ssh_sc_loc_ini(ji,jj) = tsn(ji,jj,mikt(ji,jj),jp_sal) * sshn(ji,jj) ! initial salt content in ssh 283 ENDDO 284 ENDDO 285 frc_wn_t = 0._wp ! initial heat content misfit due to free surface 286 frc_wn_s = 0._wp ! initial salt content misfit due to free surface 327 287 ENDIF 328 288 ENDIF … … 354 314 END SUBROUTINE dia_hsb_rst 355 315 316 317 SUBROUTINE dia_hsb_init 318 !!--------------------------------------------------------------------------- 319 !! *** ROUTINE dia_hsb *** 320 !! 321 !! ** Purpose: Initialization for the heat salt volume budgets 322 !! 323 !! ** Method : Compute initial heat content, salt content and volume 324 !! 325 !! ** Action : - Compute initial heat content, salt content and volume 326 !! - Initialize forcing trends 327 !! - Compute coefficients for conversion 328 !!--------------------------------------------------------------------------- 329 INTEGER :: jk ! dummy loop indice 330 INTEGER :: ierror ! local integer 331 INTEGER :: ios 332 ! 333 NAMELIST/namhsb/ ln_diahsb 334 !!---------------------------------------------------------------------- 335 336 IF(lwp) THEN 337 WRITE(numout,*) 338 WRITE(numout,*) 'dia_hsb_init : check the heat and salt budgets' 339 WRITE(numout,*) '~~~~~~~~ ' 340 ENDIF 341 342 REWIND( numnam_ref ) ! Namelist namhsb in reference namelist 343 READ ( numnam_ref, namhsb, IOSTAT = ios, ERR = 901) 344 901 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namhsb in reference namelist', lwp ) 345 346 REWIND( numnam_cfg ) ! Namelist namhsb in configuration namelist 347 READ ( numnam_cfg, namhsb, IOSTAT = ios, ERR = 902 ) 348 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namhsb in configuration namelist', lwp ) 349 IF(lwm) WRITE ( numond, namhsb ) 350 351 ! 352 IF(lwp) THEN ! Control print 353 WRITE(numout,*) 354 WRITE(numout,*) 'dia_hsb_init : check the heat and salt budgets' 355 WRITE(numout,*) '~~~~~~~~~~~~' 356 WRITE(numout,*) ' Namelist namhsb : set hsb parameters' 357 WRITE(numout,*) ' Switch for hsb diagnostic (T) or not (F) ln_diahsb = ', ln_diahsb 358 WRITE(numout,*) 359 ENDIF 360 361 IF( .NOT. ln_diahsb ) RETURN 362 ! IF( .NOT. lk_mpp_rep ) & 363 ! CALL ctl_stop (' Your global mpp_sum if performed in single precision - 64 bits -', & 364 ! & ' whereas the global sum to be precise must be done in double precision ',& 365 ! & ' please add key_mpp_rep') 366 367 ! ------------------- ! 368 ! 1 - Allocate memory ! 369 ! ------------------- ! 370 ALLOCATE( hc_loc_ini(jpi,jpj,jpk), sc_loc_ini(jpi,jpj,jpk), & 371 & e3t_ini(jpi,jpj,jpk), surf(jpi,jpj), ssh_ini(jpi,jpj), STAT=ierror ) 372 IF( ierror > 0 ) THEN 373 CALL ctl_stop( 'dia_hsb: unable to allocate hc_loc_ini' ) ; RETURN 374 ENDIF 375 376 IF(.NOT. lk_vvl ) ALLOCATE( ssh_hc_loc_ini(jpi,jpj), ssh_sc_loc_ini(jpi,jpj),STAT=ierror ) 377 IF( ierror > 0 ) THEN 378 CALL ctl_stop( 'dia_hsb: unable to allocate hc_loc_ini' ) ; RETURN 379 ENDIF 380 381 ! ----------------------------------------------- ! 382 ! 2 - Time independant variables and file opening ! 383 ! ----------------------------------------------- ! 384 IF(lwp) WRITE(numout,*) "dia_hsb: heat salt volume budgets activated" 385 IF(lwp) WRITE(numout,*) '~~~~~~~' 386 surf(:,:) = e1t(:,:) * e2t(:,:) * tmask_i(:,:) ! masked surface grid cell area 387 surf_tot = glob_sum( surf(:,:) ) ! total ocean surface area 388 389 IF( lk_bdy ) CALL ctl_warn( 'dia_hsb does not take open boundary fluxes into account' ) 390 ! 391 ! ---------------------------------- ! 392 ! 4 - initial conservation variables ! 393 ! ---------------------------------- ! 394 CALL dia_hsb_rst( nit000, 'READ' ) !* read or initialize all required files 395 ! 396 END SUBROUTINE dia_hsb_init 397 356 398 !!====================================================================== 357 399 END MODULE diahsb
Note: See TracChangeset
for help on using the changeset viewer.
