Changeset 15532
- Timestamp:
- 2021-11-24T12:47:32+01:00 (2 years ago)
- Location:
- NEMO/branches/2021/dev_r14318_RK3_stage1/src
- Files:
-
- 7 edited
-
OCE/nemogcm.F90 (modified) (1 diff)
-
OCE/stprk3_stg.F90 (modified) (2 diffs)
-
OFF/dtadyn.F90 (modified) (4 diffs)
-
TOP/PISCES/P4Z/p4zche.F90 (modified) (8 diffs)
-
TOP/PISCES/P4Z/p4zflx.F90 (modified) (3 diffs)
-
TOP/PISCES/P4Z/p4zlys.F90 (modified) (6 diffs)
-
TOP/TRP/trcsbc.F90 (modified) (9 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2021/dev_r14318_RK3_stage1/src/OCE/nemogcm.F90
r14547 r15532 489 489 #if defined key_top 490 490 ! ! Passive tracers 491 # if defined key_RK3 492 CALL trc_init( Nbb, Nbb, Naa ) 493 # else 491 494 CALL trc_init( Nbb, Nnn, Naa ) 495 # endif 492 496 #endif 493 497 IF( l_ldfslp ) CALL ldf_slp_init ! slope of lateral mixing -
NEMO/branches/2021/dev_r14318_RK3_stage1/src/OCE/stprk3_stg.F90
r15514 r15532 185 185 !===>>>>>> Modify dyn_hpg & dyn_hpg_... routines : rhd computed in dyn_hpg and pass in argument to dyn_hpg_... 186 186 187 CALL eos ( ts, Kmm, rhd ) ! Kmm in situ density anomaly for hpg computation 187 !!st IF( kstg == 3 ) THEN 188 ! CALL eos ( ts(:,:,:,:,Kmm), rhd, rhop, gdept_0 ) ! now in situ density for hpg computation 189 ! ELSE 190 CALL eos ( ts, Kmm, rhd ) ! Kmm in situ density anomaly for hpg computation 191 ! ENDIF 188 192 189 193 !!gm end … … 388 392 vv(:,:,jk,Kaa) = vv(:,:,jk,Kaa) + zvb(:,:)*vmask(:,:,jk) 389 393 END DO 394 !!st pourquoi ne pas mettre A2D(0) ici ? 390 395 391 396 -
NEMO/branches/2021/dev_r14318_RK3_stage1/src/OFF/dtadyn.F90
r14310 r15532 176 176 ENDIF 177 177 ! 178 CALL eos ( ts(:,:,:,:,Kmm), rhd, rhop, gdept_0(:,:,:) ) ! In any case, we need rhop178 CALL eos ( ts(:,:,:,:,Kmm), rhd, gdept_0(:,:,:) ) ! In any case, we need rhd 179 179 CALL eos_rab( ts(:,:,:,:,Kmm), rab_n, Kmm ) ! now local thermal/haline expension ratio at T-points 180 180 CALL bn2 ( ts(:,:,:,:,Kmm), rab_n, rn2, Kmm ) ! before Brunt-Vaisala frequency need for zdfmxl … … 193 193 ! 194 194 ! 195 CALL eos( ts(:,:,:,:,Kmm), rhd, rhop, gdept_0(:,:,:) ) ! In any case, we need rhop195 CALL eos( ts(:,:,:,:,Kmm), rhd, gdept_0(:,:,:) ) ! In any case, we need rhd 196 196 ! 197 197 IF(sn_cfctl%l_prtctl) THEN ! print control … … 666 666 ! 667 667 IF( l_ldfslp .AND. .NOT.lk_c1d ) THEN ! Computes slopes (here avt is used as workspace) 668 CALL eos ( pts, rhd, rhop,gdept_0(:,:,:) )668 CALL eos ( pts, rhd, gdept_0(:,:,:) ) 669 669 CALL eos_rab( pts, rab_n, Kmm ) ! now local thermal/haline expension ratio at T-points 670 670 CALL bn2 ( pts, rab_n, rn2, Kmm ) ! now Brunt-Vaisala … … 724 724 ts(:,:,:,jp_sal,Kmm) = sf_dyn(jf_sal)%fnow(:,:,:) * tmask(:,:,:) ! salinity 725 725 ! 726 CALL eos ( ts(:,:,:,:,Kmm), rhd, rhop, gdept_0(:,:,:) ) ! In any case, we need rhop726 CALL eos ( ts(:,:,:,:,Kmm), rhd, gdept_0(:,:,:) ) ! In any case, we need rhd 727 727 728 728 IF(sn_cfctl%l_prtctl) THEN ! print control -
NEMO/branches/2021/dev_r14318_RK3_stage1/src/TOP/PISCES/P4Z/p4zche.F90
r14086 r15532 12 12 !! ! 2011-02 (J. Simeon, J.Orr ) update O2 solubility constants 13 13 !! 3.6 ! 2016-03 (O. Aumont) Change chemistry to MOCSY standards 14 !! 4.2 ! 2020 (J. ORR ) rhop is replaced by "in situ density" rhd 14 15 !!---------------------------------------------------------------------- 15 16 !! p4z_che : Sea water chemistry computed following OCMIP protocol … … 195 196 ! ! LN(K0) OF SOLUBILITY OF CO2 (EQ. 12, WEISS, 1980) 196 197 ! ! AND FOR THE ATMOSPHERE FOR NON IDEAL GAS 197 zcek1 = 9345.17/ztkel - 60.2409 + 23.3585 * LOG(zt) + zsal*(0.023517 - 0.00023656*ztkel & 198 & + 0.0047036e-4*ztkel**2) 199 chemc(ji,jj,1) = EXP( zcek1 ) * 1E-6 * rhop(ji,jj,1) / 1000. ! mol/(L atm) 198 ! J. ORR: The previous code has been modified. It computed CO2 solubility in mol/(kg*atm), then converted that to mol/(L*atm). 199 ! But Weiss (1974) provides sets of coefficients for each of those 2 units. 200 ! Thus I have changed the code to use the coefficients for mol*L/atm. 201 ! Hence I've eliminated using the conversion (which used the variable rhop) 202 ! OLD - Coefficients for CO2 soulbility in mol/(kg*atm) (Weiss,1974, Table 1, column 2) 203 !zcek1 = 9345.17/ztkel - 60.2409 + 23.3585 * LOG(zt) + zsal*(0.023517 - 0.00023656*ztkel & 204 !& + 0.0047036e-4*ztkel**2) 205 ! NEW - Coefficients for CO2 soulbility in mol/(L*atm) (Weiss, 1974, Table 1, column 1) 206 zcek1 = 9050.69/ztkel - 58.0931 + 22.2940 * LOG(zt) + zsal*(0.027766 - 0.00025888*ztkel & 207 & + 0.0050578e-4*ztkel**2) 208 ! 209 ! OLD: chemc(ji,jj,1) = EXP( zcek1 ) * 1E-6 * rhop(ji,jj,1) / 1000. ! mol/(L atm) 210 ! The units indicated in the above line are wrong. They are actually "mol/(L*uatm)" 211 ! NEW: 212 chemc(ji,jj,1) = EXP( zcek1 ) * 1E-6 ! mol/(L * uatm) 200 213 chemc(ji,jj,2) = -1636.75 + 12.0408*ztkel - 0.0327957*ztkel**2 + 0.0000316528*ztkel**3 201 214 chemc(ji,jj,3) = 57.7 - 0.118*ztkel … … 445 458 REAL(wp) :: zca1, zba1 446 459 REAL(wp) :: zd, zsqrtd, zhmin 447 REAL(wp) :: za2, za1, za0 460 REAL(wp) :: za2, za1, za0, zrhd 448 461 REAL(wp) :: p_dictot, p_bortot, p_alkcb 449 462 !!--------------------------------------------------------------------- … … 452 465 ! 453 466 DO_3D( 1, 1, 1, 1, 1, jpk ) 454 p_alkcb = tr(ji,jj,jk,jptal,Kbb) * 1000. / (rhop(ji,jj,jk) + rtrn) 455 p_dictot = tr(ji,jj,jk,jpdic,Kbb) * 1000. / (rhop(ji,jj,jk) + rtrn) 467 zrhd = 1._wp / ( rhd(ji,jj,jk) + 1. ) 468 p_alkcb = tr(ji,jj,jk,jptal,Kbb) * zrhd 469 p_dictot = tr(ji,jj,jk,jpdic,Kbb) * zrhd 456 470 p_bortot = borat(ji,jj,jk) 457 471 IF (p_alkcb <= 0.) THEN … … 502 516 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(OUT) :: p_alknw_sup 503 517 INTEGER, INTENT(in) :: Kbb ! time level indices 504 505 p_alknw_inf(:,:,:) = -tr(:,:,:,jppo4,Kbb) * 1000. / (rhop(:,:,:) + rtrn) - sulfat(:,:,:) & 506 & - fluorid(:,:,:) 507 p_alknw_sup(:,:,:) = (2. * tr(:,:,:,jpdic,Kbb) + 2. * tr(:,:,:,jppo4,Kbb) + tr(:,:,:,jpsil,Kbb) ) & 508 & * 1000. / (rhop(:,:,:) + rtrn) + borat(:,:,:) 518 INTEGER :: ji, jj, jk 519 REAL(wp) :: zrhd 520 521 DO_3D( 1, 1, 1, 1, 1, jpk ) 522 zrhd = 1._wp / ( rhd(ji,jj,jk) + 1. ) 523 p_alknw_inf(ji,jj,jk) = -tr(ji,jj,jk,jppo4,Kbb) * zrhd - sulfat(ji,jj,jk) & 524 & - fluorid(ji,jj,jk) 525 p_alknw_sup(ji,jj,jk) = (2. * tr(ji,jj,jk,jpdic,Kbb) + 2. * tr(ji,jj,jk,jppo4,Kbb) + tr(ji,jj,jk,jpsil,Kbb) ) & 526 & * zrhd + borat(ji,jj,jk) 527 END_3D 509 528 510 529 END SUBROUTINE anw_infsup … … 536 555 REAL(wp) :: znumer_flu, zdnumer_flu, zdenom_flu, zalk_flu, zdalk_flu 537 556 REAL(wp) :: zalk_wat, zdalk_wat 538 REAL(wp) :: z fact, p_alktot, zdic, zbot, zpt, zst, zft, zsit557 REAL(wp) :: zrhd, p_alktot, zdic, zbot, zpt, zst, zft, zsit 539 558 LOGICAL :: l_exitnow 540 559 REAL(wp), PARAMETER :: pz_exp_threshold = 1.0 … … 551 570 DO_3D( 1, 1, 1, 1, 1, jpk ) 552 571 IF (rmask(ji,jj,jk) == 1.) THEN 553 p_alktot = tr(ji,jj,jk,jptal,Kbb) * 1000. / (rhop(ji,jj,jk) + rtrn) 572 zrhd = 1._wp / ( rhd(ji,jj,jk) + 1. ) 573 p_alktot = tr(ji,jj,jk,jptal,Kbb) * zrhd 554 574 aphscale = 1. + sulfat(ji,jj,jk)/aks3(ji,jj,jk) 555 575 zh_ini = p_hini(ji,jj,jk) … … 580 600 DO_3D( 1, 1, 1, 1, 1, jpk ) 581 601 IF (rmask(ji,jj,jk) == 1.) THEN 582 z fact = rhop(ji,jj,jk) / 1000. + rtrn583 p_alktot = tr(ji,jj,jk,jptal,Kbb) / zfact584 zdic = tr(ji,jj,jk,jpdic,Kbb) / zfact602 zrhd = 1._wp / ( rhd(ji,jj,jk) + 1. ) 603 p_alktot = tr(ji,jj,jk,jptal,Kbb) * zrhd 604 zdic = tr(ji,jj,jk,jpdic,Kbb) * zrhd 585 605 zbot = borat(ji,jj,jk) 586 zpt = tr(ji,jj,jk,jppo4,Kbb) / zfact* po4r587 zsit = tr(ji,jj,jk,jpsil,Kbb) / zfact606 zpt = tr(ji,jj,jk,jppo4,Kbb) * zrhd * po4r 607 zsit = tr(ji,jj,jk,jpsil,Kbb) * zrhd 588 608 zst = sulfat (ji,jj,jk) 589 609 zft = fluorid(ji,jj,jk) -
NEMO/branches/2021/dev_r14318_RK3_stage1/src/TOP/PISCES/P4Z/p4zflx.F90
r13295 r15532 10 10 !! 2.0 ! 2007-12 (C. Ethe, G. Madec) F90 11 11 !! ! 2011-02 (J. Simeon, J. Orr) Include total atm P correction 12 !! 4.2 ! 2020 (J. ORR ) rhop is replaced by "in situ density" rhd 12 13 !!---------------------------------------------------------------------- 13 14 !! p4z_flx : CALCULATES GAS EXCHANGE AND CHEMISTRY AT SEA SURFACE … … 78 79 INTEGER :: ji, jj, jm, iind, iindm1 79 80 REAL(wp) :: ztc, ztc2, ztc3, ztc4, zws, zkgwan 80 REAL(wp) :: zfld, zflu, zfld16, zflu16, z fact81 REAL(wp) :: zfld, zflu, zfld16, zflu16, zrhd 81 82 REAL(wp) :: zvapsw, zsal, zfco2, zxc2, xCO2approx, ztkel, zfugcoeff 82 83 REAL(wp) :: zph, zdic, zsch_o2, zsch_co2 … … 112 113 DO_2D( 1, 1, 1, 1 ) 113 114 ! DUMMY VARIABLES FOR DIC, H+, AND BORATE 114 z fact = rhop(ji,jj,1) / 1000. + rtrn115 zrhd = rhd(ji,jj,1) + 1._wp 115 116 zdic = tr(ji,jj,1,jpdic,Kbb) 116 zph = MAX( hi(ji,jj,1), 1.e-10 ) / zfact117 zph = MAX( hi(ji,jj,1), 1.e-10 ) / ( zrhd + rtrn ) 117 118 ! CALCULATE [H2CO3] 118 119 zh2co3(ji,jj) = zdic/(1. + ak13(ji,jj,1)/zph + ak13(ji,jj,1)*ak23(ji,jj,1)/zph**2) -
NEMO/branches/2021/dev_r14318_RK3_stage1/src/TOP/PISCES/P4Z/p4zlys.F90
r13295 r15532 12 12 !! 3.4 ! 2011-06 (O. Aumont, C. Ethe) Improvment of calcite dissolution 13 13 !! 3.6 ! 2015-05 (O. Aumont) PISCES quota 14 !! 4.2 ! 2020 (J. ORR ) rhop is replaced by "in situ density" rhd 14 15 !!---------------------------------------------------------------------- 15 16 !! p4z_lys : Compute the CaCO3 dissolution … … 33 34 34 35 INTEGER :: rmtss ! number of seconds per month 35 REAL(wp) :: calcon = 1.03E-2 ! mean calcite concentration [Ca2+] in sea water [mole/kg solution] 36 36 37 !! * Module variables 38 REAL(wp) :: calcon = 1.03E-2 !: mean calcite concentration [Ca2+] in sea water [mole/kg solution] 39 ! J. ORR: Made consistent with mocsy's choice based on literature review from Munhoven 40 ! REAL(wp) :: calcon = 1.0287E-2 !: mean calcite concentration [Ca2+] in sea water [mole/kg solution] 41 37 42 !! * Substitutions 38 43 # include "do_loop_substitute.h90" … … 59 64 ! 60 65 INTEGER :: ji, jj, jk, jn 61 REAL(wp) :: zdispot, z fact, zcalcon66 REAL(wp) :: zdispot, zrhd, zcalcon 62 67 REAL(wp) :: zomegaca, zexcess, zexcess0 63 68 CHARACTER (len=25) :: charout … … 67 72 IF( ln_timing ) CALL timing_start('p4z_lys') 68 73 ! 69 zhinit (:,:,:) = hi(:,:,:) * 1000. / ( rhop(:,:,:) + rtrn ) 74 75 zhinit (:,:,:) = hi(:,:,:) / ( rhd(:,:,:) + 1._wp ) 70 76 ! 71 77 ! ------------------------------------------- … … 78 84 zco3(ji,jj,jk) = tr(ji,jj,jk,jpdic,Kbb) * ak13(ji,jj,jk) * ak23(ji,jj,jk) / (zhi(ji,jj,jk)**2 & 79 85 & + ak13(ji,jj,jk) * zhi(ji,jj,jk) + ak13(ji,jj,jk) * ak23(ji,jj,jk) + rtrn ) 80 hi (ji,jj,jk) = zhi(ji,jj,jk) * rhop(ji,jj,jk) / 1000.86 hi (ji,jj,jk) = zhi(ji,jj,jk) * ( rhd(ji,jj,jk) + 1._wp ) 81 87 END_3D 82 88 … … 90 96 91 97 ! DEVIATION OF [CO3--] FROM SATURATION VALUE 92 ! Salinity dependance in zomegaca and divide by rh op/1000to have good units98 ! Salinity dependance in zomegaca and divide by rhd to have good units 93 99 zcalcon = calcon * ( salinprac(ji,jj,jk) / 35._wp ) 94 z fact = rhop(ji,jj,jk) / 1000._wp95 zomegaca = ( zcalcon * zco3(ji,jj,jk) ) / ( aksp(ji,jj,jk) * z fact+ rtrn )96 zco3sat(ji,jj,jk) = aksp(ji,jj,jk) * z fact/ ( zcalcon + rtrn )100 zrhd = rhd(ji,jj,jk) + 1._wp 101 zomegaca = ( zcalcon * zco3(ji,jj,jk) ) / ( aksp(ji,jj,jk) * zrhd + rtrn ) 102 zco3sat(ji,jj,jk) = aksp(ji,jj,jk) * zrhd / ( zcalcon + rtrn ) 97 103 98 104 ! SET DEGREE OF UNDER-/SUPERSATURATION -
NEMO/branches/2021/dev_r14318_RK3_stage1/src/TOP/TRP/trcsbc.F90
r15380 r15532 118 118 zsfx(:,:) = emp(:,:) 119 119 ENDIF 120 !!st DO LOOPs used to look like DO_2D( 0, 0, 0, 1 ) I changed it because it does not work in debug mode... 120 121 121 122 ! 0. initialization … … 125 126 ! 126 127 DO jn = 1, jptra 127 DO_2D( 0, 0, 0, 1)128 DO_2D( 0, 0, 0, 0 ) 128 129 sbc_trc(ji,jj,jn) = zsfx(ji,jj) * r1_rho0 * ptr(ji,jj,1,jn,Kmm) 129 130 END_2D … … 133 134 ! 134 135 DO jn = 1, jptra 135 DO_2D( 0, 0, 0, 1)136 DO_2D( 0, 0, 0, 0 ) 136 137 sbc_trc(ji,jj,jn) = ( zsfx(ji,jj) + fmmflx(ji,jj) ) * r1_rho0 * ptr(ji,jj,1,jn,Kmm) 137 138 END_2D … … 141 142 ! 142 143 DO jn = 1, jptra 143 DO_2D( 0, 0, 0, 1)144 DO_2D( 0, 0, 0, 0 ) 144 145 ! tracer flux at the ice/ocean interface (tracer/m2/s) 145 146 zftra = - trc_i(ji,jj,jn) * fmmflx(ji,jj) ! uptake of tracer in the sea ice … … 164 165 IF( l_trdtrc ) ztrtrd(:,:,:) = ptr(:,:,:,jn,Krhs) ! save trends 165 166 ! 166 DO_2D( 0, 0, 0, 1)167 DO_2D( 0, 0, 0, 0 ) 167 168 #if defined key_RK3 168 169 zse3t = 1._wp / e3t(ji,jj,1,Kmm) … … 257 258 !!not sure about trc_i case... (1) 258 259 DO jn = 1, jptra 259 DO_2D( 0, 0, 0, 1) !!st WHY 1 : exterior here ?260 DO_2D( 0, 0, 0, 0 ) !!st WHY 1 : exterior here ? 260 261 z1_rho0_e3t = r1_rho0 / e3t(ji,jj,1,Kmm) 261 262 ptr(ji,jj,1,jn,Krhs) = ptr(ji,jj,1,jn,Krhs) - emp(ji,jj) * ptr(ji,jj,1,jn,Kmm) * z1_rho0_e3t … … 281 282 ! 282 283 DO jn = 1, jptra 283 DO_2D( 0, 0, 0, 1)284 DO_2D( 0, 0, 0, 0 ) 284 285 z1_rho0_e3t = r1_rho0 / e3t(ji,jj,1,Kmm) 285 286 ptr(ji,jj,1,jn,Krhs) = ptr(ji,jj,1,jn,Krhs) + emp(ji,jj) * r1_rho0 * ptr(ji,jj,1,jn,Kmm) … … 299 300 ! 300 301 DO jn = 1, jptra 301 DO_2D( 0, 0, 0, 1)302 DO_2D( 0, 0, 0, 0 ) 302 303 z1_rho0_e3t = r1_rho0 / e3t(ji,jj,1,Kmm) 303 304 ! tracer flux at the ice/ocean interface (tracer/m2/s) … … 335 336 ! 336 337 DO jn = 1, jptra 337 DO_2D( 0, 0, 0, 1)338 DO_2D( 0, 0, 0, 0 ) 338 339 ! tracer flux at the ice/ocean interface (tracer/m2/s) 339 340 zftra = - trc_i(ji,jj,jn) * fmmflx(ji,jj) ! uptake of tracer in the sea ice
Note: See TracChangeset
for help on using the changeset viewer.
