Changeset 4077
- Timestamp:
- 2013-10-18T11:23:22+02:00 (11 years ago)
- Location:
- branches/2013/dev_r3940_CNRS4_IOCRS/NEMOGCM/NEMO
- Files:
-
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2013/dev_r3940_CNRS4_IOCRS/NEMOGCM/NEMO/OPA_SRC/DIA/diahsb.F90
r3625 r4077 21 21 USE bdy_par ! (for lk_bdy) 22 22 USE timing ! preformance summary 23 USE lib_fortran 24 USE sbcrnf 23 25 24 26 IMPLICIT NONE … … 33 35 REAL(dp) :: surf_tot , vol_tot ! 34 36 REAL(dp) :: frc_t , frc_s , frc_v ! global forcing trends 37 REAL(dp) :: frc_wn_t , frc_wn_s ! global forcing trends 35 38 REAL(dp) :: fact1 ! conversion factors 36 39 REAL(dp) :: fact21 , fact22 ! - - … … 38 41 REAL(dp), DIMENSION(:,:) , ALLOCATABLE :: surf , ssh_ini ! 39 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 40 44 41 45 !! * Substitutions … … 67 71 INTEGER :: jk ! dummy loop indice 68 72 REAL(dp) :: zdiff_hc , zdiff_sc ! heat and salt content variations 73 REAL(dp) :: zdiff_hc1 , zdiff_sc1 ! heat and salt content variations of ssh 69 74 REAL(dp) :: zdiff_v1 , zdiff_v2 ! volume variation 70 REAL(dp) :: z 1_rau0 ! local scalars75 REAL(dp) :: zerr_hc1 , zerr_sc1 ! Non conservation due to free surface 71 76 REAL(dp) :: zdeltat ! - - 72 77 REAL(dp) :: z_frc_trd_t , z_frc_trd_s ! - - 73 78 REAL(dp) :: z_frc_trd_v ! - - 79 REAL(dp) :: z_wn_trd_t , z_wn_trd_s ! - - 80 REAL(dp) :: z_ssh_hc , z_ssh_sc ! - - 74 81 !!--------------------------------------------------------------------------- 75 82 IF( nn_timing == 1 ) CALL timing_start('dia_hsb') … … 78 85 ! 1 - Trends due to forcing ! 79 86 ! ------------------------- ! 80 z1_rau0 = 1.e0 / rau0 81 z_frc_trd_v = z1_rau0 * SUM( - ( emp(:,:) - rnf(:,:) ) * surf(:,:) ) ! volume fluxes 82 z_frc_trd_t = SUM( sbc_tsc(:,:,jp_tem) * surf(:,:) ) ! heat fluxes 83 z_frc_trd_s = SUM( sbc_tsc(:,:,jp_sal) * surf(:,:) ) ! salt fluxes 87 z_frc_trd_v = r1_rau0 * glob_sum( - ( emp(:,:) - rnf(:,:) ) * surf(:,:) ) ! volume fluxes 88 z_frc_trd_t = glob_sum( sbc_tsc(:,:,jp_tem) * surf(:,:) ) ! heat fluxes 89 z_frc_trd_s = glob_sum( sbc_tsc(:,:,jp_sal) * surf(:,:) ) ! salt fluxes 90 ! Add runoff heat & salt input 91 IF( ln_rnf ) z_frc_trd_t = z_frc_trd_t + glob_sum( rnf_tsc(:,:,jp_tem) * surf(:,:) ) 92 IF( ln_rnf_sal) z_frc_trd_s = z_frc_trd_s + glob_sum( rnf_tsc(:,:,jp_sal) * surf(:,:) ) 84 93 ! Add penetrative solar radiation 85 IF( ln_traqsr ) z_frc_trd_t = z_frc_trd_t + r1_rau0_rcp * SUM( qsr (:,:) * surf(:,:) )94 IF( ln_traqsr ) z_frc_trd_t = z_frc_trd_t + r1_rau0_rcp * glob_sum( qsr (:,:) * surf(:,:) ) 86 95 ! Add geothermal heat flux 87 IF( ln_trabbc ) z_frc_trd_t = z_frc_trd_t + r1_rau0_rcp * SUM( qgh_trd0(:,:) * surf(:,:) ) 88 IF( lk_mpp ) THEN 89 CALL mpp_sum( z_frc_trd_v ) 90 CALL mpp_sum( z_frc_trd_t ) 91 ENDIF 96 IF( ln_trabbc ) z_frc_trd_t = z_frc_trd_t + glob_sum( qgh_trd0(:,:) * surf(:,:) ) 97 IF( .NOT. lk_vvl ) THEN 98 z_wn_trd_t = - glob_sum( surf(:,:) * wn(:,:,1) * tsb(:,:,1,jp_tem) ) 99 z_wn_trd_s = - glob_sum( surf(:,:) * wn(:,:,1) * tsb(:,:,1,jp_sal) ) 100 ENDIF 101 92 102 frc_v = frc_v + z_frc_trd_v * rdt 93 103 frc_t = frc_t + z_frc_trd_t * rdt 94 104 frc_s = frc_s + z_frc_trd_s * rdt 105 ! ! Advection flux through fixed surface (z=0) 106 IF( .NOT. lk_vvl ) THEN 107 frc_wn_t = frc_wn_t + z_wn_trd_t * rdt 108 frc_wn_s = frc_wn_s + z_wn_trd_s * rdt 109 ENDIF 95 110 96 111 ! ----------------------- ! … … 100 115 zdiff_hc = 0.d0 101 116 zdiff_sc = 0.d0 117 102 118 ! volume variation (calculated with ssh) 103 zdiff_v1 = SUM( surf(:,:) * tmask(:,:,1) * ( sshn(:,:) - ssh_ini(:,:) ) ) 119 zdiff_v1 = glob_sum( surf(:,:) * ( sshn(:,:) - ssh_ini(:,:) ) ) 120 121 ! heat & salt content variation (associated with ssh) 122 IF( .NOT. lk_vvl ) THEN 123 z_ssh_hc = glob_sum( surf(:,:) * ( tsn(:,:,1,jp_tem) * sshn(:,:) - ssh_hc_loc_ini(:,:) ) ) 124 z_ssh_sc = glob_sum( surf(:,:) * ( tsn(:,:,1,jp_sal) * sshn(:,:) - ssh_sc_loc_ini(:,:) ) ) 125 ENDIF 126 104 127 DO jk = 1, jpkm1 105 106 zdiff_v2 = zdiff_v2 + SUM( surf(:,:) * tmask(:,:,jk) &128 ! volume variation (calculated with scale factors) 129 zdiff_v2 = zdiff_v2 + glob_sum( surf(:,:) * tmask(:,:,jk) & 107 130 & * ( fse3t_n(:,:,jk) & 108 131 & - e3t_ini(:,:,jk) ) ) 109 132 ! heat content variation 110 zdiff_hc = zdiff_hc + SUM( surf(:,:) * tmask(:,:,jk) &133 zdiff_hc = zdiff_hc + glob_sum( surf(:,:) * tmask(:,:,jk) & 111 134 & * ( fse3t_n(:,:,jk) * tsn(:,:,jk,jp_tem) & 112 135 & - hc_loc_ini(:,:,jk) ) ) 113 136 ! salt content variation 114 zdiff_sc = zdiff_sc + SUM( surf(:,:) * tmask(:,:,jk) &137 zdiff_sc = zdiff_sc + glob_sum( surf(:,:) * tmask(:,:,jk) & 115 138 & * ( fse3t_n(:,:,jk) * tsn(:,:,jk,jp_sal) & 116 139 & - sc_loc_ini(:,:,jk) ) ) 117 140 ENDDO 118 141 119 IF( lk_mpp ) THEN120 CALL mpp_sum( zdiff_hc )121 CALL mpp_sum( zdiff_sc )122 CALL mpp_sum( zdiff_v1 )123 CALL mpp_sum( zdiff_v2 )124 ENDIF125 126 142 ! Substract forcing from heat content, salt content and volume variations 127 143 zdiff_v1 = zdiff_v1 - frc_v 128 zdiff_v2 = zdiff_v2 - frc_v144 IF( lk_vvl ) zdiff_v2 = zdiff_v2 - frc_v 129 145 zdiff_hc = zdiff_hc - frc_t 130 146 zdiff_sc = zdiff_sc - frc_s 147 IF( .NOT. lk_vvl ) THEN 148 zdiff_hc1 = zdiff_hc + z_ssh_hc 149 zdiff_sc1 = zdiff_sc + z_ssh_sc 150 zerr_hc1 = z_ssh_hc - frc_wn_t 151 zerr_sc1 = z_ssh_sc - frc_wn_s 152 ENDIF 131 153 132 154 ! ----------------------- ! … … 134 156 ! ----------------------- ! 135 157 zdeltat = 1.e0 / ( ( kt - nit000 + 1 ) * rdt ) 136 WRITE(numhsb , 9020) kt , zdiff_hc / vol_tot , zdiff_hc * fact1 * zdeltat, & 137 & zdiff_sc / vol_tot , zdiff_sc * fact21 * zdeltat, zdiff_sc * fact22 * zdeltat, & 138 & zdiff_v1 , zdiff_v1 * fact31 * zdeltat, zdiff_v1 * fact32 * zdeltat, & 139 & zdiff_v2 , zdiff_v2 * fact31 * zdeltat, zdiff_v2 * fact32 * zdeltat 158 IF( lk_vvl ) THEN 159 WRITE(numhsb , 9020) kt , zdiff_hc / vol_tot , zdiff_hc * fact1 * zdeltat, & 160 & zdiff_sc / vol_tot , zdiff_sc * fact21 * zdeltat, zdiff_sc * fact22 * zdeltat, & 161 & zdiff_v1 , zdiff_v1 * fact31 * zdeltat, zdiff_v1 * fact32 * zdeltat, & 162 & zdiff_v2 , zdiff_v2 * fact31 * zdeltat, zdiff_v2 * fact32 * zdeltat 163 ELSE 164 WRITE(numhsb , 9030) kt , zdiff_hc1 / vol_tot , zdiff_hc1 * fact1 * zdeltat, & 165 & zdiff_sc1 / vol_tot , zdiff_sc1 * fact21 * zdeltat, zdiff_sc1 * fact22 * zdeltat, & 166 & zdiff_v1 , zdiff_v1 * fact31 * zdeltat, zdiff_v1 * fact32 * zdeltat, & 167 & zerr_hc1 / vol_tot , zerr_sc1 / vol_tot 168 ENDIF 140 169 141 170 IF ( kt == nitend ) CLOSE( numhsb ) … … 144 173 145 174 9020 FORMAT(I5,11D15.7) 175 9030 FORMAT(I5,10D15.7) 146 176 ! 147 177 END SUBROUTINE dia_hsb … … 179 209 180 210 IF( .NOT. ln_diahsb ) RETURN 211 IF( .NOT. lk_mpp_rep ) & 212 CALL ctl_stop (' Your global mpp_sum if performed in single precision - 64 bits -', & 213 & ' whereas the global sum to be precise must be done in double precision ',& 214 & ' please add key_mpp_rep') 181 215 182 216 ! ------------------- ! 183 217 ! 1 - Allocate memory ! 184 218 ! ------------------- ! 185 ALLOCATE( hc_loc_ini(jpi,jpj,jpk), STAT=ierror ) 219 ALLOCATE( hc_loc_ini(jpi,jpj,jpk), sc_loc_ini(jpi,jpj,jpk), & 220 & ssh_hc_loc_ini(jpi,jpj), ssh_sc_loc_ini(jpi,jpj), & 221 & e3t_ini(jpi,jpj,jpk) , & 222 & surf(jpi,jpj), ssh_ini(jpi,jpj), STAT=ierror ) 186 223 IF( ierror > 0 ) THEN 187 224 CALL ctl_stop( 'dia_hsb: unable to allocate hc_loc_ini' ) ; RETURN 188 ENDIF189 ALLOCATE( sc_loc_ini(jpi,jpj,jpk), STAT=ierror )190 IF( ierror > 0 ) THEN191 CALL ctl_stop( 'dia_hsb: unable to allocate sc_loc_ini' ) ; RETURN192 ENDIF193 ALLOCATE( e3t_ini(jpi,jpj,jpk) , STAT=ierror )194 IF( ierror > 0 ) THEN195 CALL ctl_stop( 'dia_hsb: unable to allocate e3t_ini' ) ; RETURN196 ENDIF197 ALLOCATE( surf(jpi,jpj) , STAT=ierror )198 IF( ierror > 0 ) THEN199 CALL ctl_stop( 'dia_hsb: unable to allocate surf' ) ; RETURN200 ENDIF201 ALLOCATE( ssh_ini(jpi,jpj) , STAT=ierror )202 IF( ierror > 0 ) THEN203 CALL ctl_stop( 'dia_hsb: unable to allocate ssh_ini' ) ; RETURN204 225 ENDIF 205 226 … … 214 235 cl_name = 'heat_salt_volume_budgets.txt' ! name of output file 215 236 surf(:,:) = e1t(:,:) * e2t(:,:) * tmask(:,:,1) * tmask_i(:,:) ! masked surface grid cell area 216 surf_tot = SUM( surf(:,:) ) ! total ocean surface area237 surf_tot = glob_sum( surf(:,:) ) ! total ocean surface area 217 238 vol_tot = 0.d0 ! total ocean volume 218 239 DO jk = 1, jpkm1 219 vol_tot = vol_tot + SUM( surf(:,:) * tmask(:,:,jk) &220 & * fse3t_n(:,:,jk) )240 vol_tot = vol_tot + glob_sum( surf(:,:) * tmask(:,:,jk) & 241 & * fse3t_n(:,:,jk) ) 221 242 END DO 222 IF( lk_mpp ) THEN223 CALL mpp_sum( vol_tot )224 CALL mpp_sum( surf_tot )225 ENDIF226 243 227 244 CALL ctl_opn( numhsb , cl_name , 'UNKNOWN' , 'FORMATTED' , 'SEQUENTIAL' , 1 , numout , lwp , 1 ) 228 ! 12345678901234567890123456789012345678901234567890123456789012345678901234567890 -> 80 229 WRITE( numhsb, 9010 ) "kt | heat content budget | salt content budget ", & 230 ! 123456789012345678901234567890123456789012345 -> 45 231 & "| volume budget (ssh) ", & 232 ! 678901234567890123456789012345678901234567890 -> 45 233 & "| volume budget (e3t) " 234 WRITE( numhsb, 9010 ) " | [C] [W/m2] | [psu] [mmm/s] [SV] ", & 235 & "| [m3] [mmm/s] [SV] ", & 236 & "| [m3] [mmm/s] [SV] " 237 245 IF( lk_vvl ) THEN 246 ! 12345678901234567890123456789012345678901234567890123456789012345678901234567890 -> 80 247 WRITE( numhsb, 9010 ) "kt | heat content budget | salt content budget ", & 248 ! 123456789012345678901234567890123456789012345 -> 45 249 & "| volume budget (ssh) ", & 250 ! 678901234567890123456789012345678901234567890 -> 45 251 & "| volume budget (e3t) " 252 WRITE( numhsb, 9010 ) " | [C] [W/m2] | [psu] [mmm/s] [SV] ", & 253 & "| [m3] [mmm/s] [SV] ", & 254 & "| [m3] [mmm/s] [SV] " 255 ELSE 256 ! 12345678901234567890123456789012345678901234567890123456789012345678901234567890 -> 80 257 WRITE( numhsb, 9011 ) "kt | heat content budget | salt content budget ", & 258 ! 123456789012345678901234567890123456789012345 -> 45 259 & "| volume budget (ssh) ", & 260 ! 678901234567890123456789012345678901234567890 -> 45 261 & "| Non conservation due to free surface " 262 WRITE( numhsb, 9011 ) " | [C] [W/m2] | [psu] [mmm/s] [SV] ", & 263 & "| [m3] [mmm/s] [SV] ", & 264 & "| [heat - C] [salt - psu] " 265 ENDIF 238 266 ! --------------- ! 239 267 ! 3 - Conversions ! (factors will be multiplied by duration afterwards) … … 261 289 frc_t = 0.d0 ! heat content - - - - 262 290 frc_s = 0.d0 ! salt content - - - - 291 IF( .NOT. lk_vvl ) THEN 292 ssh_hc_loc_ini(:,:) = tsn(:,:,1,jp_tem) * ssh_ini(:,:) ! initial heat content associated with ssh 293 ssh_sc_loc_ini(:,:) = tsn(:,:,1,jp_sal) * ssh_ini(:,:) ! initial salt content associated with ssh 294 frc_wn_t = 0.d0 295 frc_wn_s = 0.d0 296 ENDIF 263 297 ! 264 298 9010 FORMAT(A80,A45,A45) 299 9011 FORMAT(A80,A45,A45) 265 300 ! 266 301 END SUBROUTINE dia_hsb_init -
branches/2013/dev_r3940_CNRS4_IOCRS/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zsms.F90
r3882 r4077 76 76 ENDIF 77 77 ! 78 IF( ln_rsttr .AND. kt == nittrc000 ) CALL p4z_rst( nittrc000, 'READ' ) !* read or initialize all required fields 78 IF( kt == nittrc000 ) THEN 79 ! 80 CALL p4z_che ! initialize the chemical constants 81 ! 82 IF( .NOT. ln_rsttr ) THEN ; CALL p4z_ph_ini ! set PH at kt=nit000 83 ELSE ; CALL p4z_rst( nittrc000, 'READ' ) !* read or initialize all required fields 84 ENDIF 85 ! 86 ENDIF 87 ! 79 88 IF( ln_pisdmp .AND. MOD( kt - nn_dttrc, nn_pisdmp ) == 0 ) CALL p4z_dmp( kt ) ! Relaxation of some tracers 80 89 ! … … 236 245 ENDIF 237 246 238 END SUBROUTINE p4z_sms_init 247 END SUBROUTINE p4z_sms_inita 248 249 SUBROUTINE p4z_ph_ini 250 !!--------------------------------------------------------------------- 251 !! *** ROUTINE p4z_ini_ph *** 252 !! 253 !! ** Purpose : Initialization of chemical variables of the carbon cycle 254 !!--------------------------------------------------------------------- 255 INTEGER :: ji, jj, jk 256 REAL(wp) :: zcaralk, zbicarb, zco3 257 REAL(wp) :: ztmas, ztmas1 258 !!--------------------------------------------------------------------- 259 260 ! Set PH from total alkalinity, borat (???), akb3 (???) and ak23 (???) 261 ! -------------------------------------------------------- 262 DO jk = 1, jpk 263 DO jj = 1, jpj 264 DO ji = 1, jpi 265 ztmas = tmask(ji,jj,jk) 266 ztmas1 = 1. - tmask(ji,jj,jk) 267 zcaralk = trn(ji,jj,jk,jptal) - borat(ji,jj,jk) / ( 1. + 1.E-8 / ( rtrn + akb3(ji,jj,jk) ) ) 268 zco3 = ( zcaralk - trn(ji,jj,jk,jpdic) ) * ztmas + 0.5e-3 * ztmas1 269 zbicarb = ( 2. * trn(ji,jj,jk,jpdic) - zcaralk ) 270 hi(ji,jj,jk) = ( ak23(ji,jj,jk) * zbicarb / zco3 ) * ztmas + 1.e-9 * ztmas1 271 END DO 272 END DO 273 END DO 274 ! 275 END SUBROUTINE p4z_ph_ini 276 239 277 240 278 SUBROUTINE p4z_rst( kt, cdrw ) … … 266 304 ELSE 267 305 ! hi(:,:,:) = 1.e-9 268 ! Set PH from total alkalinity, borat (???), akb3 (???) and ak23 (???) 269 ! -------------------------------------------------------- 270 DO jk = 1, jpk 271 DO jj = 1, jpj 272 DO ji = 1, jpi 273 ztmas = tmask(ji,jj,jk) 274 ztmas1 = 1. - tmask(ji,jj,jk) 275 zcaralk = trn(ji,jj,jk,jptal) - borat(ji,jj,jk) / ( 1. + 1.E-8 / ( rtrn + akb3(ji,jj,jk) ) ) 276 zco3 = ( zcaralk - trn(ji,jj,jk,jpdic) ) * ztmas + 0.5e-3 * ztmas1 277 zbicarb = ( 2. * trn(ji,jj,jk,jpdic) - zcaralk ) 278 hi(ji,jj,jk) = ( ak23(ji,jj,jk) * zbicarb / zco3 ) * ztmas + 1.e-9 * ztmas1 279 END DO 280 END DO 281 END DO 306 CALL p4z_ph_ini 282 307 ENDIF 283 308 CALL iom_get( numrtr, jpdom_autoglo, 'Silicalim', xksi(:,:) ) -
branches/2013/dev_r3940_CNRS4_IOCRS/NEMOGCM/NEMO/TOP_SRC/PISCES/trcini_pisces.F90
r3757 r4077 122 122 rdenita = 3._wp / 5._wp 123 123 o2ut = 131._wp / 122._wp 124 125 CALL p4z_che ! initialize the chemical constants126 124 127 125 ! Initialization of tracer concentration in case of no restart … … 162 160 xksi(:,:) = 2.e-6 163 161 xksimax(:,:) = xksi(:,:) 164 165 ! Initialization of chemical variables of the carbon cycle166 ! --------------------------------------------------------167 DO jk = 1, jpk168 DO jj = 1, jpj169 DO ji = 1, jpi170 ztmas = tmask(ji,jj,jk)171 ztmas1 = 1. - tmask(ji,jj,jk)172 zcaralk = trn(ji,jj,jk,jptal) - borat(ji,jj,jk) / ( 1. + 1.E-8 / ( rtrn + akb3(ji,jj,jk) ) )173 zco3 = ( zcaralk - trn(ji,jj,jk,jpdic) ) * ztmas + 0.5e-3 * ztmas1174 zbicarb = ( 2. * trn(ji,jj,jk,jpdic) - zcaralk )175 hi(ji,jj,jk) = ( ak23(ji,jj,jk) * zbicarb / zco3 ) * ztmas + 1.e-9 * ztmas1176 END DO177 END DO178 END DO179 162 ! 180 163 END IF
Note: See TracChangeset
for help on using the changeset viewer.