Changeset 13935 for NEMO/branches
- Timestamp:
- 2020-12-01T13:21:25+01:00 (3 years ago)
- Location:
- NEMO/branches/2020/dev_r13723_KERNEL-01_Amy_Mike_newHPGschemes_KERNEL-08_merge/src
- Files:
-
- 9 edited
-
NST/agrif_oce_update.F90 (modified) (6 diffs)
-
OCE/DYN/dynhpg.F90 (modified) (32 diffs)
-
OCE/DYN/dynspg.F90 (modified) (7 diffs)
-
OCE/DYN/dynspg_exp.F90 (modified) (2 diffs)
-
OCE/DYN/dynspg_ts.F90 (modified) (7 diffs)
-
OCE/ISF/isf_oce.F90 (modified) (9 diffs)
-
OCE/ISF/isfload.F90 (modified) (10 diffs)
-
OCE/ISF/isfstp.F90 (modified) (9 diffs)
-
OCE/oce.F90 (modified) (2 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2020/dev_r13723_KERNEL-01_Amy_Mike_newHPGschemes_KERNEL-08_merge/src/NST/agrif_oce_update.F90
r13286 r13935 886 886 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres 887 887 LOGICAL , INTENT(in ) :: before 888 !! 889 REAL(wp), DIMENSION(jpi,jpj) :: zpgu ! 2D workspace 888 890 !! 889 891 INTEGER :: ji, jj, jk … … 905 907 ! 906 908 ! Update "now" 3d velocities: 907 spgu(ji,jj) = 0._wp909 zpgu(ji,jj) = 0._wp 908 910 DO jk=1,jpkm1 909 spgu(ji,jj) = spgu(ji,jj) + e3u(ji,jj,jk,Kmm_a) * uu(ji,jj,jk,Kmm_a)911 zpgu(ji,jj) = zpgu(ji,jj) + e3u(ji,jj,jk,Kmm_a) * uu(ji,jj,jk,Kmm_a) 910 912 END DO 911 913 ! 912 zcorr = (tabres(ji,jj) - spgu(ji,jj)) * r1_hu(ji,jj,Kmm_a)914 zcorr = (tabres(ji,jj) - zpgu(ji,jj)) * r1_hu(ji,jj,Kmm_a) 913 915 DO jk=1,jpkm1 914 916 uu(ji,jj,jk,Kmm_a) = uu(ji,jj,jk,Kmm_a) + zcorr * umask(ji,jj,jk) … … 925 927 ! 926 928 ! Correct "before" velocities to hold correct bt component: 927 spgu(ji,jj) = 0.e0929 zpgu(ji,jj) = 0.e0 928 930 DO jk=1,jpkm1 929 spgu(ji,jj) = spgu(ji,jj) + e3u(ji,jj,jk,Kbb_a) * uu(ji,jj,jk,Kbb_a)931 zpgu(ji,jj) = zpgu(ji,jj) + e3u(ji,jj,jk,Kbb_a) * uu(ji,jj,jk,Kbb_a) 930 932 END DO 931 933 ! 932 zcorr = uu_b(ji,jj,Kbb_a) - spgu(ji,jj) * r1_hu(ji,jj,Kbb_a)934 zcorr = uu_b(ji,jj,Kbb_a) - zpgu(ji,jj) * r1_hu(ji,jj,Kbb_a) 933 935 DO jk=1,jpkm1 934 936 uu(ji,jj,jk,Kbb_a) = uu(ji,jj,jk,Kbb_a) + zcorr * umask(ji,jj,jk) … … 953 955 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres 954 956 LOGICAL , INTENT(in ) :: before 957 ! 958 REAL(wp), DIMENSION(jpi,jpj) :: zpgv ! 2D workspace 955 959 ! 956 960 INTEGER :: ji, jj, jk … … 971 975 ! 972 976 ! Update "now" 3d velocities: 973 spgv(ji,jj) = 0.e0977 zpgv(ji,jj) = 0.e0 974 978 DO jk=1,jpkm1 975 spgv(ji,jj) = spgv(ji,jj) + e3v(ji,jj,jk,Kmm_a) * vv(ji,jj,jk,Kmm_a)979 zpgv(ji,jj) = zpgv(ji,jj) + e3v(ji,jj,jk,Kmm_a) * vv(ji,jj,jk,Kmm_a) 976 980 END DO 977 981 ! 978 zcorr = (tabres(ji,jj) - spgv(ji,jj)) * r1_hv(ji,jj,Kmm_a)982 zcorr = (tabres(ji,jj) - zpgv(ji,jj)) * r1_hv(ji,jj,Kmm_a) 979 983 DO jk=1,jpkm1 980 984 vv(ji,jj,jk,Kmm_a) = vv(ji,jj,jk,Kmm_a) + zcorr * vmask(ji,jj,jk) … … 991 995 ! 992 996 ! Correct "before" velocities to hold correct bt component: 993 spgv(ji,jj) = 0.e0997 zpgv(ji,jj) = 0.e0 994 998 DO jk=1,jpkm1 995 spgv(ji,jj) = spgv(ji,jj) + e3v(ji,jj,jk,Kbb_a) * vv(ji,jj,jk,Kbb_a)999 zpgv(ji,jj) = zpgv(ji,jj) + e3v(ji,jj,jk,Kbb_a) * vv(ji,jj,jk,Kbb_a) 996 1000 END DO 997 1001 ! 998 zcorr = vv_b(ji,jj,Kbb_a) - spgv(ji,jj) * r1_hv(ji,jj,Kbb_a)1002 zcorr = vv_b(ji,jj,Kbb_a) - zpgv(ji,jj) * r1_hv(ji,jj,Kbb_a) 999 1003 DO jk=1,jpkm1 1000 1004 vv(ji,jj,jk,Kbb_a) = vv(ji,jj,jk,Kbb_a) + zcorr * vmask(ji,jj,jk) -
NEMO/branches/2020/dev_r13723_KERNEL-01_Amy_Mike_newHPGschemes_KERNEL-08_merge/src/OCE/DYN/dynhpg.F90
r13889 r13935 154 154 !! 155 155 INTEGER :: ji, jj, jk, ikt ! dummy loop indices ISF 156 REAL(wp) :: znad157 156 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: zts_top, zrhd ! hypothesys on isf density 158 157 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zrhdtop_isf ! density at bottom of ISF … … 265 264 INTEGER :: ji, jj, jk ! dummy loop indices 266 265 REAL(wp) :: zcoef0, zcoef1 ! temporary scalars 267 REAL(wp), DIMENSION(jpi,jpj ,jpk) :: zhpi, zhpj266 REAL(wp), DIMENSION(jpi,jpj) :: zhpi, zhpj 268 267 !!---------------------------------------------------------------------- 269 268 ! … … 273 272 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ z-coordinate case ' 274 273 ENDIF 275 276 zcoef0 = - grav * 0.5_wp ! Local constant initialization 277 278 ! Surface value 279 DO_2D( 0, 0, 0, 0 ) 274 ! 275 zcoef0 = - grav * 0.5_wp ! Local constant initialization 276 ! 277 DO_2D( 0, 0, 0, 0 ) ! Surface value 280 278 zcoef1 = zcoef0 * e3w(ji,jj,1,Kmm) 281 ! hydrostatic pressure gradient282 zhpi(ji,jj ,1) = zcoef1 * ( rhd(ji+1,jj,1) - rhd(ji,jj,1) ) * r1_e1u(ji,jj)283 zhpj(ji,jj ,1) = zcoef1 * ( rhd(ji,jj+1,1) - rhd(ji,jj,1) ) * r1_e2v(ji,jj)284 ! add to the general momentum trend285 puu(ji,jj,1,Krhs) = puu(ji,jj,1,Krhs) + zhpi(ji,jj ,1)286 pvv(ji,jj,1,Krhs) = pvv(ji,jj,1,Krhs) + zhpj(ji,jj ,1)279 ! ! hydrostatic pressure gradient 280 zhpi(ji,jj) = zcoef1 * ( rhd(ji+1,jj,1) - rhd(ji,jj,1) ) * r1_e1u(ji,jj) 281 zhpj(ji,jj) = zcoef1 * ( rhd(ji,jj+1,1) - rhd(ji,jj,1) ) * r1_e2v(ji,jj) 282 ! ! add to the general momentum trend 283 puu(ji,jj,1,Krhs) = puu(ji,jj,1,Krhs) + zhpi(ji,jj) 284 pvv(ji,jj,1,Krhs) = pvv(ji,jj,1,Krhs) + zhpj(ji,jj) 287 285 END_2D 288 289 ! 290 ! interior value (2=<jk=<jpkm1) 291 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 286 ! 287 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! interior value (2=<jk=<jpkm1) 292 288 zcoef1 = zcoef0 * e3w(ji,jj,jk,Kmm) 293 ! hydrostatic pressure gradient 294 zhpi(ji,jj,jk) = zhpi(ji,jj,jk-1) & 295 & + zcoef1 * ( ( rhd(ji+1,jj,jk)+rhd(ji+1,jj,jk-1) ) & 296 & - ( rhd(ji ,jj,jk)+rhd(ji ,jj,jk-1) ) ) * r1_e1u(ji,jj) 297 298 zhpj(ji,jj,jk) = zhpj(ji,jj,jk-1) & 299 & + zcoef1 * ( ( rhd(ji,jj+1,jk)+rhd(ji,jj+1,jk-1) ) & 300 & - ( rhd(ji,jj, jk)+rhd(ji,jj ,jk-1) ) ) * r1_e2v(ji,jj) 301 ! add to the general momentum trend 302 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + zhpi(ji,jj,jk) 303 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + zhpj(ji,jj,jk) 289 ! ! hydrostatic pressure gradient 290 zhpi(ji,jj) = zhpi(ji,jj) + zcoef1 * ( ( rhd(ji+1,jj,jk)+rhd(ji+1,jj,jk-1) ) & 291 & - ( rhd(ji ,jj,jk)+rhd(ji ,jj,jk-1) ) ) * r1_e1u(ji,jj) 292 293 zhpj(ji,jj) = zhpj(ji,jj) + zcoef1 * ( ( rhd(ji,jj+1,jk)+rhd(ji,jj+1,jk-1) ) & 294 & - ( rhd(ji,jj, jk)+rhd(ji,jj ,jk-1) ) ) * r1_e2v(ji,jj) 295 ! ! add to the general momentum trend 296 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + zhpi(ji,jj) 297 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + zhpj(ji,jj) 304 298 END_3D 305 299 ! … … 410 404 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) :: puu, pvv ! ocean velocities and RHS of momentum equation 411 405 !! 412 INTEGER :: ji, jj, jk, jii, jjj 413 REAL(wp) :: zcoef0, zuap, zvap, z nad, ztmp ! temporaryscalars414 LOGICAL :: ll_tmp1, ll_tmp2 406 INTEGER :: ji, jj, jk, jii, jjj ! dummy loop indices 407 REAL(wp) :: zcoef0, zuap, zvap, ztmp ! local scalars 408 LOGICAL :: ll_tmp1, ll_tmp2 ! local logical variables 415 409 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zhpi, zhpj 416 410 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: zcpx, zcpy !W/D pressure filter … … 426 420 ! 427 421 zcoef0 = - grav * 0.5_wp 428 IF ( ln_linssh ) THEN ; znad = 0._wp ! Fixed volume: density anomaly429 ELSE ; znad = 1._wp ! Variable volume: density430 ENDIF431 422 ! 432 423 IF( ln_wd_il ) THEN … … 470 461 CALL lbc_lnk_multi( 'dynhpg', zcpx, 'U', 1.0_wp, zcpy, 'V', 1.0_wp ) 471 462 END IF 472 473 ! Surface value 474 DO_2D( 0, 0, 0, 0 ) 475 ! hydrostatic pressure gradient along s-surfaces 476 zhpi(ji,jj,1) = & 477 & zcoef0 * ( e3w(ji+1,jj ,1,Kmm) * ( znad + rhd(ji+1,jj ,1) ) & 478 & - e3w(ji ,jj ,1,Kmm) * ( znad + rhd(ji ,jj ,1) ) ) & 479 & * r1_e1u(ji,jj) 480 zhpj(ji,jj,1) = & 481 & zcoef0 * ( e3w(ji ,jj+1,1,Kmm) * ( znad + rhd(ji ,jj+1,1) ) & 482 & - e3w(ji ,jj ,1,Kmm) * ( znad + rhd(ji ,jj ,1) ) ) & 483 & * r1_e2v(ji,jj) 484 ! s-coordinate pressure gradient correction 485 zuap = -zcoef0 * ( rhd (ji+1,jj,1) + rhd (ji,jj,1) + 2._wp * znad ) & 463 ! 464 DO_2D( 0, 0, 0, 0 ) ! Surface value 465 ! ! hydrostatic pressure gradient along s-surfaces 466 zhpi(ji,jj,1) = zcoef0 * r1_e1u(ji,jj) & 467 & * ( e3w(ji+1,jj ,1,Kmm) * rhd(ji+1,jj ,1) & 468 & - e3w(ji ,jj ,1,Kmm) * rhd(ji ,jj ,1) ) 469 zhpj(ji,jj,1) = zcoef0 * r1_e2v(ji,jj) & 470 & * ( e3w(ji ,jj+1,1,Kmm) * rhd(ji ,jj+1,1) & 471 & - e3w(ji ,jj ,1,Kmm) * rhd(ji ,jj ,1) ) 472 ! ! s-coordinate pressure gradient correction 473 zuap = -zcoef0 * ( rhd (ji+1,jj,1) + rhd (ji,jj,1) ) & 486 474 & * ( gde3w(ji+1,jj,1) - gde3w(ji,jj,1) ) * r1_e1u(ji,jj) 487 zvap = -zcoef0 * ( rhd (ji,jj+1,1) + rhd (ji,jj,1) + 2._wp * znad) &475 zvap = -zcoef0 * ( rhd (ji,jj+1,1) + rhd (ji,jj,1) ) & 488 476 & * ( gde3w(ji,jj+1,1) - gde3w(ji,jj,1) ) * r1_e2v(ji,jj) 489 477 ! … … 494 482 zvap = zvap * zcpy(ji,jj) 495 483 ENDIF 496 ! 497 ! add to the general momentum trend 484 ! ! add to the general momentum trend 498 485 puu(ji,jj,1,Krhs) = puu(ji,jj,1,Krhs) + zhpi(ji,jj,1) + zuap 499 486 pvv(ji,jj,1,Krhs) = pvv(ji,jj,1,Krhs) + zhpj(ji,jj,1) + zvap 500 487 END_2D 501 502 ! interior value (2=<jk=<jpkm1) 503 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 504 ! hydrostatic pressure gradient along s-surfaces 505 zhpi(ji,jj,jk) = zhpi(ji,jj,jk-1) + zcoef0 * r1_e1u(ji,jj) & 506 & * ( e3w(ji+1,jj,jk,Kmm) * ( rhd(ji+1,jj,jk) + rhd(ji+1,jj,jk-1) + 2*znad ) & 507 & - e3w(ji ,jj,jk,Kmm) * ( rhd(ji ,jj,jk) + rhd(ji ,jj,jk-1) + 2*znad ) ) 508 zhpj(ji,jj,jk) = zhpj(ji,jj,jk-1) + zcoef0 * r1_e2v(ji,jj) & 509 & * ( e3w(ji,jj+1,jk,Kmm) * ( rhd(ji,jj+1,jk) + rhd(ji,jj+1,jk-1) + 2*znad ) & 510 & - e3w(ji,jj ,jk,Kmm) * ( rhd(ji,jj, jk) + rhd(ji,jj ,jk-1) + 2*znad ) ) 511 ! s-coordinate pressure gradient correction 512 zuap = -zcoef0 * ( rhd (ji+1,jj ,jk) + rhd (ji,jj,jk) + 2._wp * znad ) & 488 ! 489 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! interior value (2=<jk=<jpkm1) 490 ! ! hydrostatic pressure gradient along s-surfaces 491 zhpi(ji,jj,jk) = zhpi(ji,jj,jk-1) + zcoef0 * r1_e1u(ji,jj) & 492 & * ( e3w(ji+1,jj,jk,Kmm) * ( rhd(ji+1,jj,jk) + rhd(ji+1,jj,jk-1) ) & 493 & - e3w(ji ,jj,jk,Kmm) * ( rhd(ji ,jj,jk) + rhd(ji ,jj,jk-1) ) ) 494 zhpj(ji,jj,jk) = zhpj(ji,jj,jk-1) + zcoef0 * r1_e2v(ji,jj) & 495 & * ( e3w(ji,jj+1,jk,Kmm) * ( rhd(ji,jj+1,jk) + rhd(ji,jj+1,jk-1) ) & 496 & - e3w(ji,jj ,jk,Kmm) * ( rhd(ji,jj, jk) + rhd(ji,jj ,jk-1) ) ) 497 ! ! s-coordinate pressure gradient correction 498 zuap = -zcoef0 * ( rhd (ji+1,jj ,jk) + rhd (ji,jj,jk) ) & 513 499 & * ( gde3w(ji+1,jj ,jk) - gde3w(ji,jj,jk) ) * r1_e1u(ji,jj) 514 zvap = -zcoef0 * ( rhd (ji ,jj+1,jk) + rhd (ji,jj,jk) + 2._wp * znad )&500 zvap = -zcoef0 * ( rhd (ji ,jj+1,jk) + rhd (ji,jj,jk) ) & 515 501 & * ( gde3w(ji ,jj+1,jk) - gde3w(ji,jj,jk) ) * r1_e2v(ji,jj) 516 502 ! … … 555 541 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) :: puu, pvv ! ocean velocities and RHS of momentum equation 556 542 !! 557 INTEGER :: ji, jj, jk, ikt, iktp1i, iktp1j ! dummy loop indices 558 REAL(wp) :: zcoef0, zuap, zvap, znad ! temporary scalars 543 INTEGER :: ji, jj, jk ! dummy loop indices 544 INTEGER :: ikt , ikti1, iktj1 ! local integer 545 REAL(wp) :: ze3w, ze3wi1, ze3wj1 ! local scalars 546 REAL(wp) :: zcoef0, zuap, zvap ! - - 559 547 REAL(wp), DIMENSION(jpi,jpj,jpk ) :: zhpi, zhpj 560 548 REAL(wp), DIMENSION(jpi,jpj,jpts) :: zts_top … … 563 551 ! 564 552 zcoef0 = - grav * 0.5_wp ! Local constant initialization 565 !566 znad=1._wp ! To use density and not density anomaly567 553 ! 568 554 ! ! iniitialised to 0. zhpi zhpi … … 580 566 CALL eos( zts_top, risfdep, zrhdtop_oce ) 581 567 582 !================================================================================== 583 !===== Compute surface value ===================================================== 584 !================================================================================== 568 ! !===========================! 569 ! !===== surface value =====! 570 ! !===========================! 585 571 DO_2D( 0, 0, 0, 0 ) 586 ikt = mikt(ji,jj) 587 iktp1i = mikt(ji+1,jj) 588 iktp1j = mikt(ji,jj+1) 589 ! hydrostatic pressure gradient along s-surfaces and ice shelf pressure 590 ! we assume ISF is in isostatic equilibrium 591 zhpi(ji,jj,1) = zcoef0 / e1u(ji,jj) * ( 0.5_wp * e3w(ji+1,jj,iktp1i,Kmm) & 592 & * ( 2._wp * znad + rhd(ji+1,jj,iktp1i) + zrhdtop_oce(ji+1,jj) ) & 593 & - 0.5_wp * e3w(ji,jj,ikt,Kmm) & 594 & * ( 2._wp * znad + rhd(ji,jj,ikt) + zrhdtop_oce(ji,jj) ) & 595 & + ( risfload(ji+1,jj) - risfload(ji,jj)) ) 596 zhpj(ji,jj,1) = zcoef0 / e2v(ji,jj) * ( 0.5_wp * e3w(ji,jj+1,iktp1j,Kmm) & 597 & * ( 2._wp * znad + rhd(ji,jj+1,iktp1j) + zrhdtop_oce(ji,jj+1) ) & 598 & - 0.5_wp * e3w(ji,jj,ikt,Kmm) & 599 & * ( 2._wp * znad + rhd(ji,jj,ikt) + zrhdtop_oce(ji,jj) ) & 600 & + ( risfload(ji,jj+1) - risfload(ji,jj)) ) 601 ! s-coordinate pressure gradient correction (=0 if z coordinate) 602 zuap = -zcoef0 * ( rhd (ji+1,jj,1) + rhd (ji,jj,1) + 2._wp * znad ) & 572 ikt = mikt(ji ,jj ) ; ze3w = e3w(ji ,jj ,ikt ,Kmm) 573 ikti1 = mikt(ji+1,jj ) ; ze3wi1 = e3w(ji+1,jj ,ikti1,Kmm) 574 iktj1 = mikt(ji ,jj+1) ; ze3wj1 = e3w(ji ,jj+1,iktj1,Kmm) 575 ! ! hydrostatic pressure gradient along s-surfaces and ice shelf pressure 576 ! ! we assume ISF is in isostatic equilibrium 577 zhpi(ji,jj,1) = zcoef0 * r1_e1u(ji,jj) * ( risfload(ji+1,jj) - risfload(ji,jj) & 578 & + 0.5_wp * ( ze3wi1 * ( rhd(ji+1,jj,ikti1) + zrhdtop_oce(ji+1,jj) ) & 579 & - ze3w * ( rhd(ji ,jj,ikt ) + zrhdtop_oce(ji ,jj) ) ) ) 580 zhpj(ji,jj,1) = zcoef0 * r1_e2v(ji,jj) * ( risfload(ji,jj+1) - risfload(ji,jj) & 581 & + 0.5_wp * ( ze3wj1 * ( rhd(ji,jj+1,iktj1) + zrhdtop_oce(ji,jj+1) ) & 582 & - ze3w * ( rhd(ji,jj ,ikt ) + zrhdtop_oce(ji,jj ) ) ) ) 583 ! ! s-coordinate pressure gradient correction (=0 if z coordinate) 584 zuap = -zcoef0 * ( rhd (ji+1,jj,1) + rhd (ji,jj,1) ) & 603 585 & * ( gde3w(ji+1,jj,1) - gde3w(ji,jj,1) ) * r1_e1u(ji,jj) 604 zvap = -zcoef0 * ( rhd (ji,jj+1,1) + rhd (ji,jj,1) + 2._wp * znad) &586 zvap = -zcoef0 * ( rhd (ji,jj+1,1) + rhd (ji,jj,1) ) & 605 587 & * ( gde3w(ji,jj+1,1) - gde3w(ji,jj,1) ) * r1_e2v(ji,jj) 606 ! add to the general momentum trend588 ! ! add to the general momentum trend 607 589 puu(ji,jj,1,Krhs) = puu(ji,jj,1,Krhs) + (zhpi(ji,jj,1) + zuap) * umask(ji,jj,1) 608 590 pvv(ji,jj,1,Krhs) = pvv(ji,jj,1,Krhs) + (zhpj(ji,jj,1) + zvap) * vmask(ji,jj,1) 609 591 END_2D 610 !==================================================================================611 !===== Compute interior value ===================================================== 612 !================================================================================== 613 ! interior value (2=<jk=<jpkm1)592 ! 593 ! !=============================! 594 ! !===== interior values =====! 595 ! !=============================! 614 596 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 615 ! hydrostatic pressure gradient along s-surfaces 597 ze3w = e3w(ji ,jj ,jk,Kmm) 598 ze3wi1 = e3w(ji+1,jj ,jk,Kmm) 599 ze3wj1 = e3w(ji ,jj+1,jk,Kmm) 600 ! ! hydrostatic pressure gradient along s-surfaces 616 601 zhpi(ji,jj,jk) = zhpi(ji,jj,jk-1) + zcoef0 / e1u(ji,jj) & 617 & * ( e3w(ji+1,jj,jk,Kmm) & 618 & * ( rhd(ji+1,jj,jk) + rhd(ji+1,jj,jk-1) + 2*znad ) * wmask(ji+1,jj,jk) & 619 & - e3w(ji ,jj,jk,Kmm) & 620 & * ( rhd(ji ,jj,jk) + rhd(ji ,jj,jk-1) + 2*znad ) * wmask(ji ,jj,jk) ) 602 & * ( ze3wi1 * ( rhd(ji+1,jj,jk) + rhd(ji+1,jj,jk-1) ) * wmask(ji+1,jj,jk) & 603 & - ze3w * ( rhd(ji ,jj,jk) + rhd(ji ,jj,jk-1) ) * wmask(ji ,jj,jk) ) 621 604 zhpj(ji,jj,jk) = zhpj(ji,jj,jk-1) + zcoef0 / e2v(ji,jj) & 622 & * ( e3w(ji,jj+1,jk,Kmm) & 623 & * ( rhd(ji,jj+1,jk) + rhd(ji,jj+1,jk-1) + 2*znad ) * wmask(ji,jj+1,jk) & 624 & - e3w(ji,jj ,jk,Kmm) & 625 & * ( rhd(ji,jj, jk) + rhd(ji,jj ,jk-1) + 2*znad ) * wmask(ji,jj ,jk) ) 626 ! s-coordinate pressure gradient correction 627 zuap = -zcoef0 * ( rhd (ji+1,jj ,jk) + rhd (ji,jj,jk) + 2._wp * znad ) & 605 & * ( ze3wj1 * ( rhd(ji,jj+1,jk) + rhd(ji,jj+1,jk-1) ) * wmask(ji,jj+1,jk) & 606 & - ze3w * ( rhd(ji,jj, jk) + rhd(ji,jj ,jk-1) ) * wmask(ji,jj ,jk) ) 607 ! ! s-coordinate pressure gradient correction 608 zuap = -zcoef0 * ( rhd (ji+1,jj ,jk) + rhd (ji,jj,jk) ) & 628 609 & * ( gde3w(ji+1,jj ,jk) - gde3w(ji,jj,jk) ) / e1u(ji,jj) 629 zvap = -zcoef0 * ( rhd (ji ,jj+1,jk) + rhd (ji,jj,jk) + 2._wp * znad) &610 zvap = -zcoef0 * ( rhd (ji ,jj+1,jk) + rhd (ji,jj,jk) ) & 630 611 & * ( gde3w(ji ,jj+1,jk) - gde3w(ji,jj,jk) ) / e2v(ji,jj) 631 ! add to the general momentum trend612 ! ! add to the general momentum trend 632 613 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + (zhpi(ji,jj,jk) + zuap) * umask(ji,jj,jk) 633 614 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + (zhpj(ji,jj,jk) + zvap) * vmask(ji,jj,jk) … … 657 638 LOGICAL :: ll_tmp1, ll_tmp2 ! local logical variables 658 639 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zhpi, zhpj 659 REAL(wp), DIMENSION(jpi,jpj,jpk) :: rhd_opt660 640 661 641 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zdzx, zdzy, zdzz ! Primitive grid differences ('delta_xyz') … … 720 700 z1_12 = 1.0_wp / 12._wp 721 701 722 723 !! To be removed once combined with KERNEL08. KERNEL08 action removes znad from the hpg module, spg is calculated in spg moduel instead. The final code will just have rhd(:,:,:), no rhd_opt(:,:,:)724 IF( .NOT. ln_linssh ) THEN725 rhd_opt(:,:,:) = rhd(:,:,:) + 1.0_wp ! for vvl option726 ELSE727 rhd_opt(:,:,:) = rhd(:,:,:)728 END IF729 730 702 !---------------------------------------------------------------------------------------- 731 703 ! 1. compute and store elementary vertical differences in provisional arrays … … 735 707 736 708 DO_3D( 1, 1, 1, 1, 2, jpkm1 ) 737 zdrhoz(ji,jj,jk) = rhd _opt (ji ,jj ,jk) - rhd_opt(ji,jj,jk-1)709 zdrhoz(ji,jj,jk) = rhd (ji ,jj ,jk) - rhd (ji,jj,jk-1) 738 710 zdzz (ji,jj,jk) = - gde3w(ji ,jj ,jk) + gde3w(ji,jj,jk-1) 739 711 END_3D … … 762 734 763 735 ! mb for sea-ice shelves we will need to re-write this upper boundary condition in the same form as the lower boundary condition 764 zdrho_k(:,:,1) = aco_bc_vrt * ( rhd _opt (:,:,2) - rhd_opt(:,:,1) ) - bco_bc_vrt * zdrho_k(:,:,2)736 zdrho_k(:,:,1) = aco_bc_vrt * ( rhd (:,:,2) - rhd (:,:,1) ) - bco_bc_vrt * zdrho_k(:,:,2) 765 737 zdz_k (:,:,1) = aco_bc_vrt * (-gde3w(:,:,2) + gde3w(:,:,1) ) - bco_bc_vrt * zdz_k (:,:,2) 766 738 … … 768 740 IF ( mbkt(ji,jj)>1 ) THEN 769 741 iktb = mbkt(ji,jj) 770 zdrho_k(ji,jj,iktb) = aco_bc_vrt * ( rhd _opt(ji,jj,iktb) - rhd_opt(ji,jj,iktb-1) ) - bco_bc_vrt * zdrho_k(ji,jj,iktb-1)742 zdrho_k(ji,jj,iktb) = aco_bc_vrt * ( rhd(ji,jj,iktb) - rhd(ji,jj,iktb-1) ) - bco_bc_vrt * zdrho_k(ji,jj,iktb-1) 771 743 zdz_k (ji,jj,iktb) = aco_bc_vrt * (-gde3w(ji,jj,iktb) + gde3w(ji,jj,iktb-1) ) - bco_bc_vrt * zdz_k (ji,jj,iktb-1) 772 744 END IF … … 786 758 DO_2D( 0, 1, 0, 1) 787 759 z_rho_k(ji,jj,1) = grav * ( ssh(ji,jj,Kmm) + gde3w(ji,jj,1) ) & 788 & * ( rhd _opt(ji,jj,1) &789 & + 0.5_wp * ( rhd _opt (ji,jj,2) - rhd_opt(ji,jj,1) ) &760 & * ( rhd(ji,jj,1) & 761 & + 0.5_wp * ( rhd (ji,jj,2) - rhd (ji,jj,1) ) & 790 762 & * ( ssh (ji,jj,Kmm) + gde3w(ji,jj,1) ) & 791 763 & / ( - gde3w(ji,jj,2) + gde3w(ji,jj,1) ) ) … … 797 769 798 770 DO_3D( 0, 1, 0, 1, 2, jpkm1 ) 799 z_rho_k(ji,jj,jk) = zcoef0 * ( rhd _opt (ji,jj,jk) + rhd_opt(ji,jj,jk-1) ) &771 z_rho_k(ji,jj,jk) = zcoef0 * ( rhd (ji,jj,jk) + rhd (ji,jj,jk-1) ) & 800 772 & * ( - gde3w(ji,jj,jk) + gde3w(ji,jj,jk-1) ) & 801 773 & + z_grav_10 * ( & … … 803 775 & * ( - gde3w(ji,jj,jk) + gde3w(ji,jj,jk-1) - z1_12 * ( zdz_k (ji,jj,jk) + zdz_k (ji,jj,jk-1) ) ) & 804 776 & - ( zdz_k (ji,jj,jk) - zdz_k (ji,jj,jk-1) ) & 805 & * ( rhd _opt (ji,jj,jk) - rhd_opt(ji,jj,jk-1) - z1_12 * ( zdrho_k(ji,jj,jk) + zdrho_k(ji,jj,jk-1) ) ) &777 & * ( rhd (ji,jj,jk) - rhd (ji,jj,jk-1) - z1_12 * ( zdrho_k(ji,jj,jk) + zdrho_k(ji,jj,jk-1) ) ) & 806 778 & ) 807 779 END_3D … … 812 784 813 785 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 814 zdrhox(ji,jj,jk) = rhd _opt (ji+1,jj ,jk) - rhd_opt(ji,jj,jk )786 zdrhox(ji,jj,jk) = rhd (ji+1,jj ,jk) - rhd (ji,jj,jk ) 815 787 zdzx (ji,jj,jk) = - gde3w(ji+1,jj ,jk) + gde3w(ji,jj,jk ) 816 zdrhoy(ji,jj,jk) = rhd _opt (ji ,jj+1,jk) - rhd_opt(ji,jj,jk )788 zdrhoy(ji,jj,jk) = rhd (ji ,jj+1,jk) - rhd (ji,jj,jk ) 817 789 zdzy (ji,jj,jk) = - gde3w(ji ,jj+1,jk) + gde3w(ji,jj,jk ) 818 790 END_3D … … 869 841 IF (ji < jpi) THEN 870 842 IF ( umask(ji,jj,jk) > 0.5_wp .AND. umask(ji-1,jj,jk) < 0.5_wp .AND. umask(ji+1,jj,jk) > 0.5_wp) THEN 871 zz_drho_i(ji,jj) = aco_bc_hor * ( rhd _opt (ji+1,jj,jk) - rhd_opt(ji,jj,jk) ) - bco_bc_hor * zdrho_i(ji+1,jj,jk)843 zz_drho_i(ji,jj) = aco_bc_hor * ( rhd (ji+1,jj,jk) - rhd (ji,jj,jk) ) - bco_bc_hor * zdrho_i(ji+1,jj,jk) 872 844 zz_dz_i (ji,jj) = aco_bc_hor * (-gde3w(ji+1,jj,jk) + gde3w(ji,jj,jk) ) - bco_bc_hor * zdz_i (ji+1,jj,jk) 873 845 END IF … … 876 848 IF (ji > 2) THEN 877 849 IF ( umask(ji,jj,jk) < 0.5_wp .AND. umask(ji-1,jj,jk) > 0.5_wp .AND. umask(ji-2,jj,jk) > 0.5_wp) THEN 878 zz_drho_i(ji,jj) = aco_bc_hor * ( rhd _opt (ji,jj,jk) - rhd_opt(ji-1,jj,jk) ) - bco_bc_hor * zdrho_i(ji-1,jj,jk)850 zz_drho_i(ji,jj) = aco_bc_hor * ( rhd (ji,jj,jk) - rhd (ji-1,jj,jk) ) - bco_bc_hor * zdrho_i(ji-1,jj,jk) 879 851 zz_dz_i (ji,jj) = aco_bc_hor * (-gde3w(ji,jj,jk) + gde3w(ji-1,jj,jk) ) - bco_bc_hor * zdz_i (ji-1,jj,jk) 880 852 END IF … … 883 855 IF (jj < jpj) THEN 884 856 IF ( vmask(ji,jj,jk) > 0.5_wp .AND. vmask(ji,jj-1,jk) < 0.5_wp .AND. vmask(ji,jj+1,jk) > 0.5_wp) THEN 885 zz_drho_j(ji,jj) = aco_bc_hor * ( rhd _opt (ji,jj+1,jk) - rhd_opt(ji,jj,jk) ) - bco_bc_hor * zdrho_j(ji,jj+1,jk)857 zz_drho_j(ji,jj) = aco_bc_hor * ( rhd (ji,jj+1,jk) - rhd (ji,jj,jk) ) - bco_bc_hor * zdrho_j(ji,jj+1,jk) 886 858 zz_dz_j (ji,jj) = aco_bc_hor * (-gde3w(ji,jj+1,jk) + gde3w(ji,jj,jk) ) - bco_bc_hor * zdz_j (ji,jj+1,jk) 887 859 END IF … … 890 862 IF (jj > 2) THEN 891 863 IF ( vmask(ji,jj,jk) < 0.5_wp .AND. vmask(ji,jj-1,jk) > 0.5_wp .AND. vmask(ji,jj-2,jk) > 0.5_wp) THEN 892 zz_drho_j(ji,jj) = aco_bc_hor * ( rhd _opt (ji,jj,jk) - rhd_opt(ji,jj-1,jk) ) - bco_bc_hor * zdrho_j(ji,jj-1,jk)864 zz_drho_j(ji,jj) = aco_bc_hor * ( rhd (ji,jj,jk) - rhd (ji,jj-1,jk) ) - bco_bc_hor * zdrho_j(ji,jj-1,jk) 893 865 zz_dz_j (ji,jj) = aco_bc_hor * (-gde3w(ji,jj,jk) + gde3w(ji,jj-1,jk) ) - bco_bc_hor * zdz_j (ji,jj-1,jk) 894 866 END IF … … 907 879 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 908 880 ! two -ve signs cancel in next two lines (within zcoef0 and because gde3w is a depth not a height) 909 z_rho_i(ji,jj,jk) = zcoef0 * ( rhd _opt (ji+1,jj,jk) + rhd_opt(ji,jj,jk) ) &881 z_rho_i(ji,jj,jk) = zcoef0 * ( rhd (ji+1,jj,jk) + rhd (ji,jj,jk) ) & 910 882 & * ( gde3w(ji+1,jj,jk) - gde3w(ji,jj,jk) ) 911 883 IF ( umask(ji-1, jj, jk) > 0.5 .OR. umask(ji+1, jj, jk) > 0.5 ) THEN … … 914 886 & * ( - gde3w(ji+1,jj,jk) + gde3w(ji,jj,jk) - z1_12 * ( zdz_i (ji+1,jj,jk) + zdz_i (ji,jj,jk) ) ) & 915 887 & - ( zdz_i (ji+1,jj,jk) - zdz_i (ji,jj,jk) ) & 916 & * ( rhd _opt (ji+1,jj,jk) - rhd_opt(ji,jj,jk) - z1_12 * ( zdrho_i(ji+1,jj,jk) + zdrho_i(ji,jj,jk) ) ) &888 & * ( rhd (ji+1,jj,jk) - rhd (ji,jj,jk) - z1_12 * ( zdrho_i(ji+1,jj,jk) + zdrho_i(ji,jj,jk) ) ) & 917 889 & ) 918 890 END IF 919 891 920 z_rho_j(ji,jj,jk) = zcoef0 * ( rhd _opt (ji,jj+1,jk) + rhd_opt(ji,jj,jk) ) &892 z_rho_j(ji,jj,jk) = zcoef0 * ( rhd (ji,jj+1,jk) + rhd (ji,jj,jk) ) & 921 893 & * ( gde3w(ji,jj+1,jk) - gde3w(ji,jj,jk) ) 922 894 IF ( vmask(ji, jj-1, jk) > 0.5 .OR. vmask(ji, jj+1, jk) > 0.5 ) THEN … … 925 897 & * ( - gde3w(ji,jj+1,jk) + gde3w(ji,jj,jk) - z1_12 * ( zdz_j (ji,jj+1,jk) + zdz_j (ji,jj,jk) ) ) & 926 898 & - ( zdz_j (ji,jj+1,jk) - zdz_j (ji,jj,jk) ) & 927 & * ( rhd _opt (ji,jj+1,jk) - rhd_opt(ji,jj,jk) - z1_12 * ( zdrho_j(ji,jj+1,jk) + zdrho_j(ji,jj,jk) ) ) &899 & * ( rhd (ji,jj+1,jk) - rhd (ji,jj,jk) - z1_12 * ( zdrho_j(ji,jj+1,jk) + zdrho_j(ji,jj,jk) ) ) & 928 900 & ) 929 901 END IF … … 999 971 REAL(wp) :: zrhdt1 1000 972 REAL(wp) :: zdpdx1, zdpdx2, zdpdy1, zdpdy2 973 REAL(wp), DIMENSION(jpi,jpj) :: zpgu, zpgv ! 2D workspace 974 REAL(wp), DIMENSION(jpi,jpj) :: zsshu_n, zsshv_n 1001 975 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zdept, zrhh 1002 976 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zhpi, zu, zv, fsp, xsp, asp, bsp, csp, dsp 1003 REAL(wp), DIMENSION(jpi,jpj) :: zsshu_n, zsshv_n1004 977 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: zcpx, zcpy !W/D pressure filter 1005 978 !!---------------------------------------------------------------------- … … 1014 987 zcoef0 = - grav 1015 988 znad = 1._wp 1016 IF( ln_linssh ) znad = 0._wp 1017 989 IF( ln_linssh ) znad = 1._wp 990 ! 991 ! --------------- 992 ! Surface pressure gradient to be removed 993 ! --------------- 994 DO_2D( 0, 0, 0, 0 ) 995 zpgu(ji,jj) = - grav * ( ssh(ji+1,jj,Kmm) - ssh(ji,jj,Kmm) ) * r1_e1u(ji,jj) 996 zpgv(ji,jj) = - grav * ( ssh(ji,jj+1,Kmm) - ssh(ji,jj,Kmm) ) * r1_e2v(ji,jj) 997 END_2D 998 ! 1018 999 IF( ln_wd_il ) THEN 1019 1000 ALLOCATE( zcpx(jpi,jpj) , zcpy(jpi,jpj) ) … … 1081 1062 ! Transfer the depth of "T(:,:,:)" to vertical coordinate "zdept(:,:,:)" 1082 1063 DO_2D( 1, 1, 1, 1 ) 1083 zdept(ji,jj,1) = 0.5_wp * e3w(ji,jj,1,Kmm) - ssh(ji,jj,Kmm) * znad1064 zdept(ji,jj,1) = 0.5_wp * e3w(ji,jj,1,Kmm) - ssh(ji,jj,Kmm) 1084 1065 END_2D 1085 1066 … … 1132 1113 1133 1114 DO_2D( 0, 0, 0, 0 ) 1134 zu(ji,jj,1) = - ( e3u(ji,jj,1,Kmm) - zsshu_n(ji,jj) * znad)1135 zv(ji,jj,1) = - ( e3v(ji,jj,1,Kmm) - zsshv_n(ji,jj) * znad)1115 zu(ji,jj,1) = - ( e3u(ji,jj,1,Kmm) - zsshu_n(ji,jj) ) 1116 zv(ji,jj,1) = - ( e3v(ji,jj,1,Kmm) - zsshv_n(ji,jj) ) 1136 1117 END_2D 1137 1118 … … 1215 1196 zdpdx2 = zdpdx2 * zcpx(ji,jj) * wdrampu(ji,jj) 1216 1197 ENDIF 1217 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + (zdpdx1 + zdpdx2 ) * umask(ji,jj,jk)1198 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + (zdpdx1 + zdpdx2 - zpgu(ji,jj)) * umask(ji,jj,jk) 1218 1199 ENDIF 1219 1200 … … 1275 1256 ENDIF 1276 1257 1277 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + (zdpdy1 + zdpdy2 ) * vmask(ji,jj,jk)1258 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + (zdpdy1 + zdpdy2 - zpgv(ji,jj)) * vmask(ji,jj,jk) 1278 1259 ENDIF 1279 1260 ! -
NEMO/branches/2020/dev_r13723_KERNEL-01_Amy_Mike_newHPGschemes_KERNEL-08_merge/src/OCE/DYN/dynspg.F90
r13497 r13935 79 79 INTEGER :: ji, jj, jk ! dummy loop indices 80 80 REAL(wp) :: z2dt, zg_2, zintp, zgrho0r, zld ! local scalars 81 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zpice 82 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdu, ztrdv 81 REAL(wp) , DIMENSION(jpi,jpj) :: zpgu, zpgv ! 2D workspace 82 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zpice 83 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdu, ztrdv 83 84 !!---------------------------------------------------------------------- 84 85 ! … … 96 97 ! 97 98 DO_2D( 0, 0, 0, 0 ) 98 spgu(ji,jj) = 0._wp99 spgv(ji,jj) = 0._wp99 zpgu(ji,jj) = 0._wp 100 zpgv(ji,jj) = 0._wp 100 101 END_2D 101 102 ! … … 103 104 zg_2 = grav * 0.5 104 105 DO_2D( 0, 0, 0, 0 ) ! gradient of Patm using inverse barometer ssh 105 spgu(ji,jj) = spgu(ji,jj) + zg_2 * ( ssh_ib (ji+1,jj) - ssh_ib (ji,jj) &106 zpgu(ji,jj) = zpgu(ji,jj) + zg_2 * ( ssh_ib (ji+1,jj) - ssh_ib (ji,jj) & 106 107 & + ssh_ibb(ji+1,jj) - ssh_ibb(ji,jj) ) * r1_e1u(ji,jj) 107 spgv(ji,jj) = spgv(ji,jj) + zg_2 * ( ssh_ib (ji,jj+1) - ssh_ib (ji,jj) &108 zpgv(ji,jj) = zpgv(ji,jj) + zg_2 * ( ssh_ib (ji,jj+1) - ssh_ib (ji,jj) & 108 109 & + ssh_ibb(ji,jj+1) - ssh_ibb(ji,jj) ) * r1_e2v(ji,jj) 109 110 END_2D … … 118 119 ! 119 120 DO_2D( 0, 0, 0, 0 ) ! add tide potential forcing 120 spgu(ji,jj) = spgu(ji,jj) + grav * ( pot_astro(ji+1,jj) - pot_astro(ji,jj) ) * r1_e1u(ji,jj)121 spgv(ji,jj) = spgv(ji,jj) + grav * ( pot_astro(ji,jj+1) - pot_astro(ji,jj) ) * r1_e2v(ji,jj)121 zpgu(ji,jj) = zpgu(ji,jj) + grav * ( pot_astro(ji+1,jj) - pot_astro(ji,jj) ) * r1_e1u(ji,jj) 122 zpgv(ji,jj) = zpgv(ji,jj) + grav * ( pot_astro(ji,jj+1) - pot_astro(ji,jj) ) * r1_e2v(ji,jj) 122 123 END_2D 123 124 ! … … 125 126 zld = rn_scal_load * grav 126 127 DO_2D( 0, 0, 0, 0 ) ! add scalar approximation for load potential 127 spgu(ji,jj) = spgu(ji,jj) + zld * ( pssh(ji+1,jj,Kmm) - pssh(ji,jj,Kmm) ) * r1_e1u(ji,jj)128 spgv(ji,jj) = spgv(ji,jj) + zld * ( pssh(ji,jj+1,Kmm) - pssh(ji,jj,Kmm) ) * r1_e2v(ji,jj)128 zpgu(ji,jj) = zpgu(ji,jj) + zld * ( pssh(ji+1,jj,Kmm) - pssh(ji,jj,Kmm) ) * r1_e1u(ji,jj) 129 zpgv(ji,jj) = zpgv(ji,jj) + zld * ( pssh(ji,jj+1,Kmm) - pssh(ji,jj,Kmm) ) * r1_e2v(ji,jj) 129 130 END_2D 130 131 ENDIF … … 137 138 zpice(:,:) = ( zintp * snwice_mass(:,:) + ( 1.- zintp ) * snwice_mass_b(:,:) ) * zgrho0r 138 139 DO_2D( 0, 0, 0, 0 ) 139 spgu(ji,jj) = spgu(ji,jj) + ( zpice(ji+1,jj) - zpice(ji,jj) ) * r1_e1u(ji,jj)140 spgv(ji,jj) = spgv(ji,jj) + ( zpice(ji,jj+1) - zpice(ji,jj) ) * r1_e2v(ji,jj)140 zpgu(ji,jj) = zpgu(ji,jj) + ( zpice(ji+1,jj) - zpice(ji,jj) ) * r1_e1u(ji,jj) 141 zpgv(ji,jj) = zpgv(ji,jj) + ( zpice(ji,jj+1) - zpice(ji,jj) ) * r1_e2v(ji,jj) 141 142 END_2D 142 143 DEALLOCATE( zpice ) … … 144 145 ! 145 146 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !== Add all terms to the general trend 146 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + spgu(ji,jj)147 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + spgv(ji,jj)147 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + zpgu(ji,jj) 148 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + zpgv(ji,jj) 148 149 END_3D 149 150 ! -
NEMO/branches/2020/dev_r13723_KERNEL-01_Amy_Mike_newHPGschemes_KERNEL-08_merge/src/OCE/DYN/dynspg_exp.F90
r13497 r13935 60 60 !! 61 61 INTEGER :: ji, jj, jk ! dummy loop indices 62 REAL(wp), DIMENSION(jpi,jpj) :: zpgu, zpgv ! 2D workspace 62 63 !!---------------------------------------------------------------------- 63 64 ! … … 67 68 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ (explicit free surface)' 68 69 ! 69 spgu(:,:) = 0._wp ; spgv(:,:) = 0._wp70 zpgu(:,:) = 0._wp ; zpgv(:,:) = 0._wp 70 71 ! 71 72 IF( .NOT.ln_linssh .AND. lwp ) WRITE(numout,*) ' non linear free surface: spg is included in dynhpg' 72 73 ENDIF 73 74 IF( ln_linssh ) THEN !* linear free surface : add the surface pressure gradient trend 75 ! 76 DO_2D( 0, 0, 0, 0 ) ! now surface pressure gradient 77 spgu(ji,jj) = - grav * ( ssh(ji+1,jj,Kmm) - ssh(ji,jj,Kmm) ) * r1_e1u(ji,jj) 78 spgv(ji,jj) = - grav * ( ssh(ji,jj+1,Kmm) - ssh(ji,jj,Kmm) ) * r1_e2v(ji,jj) 79 END_2D 80 ! 81 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) ! Add it to the general trend 82 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + spgu(ji,jj) 83 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + spgv(ji,jj) 84 END_3D 85 ! 86 ENDIF 74 ! 75 DO_2D( 0, 0, 0, 0 ) 76 zpgu(ji,jj) = - grav * ( ssh(ji+1,jj,Kmm) - ssh(ji,jj,Kmm) ) * r1_e1u(ji,jj) 77 zpgv(ji,jj) = - grav * ( ssh(ji,jj+1,Kmm) - ssh(ji,jj,Kmm) ) * r1_e2v(ji,jj) 78 END_2D 79 ! 80 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 81 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + zpgu(ji,jj) 82 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + zpgv(ji,jj) 83 END_3D 87 84 ! 88 85 END SUBROUTINE dyn_spg_exp -
NEMO/branches/2020/dev_r13723_KERNEL-01_Amy_Mike_newHPGschemes_KERNEL-08_merge/src/OCE/DYN/dynspg_ts.F90
r13546 r13935 162 162 REAL(wp), DIMENSION(jpi,jpj) :: zCdU_u, zCdU_v ! top/bottom stress at u- & v-points 163 163 REAL(wp), DIMENSION(jpi,jpj) :: zhU, zhV ! fluxes 164 #if defined key_qco 164 165 REAL(wp), DIMENSION(jpi, jpj, jpk) :: ze3u, ze3v 166 #endif 165 167 ! 166 168 REAL(wp) :: zwdramp ! local scalar - only used if ln_wd_dl = .True. … … 229 231 ! != zu_frc = 1/H e3*d/dt(Ua) =! (Vertical mean of Ua, the 3D trends) 230 232 ! ! --------------------------- ! 233 #if defined key_qco 231 234 DO jk = 1 , jpk 232 235 ze3u(:,:,jk) = e3u(:,:,jk,Kmm) … … 236 239 zu_frc(:,:) = SUM( ze3u(:,:,:) * uu(:,:,:,Krhs) * umask(:,:,:) , DIM=3 ) * r1_hu(:,:,Kmm) 237 240 zv_frc(:,:) = SUM( ze3v(:,:,:) * vv(:,:,:,Krhs) * vmask(:,:,:) , DIM=3 ) * r1_hv(:,:,Kmm) 241 #else 242 zu_frc(:,:) = SUM( e3u(:,:,:,Kmm) * uu(:,:,:,Krhs) * umask(:,:,:) , DIM=3 ) * r1_hu(:,:,Kmm) 243 zv_frc(:,:) = SUM( e3v(:,:,:,Kmm) * vv(:,:,:,Krhs) * vmask(:,:,:) , DIM=3 ) * r1_hv(:,:,Kmm) 244 #endif 238 245 ! 239 246 ! … … 247 254 !!gm Is it correct to do so ? I think so... 248 255 249 ! != remove 2D Coriolis and pressure gradient trends=!250 ! ! -------------------------- -----------------------!256 ! != remove 2D Coriolis trend =! 257 ! ! -------------------------- ! 251 258 ! 252 259 IF( kt == nit000 .OR. .NOT. ln_linssh ) CALL dyn_cor_2D_init( Kmm ) ! Set zwz, the barotropic Coriolis force coefficient 253 ! ! recompute zwz = f/depth at every time step for (.NOT.ln_linssh) as the water colomn height changes260 ! ! recompute zwz = f/depth at every time step for (.NOT.ln_linssh) as the water colomn height changes 254 261 ! 255 262 ! !* 2D Coriolis trends … … 258 265 ! 259 266 CALL dyn_cor_2d( ht(:,:), hu(:,:,Kmm), hv(:,:,Kmm), puu_b(:,:,Kmm), pvv_b(:,:,Kmm), zhU, zhV, & ! <<== in 260 & zu_trd, zv_trd ) ! ==>> out 261 ! 262 IF( .NOT.ln_linssh ) THEN !* surface pressure gradient (variable volume only) 263 ! 264 IF( ln_wd_il ) THEN ! W/D : limiter applied to spgspg 265 CALL wad_spg( pssh(:,:,Kmm), zcpx, zcpy ) ! Calculating W/D gravity filters, zcpx and zcpy 266 DO_2D( 0, 0, 0, 0 ) ! SPG with the application of W/D gravity filters 267 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( pssh(ji+1,jj ,Kmm) - pssh(ji ,jj ,Kmm) ) & 268 & * r1_e1u(ji,jj) * zcpx(ji,jj) * wdrampu(ji,jj) !jth 269 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( pssh(ji ,jj+1,Kmm) - pssh(ji ,jj ,Kmm) ) & 270 & * r1_e2v(ji,jj) * zcpy(ji,jj) * wdrampv(ji,jj) !jth 271 END_2D 272 ELSE ! now suface pressure gradient 273 DO_2D( 0, 0, 0, 0 ) 274 zu_trd(ji,jj) = zu_trd(ji,jj) - grav * ( pssh(ji+1,jj ,Kmm) - pssh(ji ,jj ,Kmm) ) * r1_e1u(ji,jj) 275 zv_trd(ji,jj) = zv_trd(ji,jj) - grav * ( pssh(ji ,jj+1,Kmm) - pssh(ji ,jj ,Kmm) ) * r1_e2v(ji,jj) 276 END_2D 277 ENDIF 278 ! 279 ENDIF 267 & zu_trd, zv_trd ) ! ==>> out 280 268 ! 281 269 DO_2D( 0, 0, 0, 0 ) ! Remove coriolis term (and possibly spg) from barotropic trend … … 287 275 ! ! ----------------------------------------------------------- ! 288 276 CALL dyn_drg_init( Kbb, Kmm, puu, pvv, puu_b ,pvv_b, zu_frc, zv_frc, zCdU_u, zCdU_v ) ! also provide the barotropic drag coefficients 277 ! 289 278 ! != Add atmospheric pressure forcing =! 290 279 ! ! ---------------------------------- ! … … 330 319 ELSE ! CENTRED integration: use kt-1/2 + kt+1/2 fluxes (NOW) 331 320 zztmp = r1_rho0 * r1_2 332 zssh_frc(:,:) = zztmp * ( emp(:,:) + emp_b(:,:)&333 & - rnf(:,:) - rnf_b(:,:)&334 & + fwfisf_cav(:,:) + fwfisf_cav_b(:,:)&335 & + fwfisf_par(:,:) + fwfisf_par_b(:,:))321 zssh_frc(:,:) = zztmp * ( emp(:,:) + emp_b(:,:) & 322 & - rnf(:,:) - rnf_b(:,:) & 323 & + fwfisf_cav(:,:) + fwfisf_cav_b(:,:) & 324 & + fwfisf_par(:,:) + fwfisf_par_b(:,:) ) 336 325 ENDIF 337 326 ! != Add Stokes drift divergence =! (if exist) -
NEMO/branches/2020/dev_r13723_KERNEL-01_Amy_Mike_newHPGschemes_KERNEL-08_merge/src/OCE/ISF/isf_oce.F90
r13558 r13935 1 1 MODULE isf_oce 2 2 !!====================================================================== 3 !! *** MODULE sbcisf***4 !! Surface module : compute iceshelf melt and heat flux3 !! *** MODULE isf_oce *** 4 !! Ice shelves : ice shelves variables defined in memory 5 5 !!====================================================================== 6 6 !! History : 3.2 ! 2011-02 (C.Harris ) Original code isf cav … … 146 146 END SUBROUTINE isf_alloc_par 147 147 148 148 149 SUBROUTINE isf_alloc_cav() 149 150 !!--------------------------------------------------------------------- … … 173 174 END SUBROUTINE isf_alloc_cav 174 175 176 175 177 SUBROUTINE isf_alloc_cpl() 176 178 !!--------------------------------------------------------------------- … … 184 186 ierr = 0 185 187 ! 186 ALLOCATE( risfcpl_ssh(jpi,jpj) , risfcpl_tsc(jpi,jpj,jpk,jpts), risfcpl_vol(jpi,jpj,jpk), STAT=ialloc )187 ierr = ierr + ialloc 188 ! 189 risfcpl_tsc(:,:,:,:) = 0. 0 ; risfcpl_vol(:,:,:) = 0.0 ; risfcpl_ssh(:,:) = 0.0190 191 IF ( ln_isfcpl_cons ) THEN192 ALLOCATE( risfcpl_cons_tsc(jpi,jpj,jpk,jpts) , risfcpl_cons_vol(jpi,jpj,jpk) , risfcpl_cons_ssh(jpi,jpj), STAT=ialloc )188 ALLOCATE( risfcpl_ssh(jpi,jpj) , risfcpl_tsc(jpi,jpj,jpk,jpts) , risfcpl_vol(jpi,jpj,jpk) , STAT=ialloc ) 189 ierr = ierr + ialloc 190 ! 191 risfcpl_tsc(:,:,:,:) = 0._wp ; risfcpl_vol(:,:,:) = 0._wp ; risfcpl_ssh(:,:) = 0._wp 192 193 IF ( ln_isfcpl_cons ) THEN 194 ALLOCATE( risfcpl_cons_tsc(jpi,jpj,jpk,jpts) , risfcpl_cons_vol(jpi,jpj,jpk) , risfcpl_cons_ssh(jpi,jpj) , STAT=ialloc ) 193 195 ierr = ierr + ialloc 194 196 ! 195 risfcpl_cons_tsc(:,:,:,:) = 0. 0 ; risfcpl_cons_vol(:,:,:) = 0.0 ; risfcpl_cons_ssh(:,:) = 0.0197 risfcpl_cons_tsc(:,:,:,:) = 0._wp ; risfcpl_cons_vol(:,:,:) = 0._wp ; risfcpl_cons_ssh(:,:) = 0._wp 196 198 ! 197 199 END IF … … 202 204 END SUBROUTINE isf_alloc_cpl 203 205 206 204 207 SUBROUTINE isf_dealloc_cpl() 205 208 !!--------------------------------------------------------------------- … … 213 216 ierr = 0 214 217 ! 215 DEALLOCATE( risfcpl_ssh , risfcpl_tsc, risfcpl_vol, STAT=ialloc )218 DEALLOCATE( risfcpl_ssh , risfcpl_tsc , risfcpl_vol , STAT=ialloc ) 216 219 ierr = ierr + ialloc 217 220 ! … … 221 224 END SUBROUTINE isf_dealloc_cpl 222 225 226 223 227 SUBROUTINE isf_alloc() 224 228 !!--------------------------------------------------------------------- … … 233 237 ierr = 0 ! set to zero if no array to be allocated 234 238 ! 235 ALLOCATE( fwfisf_par(jpi,jpj) , fwfisf_par_b(jpi,jpj),&236 & fwfisf_cav(jpi,jpj) , fwfisf_cav_b(jpi,jpj),&237 & fwfisf_oasis(jpi,jpj),STAT=ialloc )238 ierr = ierr + ialloc 239 ! 240 ALLOCATE( risf_par_tsc(jpi,jpj,jpts), risf_par_tsc_b(jpi,jpj,jpts), STAT=ialloc )241 ierr = ierr + ialloc 242 ! 243 ALLOCATE( risf_cav_tsc(jpi,jpj,jpts), risf_cav_tsc_b(jpi,jpj,jpts), STAT=ialloc )244 ierr = ierr + ialloc 245 ! 246 ALLOCATE( risfload(jpi,jpj), STAT=ialloc)247 ierr = ierr + ialloc 248 ! 249 ALLOCATE( mskisf_cav(jpi,jpj) , STAT=ialloc)239 ALLOCATE( fwfisf_par (jpi,jpj) , fwfisf_par_b(jpi,jpj) , & 240 & fwfisf_cav (jpi,jpj) , fwfisf_cav_b(jpi,jpj) , & 241 & fwfisf_oasis(jpi,jpj) , STAT=ialloc ) 242 ierr = ierr + ialloc 243 ! 244 ALLOCATE( risf_par_tsc(jpi,jpj,jpts) , risf_par_tsc_b(jpi,jpj,jpts) , STAT=ialloc ) 245 ierr = ierr + ialloc 246 ! 247 ALLOCATE( risf_cav_tsc(jpi,jpj,jpts) , risf_cav_tsc_b(jpi,jpj,jpts) , STAT=ialloc ) 248 ierr = ierr + ialloc 249 ! 250 ALLOCATE( risfload(jpi,jpj) , STAT=ialloc ) 251 ierr = ierr + ialloc 252 ! 253 ALLOCATE( mskisf_cav(jpi,jpj) , STAT=ialloc ) 250 254 ierr = ierr + ialloc 251 255 ! … … 254 258 ! 255 259 ! initalisation of fwf and tsc array to 0 256 risfload(:,:) = 0.0_wp 257 fwfisf_oasis(:,:) = 0.0_wp 258 fwfisf_par(:,:) = 0.0_wp ; fwfisf_par_b(:,:) = 0.0_wp 259 fwfisf_cav(:,:) = 0.0_wp ; fwfisf_cav_b(:,:) = 0.0_wp 260 risf_cav_tsc(:,:,:) = 0.0_wp ; risf_cav_tsc_b(:,:,:) = 0.0_wp 261 risf_par_tsc(:,:,:) = 0.0_wp ; risf_par_tsc_b(:,:,:) = 0.0_wp 262 ! 263 260 risfload (:,:) = 0._wp 261 fwfisf_oasis(:,:) = 0._wp 262 fwfisf_par (:,:) = 0._wp ; fwfisf_par_b (:,:) = 0._wp 263 fwfisf_cav (:,:) = 0._wp ; fwfisf_cav_b (:,:) = 0._wp 264 risf_cav_tsc(:,:,:) = 0._wp ; risf_cav_tsc_b(:,:,:) = 0._wp 265 risf_par_tsc(:,:,:) = 0._wp ; risf_par_tsc_b(:,:,:) = 0._wp 266 ! 264 267 END SUBROUTINE isf_alloc 265 268 269 !!====================================================================== 266 270 END MODULE isf_oce -
NEMO/branches/2020/dev_r13723_KERNEL-01_Amy_Mike_newHPGschemes_KERNEL-08_merge/src/OCE/ISF/isfload.F90
r13295 r13935 2 2 !!====================================================================== 3 3 !! *** MODULE isfload *** 4 !! isfload module :compute ice shelf load (needed for the hpg)4 !! Ice Shelves : compute ice shelf load (needed for the hpg) 5 5 !!====================================================================== 6 6 !! History : 4.1 ! 2019-09 (P. Mathiot) original code … … 8 8 9 9 !!---------------------------------------------------------------------- 10 !! isf load : compute ice shelf load10 !! isf_load : compute ice shelf load 11 11 !!---------------------------------------------------------------------- 12 12 … … 23 23 PRIVATE 24 24 25 PUBLIC isf_load 25 PUBLIC isf_load ! called by isfstp.F90 26 ! 26 27 !! * Substitutions 27 28 # include "do_loop_substitute.h90" 28 29 # include "domzgr_substitute.h90" 29 30 !!---------------------------------------------------------------------- 31 !! NEMO/OCE 4.0 , NEMO Consortium (2018) 32 !! $Id: sbcisf.F90 10536 2019-01-16 19:21:09Z mathiot $ 33 !! Software governed by the CeCILL license (see ./LICENSE) 34 !!---------------------------------------------------------------------- 30 35 CONTAINS 31 36 … … 37 42 !! 38 43 !!-------------------------------------------------------------------- 39 !!-------------------------- OUT ------------------------------------- 40 REAL(wp), DIMENSION(jpi,jpj), INTENT( out) :: pisfload 41 !!-------------------------- IN ------------------------------------- 42 INTEGER, INTENT(in) :: Kmm ! ocean time level index 44 INTEGER, INTENT(in ) :: Kmm ! ocean time level index 45 REAL(wp), DIMENSION(jpi,jpj), INTENT( out) :: pisfload ! ice shelf load 43 46 !!---------------------------------------------------------------------- 44 47 ! … … 46 49 ! the smaller the residual flow is, the better it is. 47 50 ! 48 ! ice shelf cavity51 ! type of ice shelf cavity 49 52 SELECT CASE ( cn_isfload ) 50 53 CASE ( 'uniform' ) … … 56 59 END SUBROUTINE isf_load 57 60 58 SUBROUTINE isf_load_uniform( Kmm, pisfload ) 61 62 SUBROUTINE isf_load_uniform( Kmm, pload ) 59 63 !!-------------------------------------------------------------------- 60 64 !! *** SUBROUTINE isf_load *** … … 67 71 !! 68 72 !!-------------------------------------------------------------------- 69 !!-------------------------- OUT ------------------------------------- 70 REAL(wp), DIMENSION(jpi,jpj), INTENT( out) :: pisfload 71 !!-------------------------- IN ------------------------------------- 72 INTEGER, INTENT(in) :: Kmm ! ocean time level index 73 !!-------------------------------------------------------------------- 73 INTEGER, INTENT(in ) :: Kmm ! ocean time level index 74 REAL(wp), DIMENSION(jpi,jpj), INTENT( out) :: pload ! ice shelf load 75 ! 74 76 INTEGER :: ji, jj, jk 75 77 INTEGER :: ikt 76 REAL(wp) :: znad !77 78 REAL(wp), DIMENSION(jpi,jpj) :: zrhdtop_isf ! water density displaced by the ice shelf (at the interface) 78 79 REAL(wp), DIMENSION(jpi,jpj,jpts) :: zts_top ! water properties displaced by the ice shelf 79 80 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zrhd ! water density displaced by the ice shelf 80 81 !!---------------------------------------------------------------------- 81 !82 znad = 1._wp !- To use density and not density anomaly83 82 ! 84 83 ! !- assume water displaced by the ice shelf is at T=rn_isfload_T and S=rn_isfload_S (rude) … … 87 86 DO jk = 1, jpk !- compute density of the water displaced by the ice shelf 88 87 CALL eos( zts_top(:,:,:), gdept(:,:,jk,Kmm), zrhd(:,:,jk) ) 88 !!st ==>> CALL eos( zts_top(:,:,:), gdept_0(:,:,jk), zrhd(:,:,jk) ) 89 89 END DO 90 90 ! … … 93 93 ! 94 94 ! !- Surface value + ice shelf gradient 95 p isfload(:,:) = 0._wp! compute pressure due to ice shelf load95 pload(:,:) = 0._wp ! compute pressure due to ice shelf load 96 96 DO_2D( 1, 1, 1, 1 ) 97 97 ikt = mikt(ji,jj) 98 98 ! 99 99 IF ( ikt > 1 ) THEN 100 ! ! top layer of the ice shelf 101 pload(ji,jj) = pload(ji,jj) & 102 & + zrhd (ji,jj,1) * e3w(ji,jj,1,Kmm) 100 103 ! 101 ! top layer of the ice shelf 102 pisfload(ji,jj) = pisfload(ji,jj) + (znad + zrhd(ji,jj,1) ) & 103 & * e3w(ji,jj,1,Kmm) 104 ! 105 ! core layers of the ice shelf 106 DO jk = 2, ikt-1 107 pisfload(ji,jj) = pisfload(ji,jj) + (2._wp * znad + zrhd(ji,jj,jk-1) + zrhd(ji,jj,jk)) & 108 & * e3w(ji,jj,jk,Kmm) 104 DO jk = 2, ikt-1 ! core layers of the ice shelf 105 pload(ji,jj) = pload(ji,jj) + (zrhd(ji,jj,jk-1) + zrhd(ji,jj,jk)) & 106 & * e3w(ji,jj,jk,Kmm) 109 107 END DO 110 ! 111 ! deepest part of the ice shelf (between deepest T point and ice/ocean interface112 pisfload(ji,jj) = pisfload(ji,jj) + (2._wp * znad + zrhdtop_isf(ji,jj) + zrhd(ji,jj,ikt-1)) &113 & * ( risfdep(ji,jj) - gdept(ji,jj,ikt-1,Kmm) )108 ! ! deepest part of the ice shelf (between deepest T point and ice/ocean interface 109 pload(ji,jj) = pload(ji,jj) + ( zrhdtop_isf(ji,jj) + zrhd(ji,jj,ikt-1) ) & 110 & * ( risfdep(ji,jj) - gdept(ji,jj,ikt-1,Kmm) ) 111 !!st ==>> & * ( risfdep(ji,jj) - gdept_0(ji,jj,ikt-1) ) 114 112 ! 115 113 END IF … … 117 115 ! 118 116 END SUBROUTINE isf_load_uniform 119 117 118 !!====================================================================== 120 119 END MODULE isfload -
NEMO/branches/2020/dev_r13723_KERNEL-01_Amy_Mike_newHPGschemes_KERNEL-08_merge/src/OCE/ISF/isfstp.F90
r13237 r13935 2 2 !!====================================================================== 3 3 !! *** MODULE isfstp *** 4 !! Surface module: compute iceshelf load, melt and heat flux4 !! Ice Shelves : compute iceshelf load, melt and heat flux 5 5 !!====================================================================== 6 6 !! History : 3.2 ! 2011-02 (C.Harris ) Original code isf cav … … 42 42 !! Software governed by the CeCILL license (see ./LICENSE) 43 43 !!---------------------------------------------------------------------- 44 45 44 CONTAINS 46 45 47 SUBROUTINE isf_stp( kt, Kmm )46 SUBROUTINE isf_stp( kt, Kmm ) 48 47 !!--------------------------------------------------------------------- 49 48 !! *** ROUTINE isf_stp *** … … 58 57 !! - compute fluxes 59 58 !! - write restart variables 60 !! 61 !!---------------------------------------------------------------------- 62 INTEGER, INTENT(in) :: kt ! ocean time step 63 INTEGER, INTENT(in) :: Kmm ! ocean time level index 64 !!---------------------------------------------------------------------- 65 INTEGER :: jk ! loop index 66 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3t ! e3t 59 !!---------------------------------------------------------------------- 60 INTEGER, INTENT(in) :: kt ! ocean time step 61 INTEGER, INTENT(in) :: Kmm ! ocean time level index 62 ! 63 INTEGER :: jk ! loop index 64 #if defined key_qco 65 REAL(wp), DIMENSION(jpi,jpj,jpk) :: ze3t ! 3D workspace 66 #endif 67 67 !!--------------------------------------------------------------------- 68 68 ! … … 83 83 ! 1.2: compute misfkb, rhisf_tbl, rfrac (deepest level, thickness, fraction of deepest cell affected by tbl) 84 84 rhisf_tbl_cav(:,:) = rn_htbl * mskisf_cav(:,:) 85 #if defined key_qco 85 86 DO jk = 1, jpk 86 87 ze3t(:,:,jk) = e3t(:,:,jk,Kmm) 87 88 END DO 88 CALL isf_tbl_lvl(ht(:,:), ze3t, misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav) 89 CALL isf_tbl_lvl( ht(:,:), ze3t, misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav ) 90 #else 91 CALL isf_tbl_lvl( ht(:,:), e3t, misfkt_cav, misfkb_cav, rhisf_tbl_cav, rfrac_tbl_cav ) 92 #endif 89 93 ! 90 94 ! 1.3: compute ice shelf melt 91 CALL isf_cav( kt, Kmm, risf_cav_tsc, fwfisf_cav )95 CALL isf_cav( kt, Kmm, risf_cav_tsc, fwfisf_cav ) 92 96 ! 93 97 END IF … … 108 112 ! by simplicity, we assume the top level where param applied do not change with time (done in init part) 109 113 rhisf_tbl_par(:,:) = rhisf0_tbl_par(:,:) 114 #if defined key_qco 110 115 DO jk = 1, jpk 111 116 ze3t(:,:,jk) = e3t(:,:,jk,Kmm) 112 117 END DO 113 CALL isf_tbl_lvl(ht(:,:), ze3t, misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par) 118 CALL isf_tbl_lvl( ht(:,:), ze3t, misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par ) 119 #else 120 CALL isf_tbl_lvl( ht(:,:), e3t, misfkt_par, misfkb_par, rhisf_tbl_par, rfrac_tbl_par ) 121 #endif 114 122 ! 115 123 ! 2.3: compute ice shelf melt 116 CALL isf_par( kt, Kmm, risf_par_tsc, fwfisf_par )124 CALL isf_par( kt, Kmm, risf_par_tsc, fwfisf_par ) 117 125 ! 118 126 END IF … … 128 136 END SUBROUTINE isf_stp 129 137 130 SUBROUTINE isf_init(Kbb, Kmm, Kaa) 138 139 SUBROUTINE isf_init( Kbb, Kmm, Kaa ) 131 140 !!--------------------------------------------------------------------- 132 141 !! *** ROUTINE isfstp_init *** … … 142 151 !! - call cav/param/isfcpl init routine 143 152 !!---------------------------------------------------------------------- 144 INTEGER, INTENT(in) :: Kbb, Kmm, Kaa ! ocean time level indices 153 INTEGER, INTENT(in) :: Kbb, Kmm, Kaa ! ocean time level indices 154 !!---------------------------------------------------------------------- 145 155 ! 146 156 ! constrain: l_isfoasis need to be known 147 157 ! 148 ! Read namelist 149 CALL isf_nam() 150 ! 151 ! Allocate public array 152 CALL isf_alloc() 153 ! 154 ! check option compatibility 155 CALL isf_ctl() 156 ! 157 ! compute ice shelf load 158 IF ( ln_isfcav ) CALL isf_load( Kmm, risfload ) 158 CALL isf_nam() ! Read namelist 159 ! 160 CALL isf_alloc() ! Allocate public array 161 ! 162 CALL isf_ctl() ! check option compatibility 163 ! 164 IF( ln_isfcav ) CALL isf_load( Kmm, risfload ) ! compute ice shelf load 159 165 ! 160 166 ! terminate routine now if no ice shelf melt formulation specify 161 IF ( ln_isf ) THEN 162 ! 163 !--------------------------------------------------------------------------------------------------------------------- 164 ! initialisation melt in the cavity 165 IF ( ln_isfcav_mlt ) CALL isf_cav_init() 166 ! 167 !--------------------------------------------------------------------------------------------------------------------- 168 ! initialisation parametrised melt 169 IF ( ln_isfpar_mlt ) CALL isf_par_init() 170 ! 171 !--------------------------------------------------------------------------------------------------------------------- 172 ! initialisation ice sheet coupling 173 IF( ln_isfcpl ) CALL isfcpl_init(Kbb, Kmm, Kaa) 167 IF( ln_isf ) THEN 168 ! 169 IF( ln_isfcav_mlt ) CALL isf_cav_init() ! initialisation melt in the cavity 170 ! 171 IF( ln_isfpar_mlt ) CALL isf_par_init() ! initialisation parametrised melt 172 ! 173 IF( ln_isfcpl ) CALL isfcpl_init( Kbb, Kmm, Kaa ) ! initialisation ice sheet coupling 174 174 ! 175 175 END IF … … 177 177 END SUBROUTINE isf_init 178 178 179 179 180 SUBROUTINE isf_ctl() 180 181 !!--------------------------------------------------------------------- … … 283 284 END IF 284 285 END SUBROUTINE isf_ctl 285 ! 286 287 286 288 SUBROUTINE isf_nam 287 289 !!--------------------------------------------------------------------- -
NEMO/branches/2020/dev_r13723_KERNEL-01_Amy_Mike_newHPGschemes_KERNEL-08_merge/src/OCE/oce.F90
r13237 r13935 49 49 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ub2_i_b, vb2_i_b !: Half step time integrated fluxes 50 50 #endif 51 !52 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: spgu, spgv !: horizontal surface pressure gradient53 51 54 52 !! interpolated gradient (only used in zps case) … … 83 81 ! 84 82 ierr(:) = 0 85 ALLOCATE( uu (jpi,jpj,jpk,jpt) , vv (jpi,jpj,jpk,jpt) ,&86 & ww (jpi,jpj,jpk) , hdiv(jpi,jpj,jpk) ,&87 & ts (jpi,jpj,jpk,jpts,jpt) ,&88 & rab_b(jpi,jpj,jpk,jpts) , rab_n(jpi,jpj,jpk,jpts) ,&89 & rn2b (jpi,jpj,jpk) , rn2 (jpi,jpj,jpk) ,&90 & rhd (jpi,jpj,jpk) , rhop (jpi,jpj,jpk) 83 ALLOCATE( uu (jpi,jpj,jpk,jpt) , vv (jpi,jpj,jpk,jpt) , & 84 & ww (jpi,jpj,jpk) , hdiv(jpi,jpj,jpk) , & 85 & ts (jpi,jpj,jpk,jpts,jpt) , & 86 & rab_b(jpi,jpj,jpk,jpts) , rab_n(jpi,jpj,jpk,jpts) , & 87 & rn2b (jpi,jpj,jpk) , rn2 (jpi,jpj,jpk) , & 88 & rhd (jpi,jpj,jpk) , rhop (jpi,jpj,jpk) , STAT=ierr(1) ) 91 89 ! 92 ALLOCATE( ssh(jpi,jpj,jpt) , uu_b(jpi,jpj,jpt) , vv_b(jpi,jpj,jpt) , & 93 & spgu (jpi,jpj) , spgv(jpi,jpj) , & 94 & gtsu(jpi,jpj,jpts), gtsv(jpi,jpj,jpts) , & 95 & gru(jpi,jpj) , grv(jpi,jpj) , & 96 & gtui(jpi,jpj,jpts), gtvi(jpi,jpj,jpts) , & 97 & grui(jpi,jpj) , grvi(jpi,jpj) , & 98 & riceload(jpi,jpj) , STAT=ierr(2) ) 90 ALLOCATE( ssh (jpi,jpj,jpt) , uu_b(jpi,jpj,jpt) , vv_b(jpi,jpj,jpt) , & 91 & gtsu(jpi,jpj,jpts) , gtsv(jpi,jpj,jpts) , & 92 & gru (jpi,jpj) , grv (jpi,jpj) , & 93 & gtui(jpi,jpj,jpts) , gtvi(jpi,jpj,jpts) , & 94 & grui(jpi,jpj) , grvi(jpi,jpj) , & 95 & riceload(jpi,jpj) , STAT=ierr(2) ) 99 96 ! 100 97 ALLOCATE( fraqsr_1lev(jpi,jpj) , STAT=ierr(3) )
Note: See TracChangeset
for help on using the changeset viewer.
