Changeset 11080
- Timestamp:
- 2019-06-06T15:32:24+02:00 (6 years ago)
- Location:
- branches/UKMO/r8395_India_uncoupled/NEMOGCM/NEMO/OPA_SRC
- Files:
-
- 7 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/r8395_India_uncoupled/NEMOGCM/NEMO/OPA_SRC/BDY/bdydta.F90
r10684 r11080 759 759 ENDDO 760 760 761 DO jfld = 1, nb_bdy_fld_sum 762 bf(jfld)%igrd = igrid(jfld) 763 bf(jfld)%ibdy = ibdy(jfld) 764 ENDDO 765 761 766 ! Initialise local boundary data arrays 762 767 ! nn_xxx_dta=0 : allocate space - will be filled from initial conditions later -
branches/UKMO/r8395_India_uncoupled/NEMOGCM/NEMO/OPA_SRC/DYN/dynnxt.F90
r10684 r11080 31 31 USE dynadv ! dynamics: vector invariant versus flux form 32 32 USE dynspg_ts ! surface pressure gradient: split-explicit scheme 33 USE dynspg 33 34 USE domvvl ! variable volume 34 35 USE bdy_oce , ONLY: ln_bdy … … 179 180 vn(:,:,jk) = va(:,:,jk) 180 181 END DO 182 ! limit velocities 183 IF (ln_ulimit) THEN 184 call dyn_limit_velocity (kt) 185 ENDIF 181 186 IF(.NOT.ln_linssh ) THEN 182 187 DO jk = 1, jpkm1 … … 203 208 END DO 204 209 END DO 210 ! limit velocities 211 IF (ln_ulimit) THEN 212 call dyn_limit_velocity (kt) 213 ENDIF 205 214 ! ! ================! 206 215 ELSE ! Variable volume ! … … 250 259 END DO 251 260 END DO 261 ! limit velocities 262 IF (ln_ulimit) THEN 263 call dyn_limit_velocity (kt) 264 ENDIF 252 265 ! 253 266 ELSE ! Asselin filter applied on thickness weighted velocity … … 277 290 END DO 278 291 END DO 292 ! limit velocities 293 IF (ln_ulimit) THEN 294 call dyn_limit_velocity (kt) 295 ENDIF 279 296 e3u_b(:,:,1:jpkm1) = ze3u_f(:,:,1:jpkm1) ! e3u_b <-- filtered scale factor 280 297 e3v_b(:,:,1:jpkm1) = ze3v_f(:,:,1:jpkm1) … … 353 370 END SUBROUTINE dyn_nxt 354 371 372 SUBROUTINE dyn_limit_velocity (kt) 373 ! limits maxming vlaues of un and vn by volume courant number 374 INTEGER, INTENT( in ) :: kt ! ocean time-step index 375 ! 376 INTEGER :: ji, jj, jk ! dummy loop indices 377 REAL(wp) :: zzu,zplim,zmlim,isp,ism,zcn,ze3e1,zzcn,zcnn,idivp,idivm 378 379 ! limit fluxes 380 zcn =cn_ulimit !0.9 ! maximum velocity inverse courant number 381 zcnn = cnn_ulimit !0.54 ! how much to reduce cn by in divergen flow 382 383 DO jk = 1, jpkm1 384 DO jj = 1, jpjm1 385 DO ji = 1, jpim1 386 ! U direction 387 zzu = un(ji,jj,jk) 388 ze3e1 = e3u_n(ji ,jj,jk) * e2u(ji ,jj) 389 ! ips is 1 if flow is positive othersize zero 390 isp = 0.5 * (sign(1.0_wp,zzu) + 1.0_wp ) 391 ism = -0.5 * (sign(1.0_wp,zzu) - 1.0_wp ) 392 ! idev is 1 if divergent flow otherwise zero 393 idivp = -isp * 0.5 * (sign(1.0_wp, un(ji-1,jj,jk)) - 1.0_wp ) 394 idivm = ism * 0.5 * (sign(1.0_wp, un(ji+1,jj,jk)) + 1.0_wp ) 395 zzcn = (idivp+idivm)*(zcnn-1.0_wp)+1.0_wp 396 zzcn = zzcn * zcn 397 zplim = zzcn * (e3t_n(ji ,jj,jk) * e1t(ji ,jj) * e2t(ji ,jj)) / & 398 (2.0*rdt * ze3e1)*umask(ji,jj,jk) 399 zmlim = -zzcn * (e3t_n(ji+1,jj,jk) * e1t(ji+1,jj) * e2t(ji+1,jj)) / & 400 (2.0*rdt * ze3e1)*umask(ji,jj,jk) 401 ! limit currents 402 un(ji,jj,jk) = min ( zzu,zplim) * isp + max (zzu,zmlim) *ism 403 ! V direction 404 zzu = vn(ji,jj,jk) 405 ze3e1 = e3v_n(ji ,jj,jk) * e1v(ji ,jj) 406 isp = 0.5 * (sign(1.0_wp,zzu) + 1.0_wp ) 407 ism = -0.5 * (sign(1.0_wp,zzu) - 1.0_wp ) 408 ! idev is 1 if divergent flow otherwise zero 409 idivp = -isp * 0.5 * (sign(1.0_wp, vn(ji,jj-1,jk)) - 1.0_wp ) 410 idivm = ism * 0.5 * (sign(1.0_wp, vn(ji,jj+1,jk)) + 1.0_wp ) 411 zzcn = (idivp+idivm)*(zcnn-1.0_wp)+1.0_wp 412 zzcn = zzcn * zcn 413 zplim = zzcn * (e3t_n(ji,jj ,jk) * e1t(ji,jj ) * e2t(ji,jj )) / & 414 (2.0*rdt * ze3e1)*vmask(ji,jj,jk) 415 zmlim = -zzcn * (e3t_n(ji,jj+1,jk) * e1t(ji,jj+1) * e2t(ji,jj+1)) / & 416 (2.0*rdt * ze3e1)*vmask(ji,jj,jk) 417 vn(ji,jj,jk) = min ( zzu,zplim) * isp + max (zzu,zmlim) *ism 418 ENDDO 419 ENDDO 420 ENDDO 421 422 END SUBROUTINE dyn_limit_velocity 423 355 424 !!========================================================================= 356 425 END MODULE dynnxt -
branches/UKMO/r8395_India_uncoupled/NEMOGCM/NEMO/OPA_SRC/DYN/dynspg.F90
r10684 r11080 39 39 INTEGER :: nspg = 0 ! type of surface pressure gradient scheme defined from lk_dynspg_... 40 40 41 LOGICAL, PUBLIC :: ln_ulimit 42 REAL(wp), PUBLIC :: cn_ulimit,cnn_ulimit 43 41 44 ! ! Parameter to control the surface pressure gradient scheme 42 45 INTEGER, PARAMETER :: np_TS = 1 ! split-explicit time stepping (Time-Splitting) … … 183 186 NAMELIST/namdyn_spg/ ln_dynspg_exp , ln_dynspg_ts, & 184 187 & ln_bt_fw, ln_bt_av , ln_bt_auto , & 185 & nn_baro , rn_bt_cmax, nn_bt_flt 188 & nn_baro , rn_bt_cmax, nn_bt_flt , & 189 ln_ulimit, cn_ulimit, cnn_ulimit 186 190 !!---------------------------------------------------------------------- 187 191 ! … … 203 207 WRITE(numout,*) ' Explicit free surface ln_dynspg_exp = ', ln_dynspg_exp 204 208 WRITE(numout,*) ' Free surface with time splitting ln_dynspg_ts = ', ln_dynspg_ts 209 WRITE(numout,*) ' Limit velocities ln_ulimit = ', ln_ulimit 210 WRITE(numout,*) ' Limit velocities max inverse Courant number = ', cn_ulimit 211 WRITE(numout,*) ' Limit velocities multiplier for divergant flow = ', cnn_ulimit 205 212 ENDIF 206 213 ! ! Control of surface pressure gradient scheme options -
branches/UKMO/r8395_India_uncoupled/NEMOGCM/NEMO/OPA_SRC/SBC/fldread.F90
r10684 r11080 810 810 CALL iom_get ( num, jpdom_unknown, 'gdepv', dta_read_z(1:ilendta,1:ipj,1:jpk_bdy) ) 811 811 CALL iom_get ( num, jpdom_unknown, 'e3v', dta_read_dz(1:ilendta,1:ipj,1:jpk_bdy) ) 812 CASE DEFAULT ; CALL ctl_stop( 'unrecognised value for igrd in fld_map' ) 812 813 END SELECT 813 814 … … 909 910 WRITE(ibstr,"(I10.10)") map%ptr(ib) 910 911 CALL ctl_warn('fld_bdy_interp: U depths differ between grids at BDY point '//TRIM(ibstr)//' by more than 1%') 911 IF(lwp) WRITE(*,*) 'DEPTHU', zh, sum(e3u_n(zij,zjj,:), mask=umask(zij,zjj,:)==1), sum(umask(zij,zjj,:)), &912 & hu_n(zij,zjj), map%ptr(ib), ib, zij, zjj, narea-1 , &913 & dta_read(map%ptr(ib),1,:)912 !IF(lwp) WRITE(*,*) 'DEPTHU', zh, sum(e3u_n(zij,zjj,:), mask=umask(zij,zjj,:)==1), sum(umask(zij,zjj,:)), & 913 ! & hu_n(zij,zjj), map%ptr(ib), ib, zij, zjj, narea-1 , & 914 ! & dta_read(map%ptr(ib),1,:) 914 915 ENDIF 915 916 CASE(3) -
branches/UKMO/r8395_India_uncoupled/NEMOGCM/NEMO/OPA_SRC/USR/usrdef_istate.F90
r10684 r11080 62 62 DO jj = 1, jpj 63 63 DO ji = 1, jpi 64 pts(ji,jj,jk,jp_tem) = ( ( 16. - 12. * TANH( (pdept(ji,jj,jk) - 400) / 700 ) ) & 65 & * (-TANH( (500. - pdept(ji,jj,jk)) / 150. ) + 1.) / 2. & 66 & + ( 15. * ( 1. - TANH( (pdept(ji,jj,jk)-50.) / 1500.) ) & 67 & - 1.4 * TANH((pdept(ji,jj,jk)-100.) / 100.) & 68 & + 7. * (1500. - pdept(ji,jj,jk) ) / 1500.) & 69 & * (-TANH( (pdept(ji,jj,jk) - 500.) / 150.) + 1.) / 2. ) * ptmask(ji,jj,jk) 70 71 pts(ji,jj,jk,jp_sal) = ( ( 36.25 - 1.13 * TANH( (pdept(ji,jj,jk) - 305) / 460 ) ) & 72 & * (-TANH((500. - pdept(ji,jj,jk)) / 150.) + 1.) / 2 & 73 & + ( 35.55 + 1.25 * (5000. - pdept(ji,jj,jk)) / 5000. & 74 & - 1.62 * TANH( (pdept(ji,jj,jk) - 60. ) / 650. ) & 75 & + 0.2 * TANH( (pdept(ji,jj,jk) - 35. ) / 100. ) & 76 & + 0.2 * TANH( (pdept(ji,jj,jk) - 1000.) / 5000.) ) & 77 & * (-TANH( (pdept(ji,jj,jk) - 500.) / 150.) + 1.) / 2 ) * ptmask(ji,jj,jk) 64 pts(ji,jj,jk,jp_tem) = 20._wp * ptmask(ji,jj,jk) 65 pts(ji,jj,jk,jp_sal) = 36.25_wp * ptmask(ji,jj,jk) 78 66 END DO 79 67 END DO -
branches/UKMO/r8395_India_uncoupled/NEMOGCM/NEMO/OPA_SRC/USR/usrdef_sbc.F90
r7753 r11080 1 1 MODULE usrdef_sbc 2 2 !!====================================================================== 3 !! *** MODULEusrdef_sbc ***4 !! 5 !! === GYREconfiguration ===3 !! *** MODULE usrdef_sbc *** 4 !! 5 !! === WAD_TEST_CASES configuration === 6 6 !! 7 7 !! User defined : surface forcing of a user configuration … … 11 11 12 12 !!---------------------------------------------------------------------- 13 !! usrdef_sbc : user defined surface bounday conditions in GYREcase13 !! usrdef_sbc : user defined surface bounday conditions in WAD_TEST_CASES case 14 14 !!---------------------------------------------------------------------- 15 15 USE oce ! ocean dynamics and tracers … … 21 21 USE lib_mpp ! distribued memory computing library 22 22 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 23 USE lib_fortran ! 23 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 24 24 25 25 IMPLICIT NONE … … 41 41 SUBROUTINE usrdef_sbc_oce( kt ) 42 42 !!--------------------------------------------------------------------- 43 !! *** ROUTINE usr def_sbc ***43 !! *** ROUTINE usr_def_sbc *** 44 44 !! 45 !! ** Purpose : provide at each time-step the GYREsurface boundary45 !! ** Purpose : provide at each time-step the surface boundary 46 46 !! condition, i.e. the momentum, heat and freshwater fluxes. 47 47 !! 48 !! ** Method : a nalytical seasonal cycle for GYRE configuration.48 !! ** Method : all 0 fields, for WAD_TEST_CASES case 49 49 !! CAUTION : never mask the surface stress field ! 50 50 !! 51 !! ** Action : - set t he ocean surface boundary condition, i.e.51 !! ** Action : - set to ZERO all the ocean surface boundary condition, i.e. 52 52 !! utau, vtau, taum, wndm, qns, qsr, emp, sfx 53 53 !! 54 !! Reference : Hazeleger, W., and S. Drijfhout, JPO, 30, 677-695, 2000.55 54 !!---------------------------------------------------------------------- 56 55 INTEGER, INTENT(in) :: kt ! ocean time step 57 !!58 INTEGER :: ji, jj ! dummy loop indices59 INTEGER :: zyear0 ! initial year60 INTEGER :: zmonth0 ! initial month61 INTEGER :: zday0 ! initial day62 INTEGER :: zday_year0 ! initial day since january 1st63 REAL(wp) :: ztau , ztau_sais ! wind intensity and of the seasonal cycle64 REAL(wp) :: ztime ! time in hour65 REAL(wp) :: ztimemax , ztimemin ! 21th June, and 21th decem. if date0 = 1st january66 REAL(wp) :: ztimemax1, ztimemin1 ! 21th June, and 21th decem. if date0 = 1st january67 REAL(wp) :: ztimemax2, ztimemin2 ! 21th June, and 21th decem. if date0 = 1st january68 REAL(wp) :: ztaun ! intensity69 REAL(wp) :: zemp_s, zemp_n, zemp_sais, ztstar70 REAL(wp) :: zcos_sais1, zcos_sais2, ztrp, zconv, t_star71 REAL(wp) :: zsumemp, zsurf72 REAL(wp) :: zrhoa = 1.22 ! Air density kg/m373 REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient74 REAL(wp) :: ztx, zty, zmod, zcoef ! temporary variables75 REAL(wp) :: zyydd ! number of days in one year76 56 !!--------------------------------------------------------------------- 77 zyydd = REAL(nyear_len(1),wp) 78 79 ! ---------------------------- ! 80 ! heat and freshwater fluxes ! 81 ! ---------------------------- ! 82 !same temperature, E-P as in HAZELEGER 2000 83 84 zyear0 = ndate0 / 10000 ! initial year 85 zmonth0 = ( ndate0 - zyear0 * 10000 ) / 100 ! initial month 86 zday0 = ndate0 - zyear0 * 10000 - zmonth0 * 100 ! initial day betwen 1 and 30 87 zday_year0 = ( zmonth0 - 1 ) * 30.+zday0 ! initial day betwen 1 and 360 88 89 ! current day (in hours) since january the 1st of the current year 90 ztime = REAL( kt ) * rdt / (rmmss * rhhmm) & ! total incrementation (in hours) 91 & - (nyear - 1) * rjjhh * zyydd ! minus years since beginning of experiment (in hours) 92 93 ztimemax1 = ((5.*30.)+21.)* 24. ! 21th june at 24h in hours 94 ztimemin1 = ztimemax1 + rjjhh * zyydd / 2 ! 21th december in hours 95 ztimemax2 = ((6.*30.)+21.)* 24. ! 21th july at 24h in hours 96 ztimemin2 = ztimemax2 - rjjhh * zyydd / 2 ! 21th january in hours 97 ! ! NB: rjjhh * zyydd / 4 = one seasonal cycle in hours 98 99 ! amplitudes 100 zemp_S = 0.7 ! intensity of COS in the South 101 zemp_N = 0.8 ! intensity of COS in the North 102 zemp_sais = 0.1 103 zTstar = 28.3 ! intemsity from 28.3 a -5 deg 104 105 ! 1/2 period between 21th June and 21th December and between 21th July and 21th January 106 zcos_sais1 = COS( (ztime - ztimemax1) / (ztimemin1 - ztimemax1) * rpi ) 107 zcos_sais2 = COS( (ztime - ztimemax2) / (ztimemax2 - ztimemin2) * rpi ) 108 109 ztrp= - 40.e0 ! retroaction term on heat fluxes (W/m2/K) 110 zconv = 3.16e-5 ! convertion factor: 1 m/yr => 3.16e-5 mm/s 111 DO jj = 1, jpj 112 DO ji = 1, jpi 113 ! domain from 15 deg to 50 deg between 27 and 28 degC at 15N, -3 114 ! and 13 degC at 50N 53.5 + or - 11 = 1/4 period : 115 ! 64.5 in summer, 42.5 in winter 116 t_star = zTstar * ( 1. + 1. / 50. * zcos_sais2 ) & 117 & * COS( rpi * (gphit(ji,jj) - 5.) & 118 & / ( 53.5 * ( 1 + 11 / 53.5 * zcos_sais2 ) * 2.) ) 119 ! 23.5 deg : tropics 120 qsr (ji,jj) = 230 * COS( 3.1415 * ( gphit(ji,jj) - 23.5 * zcos_sais1 ) / ( 0.9 * 180 ) ) 121 qns (ji,jj) = ztrp * ( tsb(ji,jj,1,jp_tem) - t_star ) - qsr(ji,jj) 122 IF( gphit(ji,jj) >= 14.845 .AND. 37.2 >= gphit(ji,jj) ) THEN ! zero at 37.8 deg, max at 24.6 deg 123 emp (ji,jj) = zemp_S * zconv & 124 & * SIN( rpi / 2 * (gphit(ji,jj) - 37.2) / (24.6 - 37.2) ) & 125 & * ( 1 - zemp_sais / zemp_S * zcos_sais1) 126 ELSE 127 emp (ji,jj) = - zemp_N * zconv & 128 & * SIN( rpi / 2 * (gphit(ji,jj) - 37.2) / (46.8 - 37.2) ) & 129 & * ( 1 - zemp_sais / zemp_N * zcos_sais1 ) 130 ENDIF 131 END DO 132 END DO 133 134 zsumemp = GLOB_SUM( emp (:,:) ) 135 zsurf = GLOB_SUM( tmask(:,:,1) ) 136 zsumemp = zsumemp / zsurf ! Default GYRE configuration 137 138 ! freshwater (mass flux) and update of qns with heat content of emp 139 emp (:,:) = emp(:,:) - zsumemp * tmask(:,:,1) ! freshwater flux (=0 in domain average) 140 sfx (:,:) = 0.0_wp ! no salt flux 141 qns (:,:) = qns(:,:) - emp(:,:) * sst_m(:,:) * rcp ! evap and precip are at SST 142 143 144 ! ---------------------------- ! 145 ! momentum fluxes ! 146 ! ---------------------------- ! 147 ! same wind as in Wico 148 !test date0 : ndate0 = 010203 149 zyear0 = ndate0 / 10000 150 zmonth0 = ( ndate0 - zyear0 * 10000 ) / 100 151 zday0 = ndate0 - zyear0 * 10000 - zmonth0 * 100 152 !Calculates nday_year, day since january 1st 153 zday_year0 = (zmonth0-1)*30.+zday0 154 155 !accumulates days of previous months of this year 156 ! day (in hours) since january the 1st 157 ztime = FLOAT( kt ) * rdt / (rmmss * rhhmm) & ! incrementation in hour 158 & - (nyear - 1) * rjjhh * zyydd ! - nber of hours the precedent years 159 ztimemax = ((5.*30.)+21.)* 24. ! 21th june in hours 160 ztimemin = ztimemax + rjjhh * zyydd / 2 ! 21th december in hours 161 ! ! NB: rjjhh * zyydd / 4 = 1 seasonal cycle in hours 162 163 ! mean intensity at 0.105 ; srqt(2) because projected with 45deg angle 164 ztau = 0.105 / SQRT( 2. ) 165 ! seasonal oscillation intensity 166 ztau_sais = 0.015 167 ztaun = ztau - ztau_sais * COS( (ztime - ztimemax) / (ztimemin - ztimemax) * rpi ) 168 DO jj = 1, jpj 169 DO ji = 1, jpi 170 ! domain from 15deg to 50deg and 1/2 period along 14deg 171 ! so 5/4 of half period with seasonal cycle 172 utau(ji,jj) = - ztaun * SIN( rpi * (gphiu(ji,jj) - 15.) / (29.-15.) ) 173 vtau(ji,jj) = ztaun * SIN( rpi * (gphiv(ji,jj) - 15.) / (29.-15.) ) 174 END DO 175 END DO 176 177 ! module of wind stress and wind speed at T-point 178 zcoef = 1. / ( zrhoa * zcdrag ) 179 DO jj = 2, jpjm1 180 DO ji = fs_2, fs_jpim1 ! vect. opt. 181 ztx = utau(ji-1,jj ) + utau(ji,jj) 182 zty = vtau(ji ,jj-1) + vtau(ji,jj) 183 zmod = 0.5 * SQRT( ztx * ztx + zty * zty ) 184 taum(ji,jj) = zmod 185 wndm(ji,jj) = SQRT( zmod * zcoef ) 186 END DO 187 END DO 188 CALL lbc_lnk( taum(:,:), 'T', 1. ) ; CALL lbc_lnk( wndm(:,:), 'T', 1. ) 189 190 ! ---------------------------------- ! 191 ! control print at first time-step ! 192 ! ---------------------------------- ! 193 IF( kt == nit000 .AND. lwp ) THEN 194 WRITE(numout,*) 195 WRITE(numout,*)'usrdef_sbc_oce : analytical surface fluxes for GYRE configuration' 196 WRITE(numout,*)'~~~~~~~~~~~ ' 197 WRITE(numout,*)' nyear = ', nyear 198 WRITE(numout,*)' nmonth = ', nmonth 199 WRITE(numout,*)' nday = ', nday 200 WRITE(numout,*)' nday_year = ', nday_year 201 WRITE(numout,*)' ztime = ', ztime 202 WRITE(numout,*)' ztimemax = ', ztimemax 203 WRITE(numout,*)' ztimemin = ', ztimemin 204 WRITE(numout,*)' ztimemax1 = ', ztimemax1 205 WRITE(numout,*)' ztimemin1 = ', ztimemin1 206 WRITE(numout,*)' ztimemax2 = ', ztimemax2 207 WRITE(numout,*)' ztimemin2 = ', ztimemin2 208 WRITE(numout,*)' zyear0 = ', zyear0 209 WRITE(numout,*)' zmonth0 = ', zmonth0 210 WRITE(numout,*)' zday0 = ', zday0 211 WRITE(numout,*)' zday_year0 = ', zday_year0 212 WRITE(numout,*)' zyydd = ', zyydd 213 WRITE(numout,*)' zemp_S = ', zemp_S 214 WRITE(numout,*)' zemp_N = ', zemp_N 215 WRITE(numout,*)' zemp_sais = ', zemp_sais 216 WRITE(numout,*)' zTstar = ', zTstar 217 WRITE(numout,*)' zsumemp = ', zsumemp 218 WRITE(numout,*)' zsurf = ', zsurf 219 WRITE(numout,*)' ztrp = ', ztrp 220 WRITE(numout,*)' zconv = ', zconv 221 WRITE(numout,*)' ndastp = ', ndastp 222 WRITE(numout,*)' adatrj = ', adatrj 57 ! 58 IF( kt == nit000 ) THEN 59 ! 60 IF(lwp) WRITE(numout,*)' usr_sbc : WAD_TEST_CASES case: NO surface forcing' 61 IF(lwp) WRITE(numout,*)' ~~~~~~~~~~~ utau = vtau = taum = wndm = qns = qsr = emp = sfx = 0' 62 ! 63 utau(:,:) = 0._wp 64 vtau(:,:) = 0._wp 65 taum(:,:) = 0._wp 66 wndm(:,:) = 0._wp 67 ! 68 emp (:,:) = 0._wp 69 sfx (:,:) = 0._wp 70 qns (:,:) = 0._wp 71 qsr (:,:) = 0._wp 72 ! 223 73 ENDIF 224 74 ! -
branches/UKMO/r8395_India_uncoupled/NEMOGCM/NEMO/OPA_SRC/stpctl.F90
r10684 r11080 52 52 INTEGER :: ji, jj, jk ! dummy loop indices 53 53 INTEGER :: ii, ij, ik ! local integers 54 REAL(wp) :: zumax, zsmin, zssh2, zsshmax ! local scalars54 REAL(wp) :: velmax2, zsmin, zssh2, zsshmax ! local scalars 55 55 INTEGER, DIMENSION(3) :: ilocu ! 56 56 INTEGER, DIMENSION(2) :: ilocs ! … … 70 70 ! !* Test maximum of velocity (zonal only) 71 71 ! ! ------------------------ 72 !! zumax= MAXVAL( ABS( un(:,:,:) ) ) ! slower than the following loop on NEC SX573 zumax= 0.e072 !! velmax2 = MAXVAL( ABS( un(:,:,:) ) ) ! slower than the following loop on NEC SX5 73 velmax2 = 0.e0 74 74 DO jk = 1, jpk 75 75 DO jj = 1, jpj 76 76 DO ji = 1, jpi 77 zumax = MAX(zumax,ABS(un(ji,jj,jk)))77 velmax2 = MAX( velmax2,un(ji,jj,jk)**2 + vn(ji,jj,jk)**2 ) 78 78 END DO 79 79 END DO 80 80 END DO 81 IF( lk_mpp ) CALL mpp_max( zumax )! max over the global domain81 IF( lk_mpp ) CALL mpp_max( velmax2 ) ! max over the global domain 82 82 ! 83 IF( MOD( kt, nwrite ) == 1 .AND. lwp ) WRITE(numout,*) ' ==>> time-step= ',kt,' abs(U) max: ', zumax83 IF( MOD( kt, nwrite ) == 1 .AND. lwp ) WRITE(numout,*) ' ==>> time-step= ',kt,' 3d speed max: ', velmax2 84 84 ! 85 IF( zumax > 20.e0) THEN85 IF( velmax2 > 20.e0**2 ) THEN 86 86 IF( lk_mpp ) THEN 87 CALL mpp_maxloc( ABS(un),umask,zumax,ii,ij,ik)87 CALL mpp_maxloc( un(:,:,:)**2+vn(:,:,:)**2,umask,velmax2,ii,ij,ik) 88 88 ELSE 89 ilocu = MAXLOC( ABS( un(:,:,:) ))89 ilocu = MAXLOC( un(:,:,:)**2 + vn(:,:,:)**2 ) 90 90 ii = ilocu(1) + nimpp - 1 91 91 ij = ilocu(2) + njmpp - 1 … … 94 94 IF(lwp) THEN 95 95 WRITE(numout,cform_err) 96 WRITE(numout,*) ' stpctl: the zonal velocityis larger than 20 m/s'96 WRITE(numout,*) ' stpctl: the speed is larger than 20 m/s' 97 97 WRITE(numout,*) ' ====== ' 98 WRITE(numout,9400) kt, zumax, ii, ij, ik98 WRITE(numout,9400) kt, velmax2, ii, ij, ik 99 99 WRITE(numout,*) 100 100 WRITE(numout,*) ' output of last fields in numwso' … … 102 102 kindic = -3 103 103 ENDIF 104 9400 FORMAT (' kt=',i6,' max abs( U): ',1pg11.4,', i j k: ',3i5)104 9400 FORMAT (' kt=',i6,' max abs(vel)**2: ',1pg11.4,', i j k: ',3i5) 105 105 ! 106 106 ! !* Test minimum of salinity … … 180 180 zssh2 = glob_sum( sshn(:,:) * sshn(:,:) ) 181 181 ! 182 IF(lwp) WRITE(numsol,9700) kt, zssh2, zumax, zsmin ! ssh statistics182 IF(lwp) WRITE(numsol,9700) kt, zssh2, velmax2, zsmin ! ssh statistics 183 183 ! 184 9700 FORMAT(' it :', i8, ' ssh2: ', d23.16, ' Umax: ',d23.16,'Smin: ',d23.16)184 9700 FORMAT(' it :', i8, ' ssh2: ', d23.16, ' velmax2: ',d23.16,' SSSmin: ',d23.16) 185 185 ! 186 186 END SUBROUTINE stp_ctl
Note: See TracChangeset
for help on using the changeset viewer.