Changeset 14064
- Timestamp:
- 2020-12-03T18:01:12+01:00 (3 years ago)
- Location:
- NEMO/trunk/src
- Files:
-
- 9 edited
-
NST/agrif_oce_update.F90 (modified) (6 diffs)
-
OCE/DYN/dynhpg.F90 (modified) (33 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/trunk/src/NST/agrif_oce_update.F90
r14053 r14064 915 915 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres 916 916 LOGICAL , INTENT(in ) :: before 917 !! 918 REAL(wp), DIMENSION(jpi,jpj) :: zpgu ! 2D workspace 917 919 !! 918 920 INTEGER :: ji, jj, jk … … 934 936 ! 935 937 ! Update "now" 3d velocities: 936 spgu(ji,jj) = 0._wp938 zpgu(ji,jj) = 0._wp 937 939 DO jk=1,jpkm1 938 spgu(ji,jj) = spgu(ji,jj) + e3u(ji,jj,jk,Kmm_a) * uu(ji,jj,jk,Kmm_a)940 zpgu(ji,jj) = zpgu(ji,jj) + e3u(ji,jj,jk,Kmm_a) * uu(ji,jj,jk,Kmm_a) 939 941 END DO 940 942 ! 941 zcorr = (tabres(ji,jj) - spgu(ji,jj)) * r1_hu(ji,jj,Kmm_a)943 zcorr = (tabres(ji,jj) - zpgu(ji,jj)) * r1_hu(ji,jj,Kmm_a) 942 944 DO jk=1,jpkm1 943 945 uu(ji,jj,jk,Kmm_a) = uu(ji,jj,jk,Kmm_a) + zcorr * umask(ji,jj,jk) … … 954 956 ! 955 957 ! Correct "before" velocities to hold correct bt component: 956 spgu(ji,jj) = 0.e0958 zpgu(ji,jj) = 0.e0 957 959 DO jk=1,jpkm1 958 spgu(ji,jj) = spgu(ji,jj) + e3u(ji,jj,jk,Kbb_a) * uu(ji,jj,jk,Kbb_a)960 zpgu(ji,jj) = zpgu(ji,jj) + e3u(ji,jj,jk,Kbb_a) * uu(ji,jj,jk,Kbb_a) 959 961 END DO 960 962 ! 961 zcorr = uu_b(ji,jj,Kbb_a) - spgu(ji,jj) * r1_hu(ji,jj,Kbb_a)963 zcorr = uu_b(ji,jj,Kbb_a) - zpgu(ji,jj) * r1_hu(ji,jj,Kbb_a) 962 964 DO jk=1,jpkm1 963 965 uu(ji,jj,jk,Kbb_a) = uu(ji,jj,jk,Kbb_a) + zcorr * umask(ji,jj,jk) … … 982 984 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: tabres 983 985 LOGICAL , INTENT(in ) :: before 986 ! 987 REAL(wp), DIMENSION(jpi,jpj) :: zpgv ! 2D workspace 984 988 ! 985 989 INTEGER :: ji, jj, jk … … 1000 1004 ! 1001 1005 ! Update "now" 3d velocities: 1002 spgv(ji,jj) = 0.e01006 zpgv(ji,jj) = 0.e0 1003 1007 DO jk=1,jpkm1 1004 spgv(ji,jj) = spgv(ji,jj) + e3v(ji,jj,jk,Kmm_a) * vv(ji,jj,jk,Kmm_a)1008 zpgv(ji,jj) = zpgv(ji,jj) + e3v(ji,jj,jk,Kmm_a) * vv(ji,jj,jk,Kmm_a) 1005 1009 END DO 1006 1010 ! 1007 zcorr = (tabres(ji,jj) - spgv(ji,jj)) * r1_hv(ji,jj,Kmm_a)1011 zcorr = (tabres(ji,jj) - zpgv(ji,jj)) * r1_hv(ji,jj,Kmm_a) 1008 1012 DO jk=1,jpkm1 1009 1013 vv(ji,jj,jk,Kmm_a) = vv(ji,jj,jk,Kmm_a) + zcorr * vmask(ji,jj,jk) … … 1020 1024 ! 1021 1025 ! Correct "before" velocities to hold correct bt component: 1022 spgv(ji,jj) = 0.e01026 zpgv(ji,jj) = 0.e0 1023 1027 DO jk=1,jpkm1 1024 spgv(ji,jj) = spgv(ji,jj) + e3v(ji,jj,jk,Kbb_a) * vv(ji,jj,jk,Kbb_a)1028 zpgv(ji,jj) = zpgv(ji,jj) + e3v(ji,jj,jk,Kbb_a) * vv(ji,jj,jk,Kbb_a) 1025 1029 END DO 1026 1030 ! 1027 zcorr = vv_b(ji,jj,Kbb_a) - spgv(ji,jj) * r1_hv(ji,jj,Kbb_a)1031 zcorr = vv_b(ji,jj,Kbb_a) - zpgv(ji,jj) * r1_hv(ji,jj,Kbb_a) 1028 1032 DO jk=1,jpkm1 1029 1033 vv(ji,jj,jk,Kbb_a) = vv(ji,jj,jk,Kbb_a) + zcorr * vmask(ji,jj,jk) -
NEMO/trunk/src/OCE/DYN/dynhpg.F90
r14060 r14064 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 ! … … 411 405 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) :: puu, pvv ! ocean velocities and RHS of momentum equation 412 406 !! 413 INTEGER :: ji, jj, jk, jii, jjj 414 REAL(wp) :: zcoef0, zuap, zvap, z nad, ztmp ! temporaryscalars415 LOGICAL :: ll_tmp1, ll_tmp2 407 INTEGER :: ji, jj, jk, jii, jjj ! dummy loop indices 408 REAL(wp) :: zcoef0, zuap, zvap, ztmp ! local scalars 409 LOGICAL :: ll_tmp1, ll_tmp2 ! local logical variables 416 410 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zhpi, zhpj 417 411 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: zcpx, zcpy !W/D pressure filter … … 427 421 ! 428 422 zcoef0 = - grav * 0.5_wp 429 IF ( ln_linssh ) THEN ; znad = 0._wp ! Fixed volume: density anomaly430 ELSE ; znad = 1._wp ! Variable volume: density431 ENDIF432 423 ! 433 424 IF( ln_wd_il ) THEN … … 471 462 CALL lbc_lnk_multi( 'dynhpg', zcpx, 'U', 1.0_wp, zcpy, 'V', 1.0_wp ) 472 463 END IF 473 474 ! Surface value 475 DO_2D( 0, 0, 0, 0 ) 476 ! hydrostatic pressure gradient along s-surfaces 477 zhpi(ji,jj,1) = & 478 & zcoef0 * ( e3w(ji+1,jj ,1,Kmm) * ( znad + rhd(ji+1,jj ,1) ) & 479 & - e3w(ji ,jj ,1,Kmm) * ( znad + rhd(ji ,jj ,1) ) ) & 480 & * r1_e1u(ji,jj) 481 zhpj(ji,jj,1) = & 482 & zcoef0 * ( e3w(ji ,jj+1,1,Kmm) * ( znad + rhd(ji ,jj+1,1) ) & 483 & - e3w(ji ,jj ,1,Kmm) * ( znad + rhd(ji ,jj ,1) ) ) & 484 & * r1_e2v(ji,jj) 485 ! s-coordinate pressure gradient correction 486 zuap = -zcoef0 * ( rhd (ji+1,jj,1) + rhd (ji,jj,1) + 2._wp * znad ) & 464 ! 465 DO_2D( 0, 0, 0, 0 ) ! Surface value 466 ! ! hydrostatic pressure gradient along s-surfaces 467 zhpi(ji,jj,1) = zcoef0 * r1_e1u(ji,jj) & 468 & * ( e3w(ji+1,jj ,1,Kmm) * rhd(ji+1,jj ,1) & 469 & - e3w(ji ,jj ,1,Kmm) * rhd(ji ,jj ,1) ) 470 zhpj(ji,jj,1) = zcoef0 * r1_e2v(ji,jj) & 471 & * ( e3w(ji ,jj+1,1,Kmm) * rhd(ji ,jj+1,1) & 472 & - e3w(ji ,jj ,1,Kmm) * rhd(ji ,jj ,1) ) 473 ! ! s-coordinate pressure gradient correction 474 zuap = -zcoef0 * ( rhd (ji+1,jj,1) + rhd (ji,jj,1) ) & 487 475 & * ( gde3w(ji+1,jj,1) - gde3w(ji,jj,1) ) * r1_e1u(ji,jj) 488 zvap = -zcoef0 * ( rhd (ji,jj+1,1) + rhd (ji,jj,1) + 2._wp * znad) &476 zvap = -zcoef0 * ( rhd (ji,jj+1,1) + rhd (ji,jj,1) ) & 489 477 & * ( gde3w(ji,jj+1,1) - gde3w(ji,jj,1) ) * r1_e2v(ji,jj) 490 478 ! … … 495 483 zvap = zvap * zcpy(ji,jj) 496 484 ENDIF 497 ! 498 ! add to the general momentum trend 485 ! ! add to the general momentum trend 499 486 puu(ji,jj,1,Krhs) = puu(ji,jj,1,Krhs) + zhpi(ji,jj,1) + zuap 500 487 pvv(ji,jj,1,Krhs) = pvv(ji,jj,1,Krhs) + zhpj(ji,jj,1) + zvap 501 488 END_2D 502 503 ! interior value (2=<jk=<jpkm1) 504 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 505 ! hydrostatic pressure gradient along s-surfaces 506 zhpi(ji,jj,jk) = zhpi(ji,jj,jk-1) + zcoef0 * r1_e1u(ji,jj) & 507 & * ( e3w(ji+1,jj,jk,Kmm) * ( rhd(ji+1,jj,jk) + rhd(ji+1,jj,jk-1) + 2*znad ) & 508 & - e3w(ji ,jj,jk,Kmm) * ( rhd(ji ,jj,jk) + rhd(ji ,jj,jk-1) + 2*znad ) ) 509 zhpj(ji,jj,jk) = zhpj(ji,jj,jk-1) + zcoef0 * r1_e2v(ji,jj) & 510 & * ( e3w(ji,jj+1,jk,Kmm) * ( rhd(ji,jj+1,jk) + rhd(ji,jj+1,jk-1) + 2*znad ) & 511 & - e3w(ji,jj ,jk,Kmm) * ( rhd(ji,jj, jk) + rhd(ji,jj ,jk-1) + 2*znad ) ) 512 ! s-coordinate pressure gradient correction 513 zuap = -zcoef0 * ( rhd (ji+1,jj ,jk) + rhd (ji,jj,jk) + 2._wp * znad ) & 489 ! 490 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) ! interior value (2=<jk=<jpkm1) 491 ! ! hydrostatic pressure gradient along s-surfaces 492 zhpi(ji,jj,jk) = zhpi(ji,jj,jk-1) + zcoef0 * r1_e1u(ji,jj) & 493 & * ( e3w(ji+1,jj,jk,Kmm) * ( rhd(ji+1,jj,jk) + rhd(ji+1,jj,jk-1) ) & 494 & - e3w(ji ,jj,jk,Kmm) * ( rhd(ji ,jj,jk) + rhd(ji ,jj,jk-1) ) ) 495 zhpj(ji,jj,jk) = zhpj(ji,jj,jk-1) + zcoef0 * r1_e2v(ji,jj) & 496 & * ( e3w(ji,jj+1,jk,Kmm) * ( rhd(ji,jj+1,jk) + rhd(ji,jj+1,jk-1) ) & 497 & - e3w(ji,jj ,jk,Kmm) * ( rhd(ji,jj, jk) + rhd(ji,jj ,jk-1) ) ) 498 ! ! s-coordinate pressure gradient correction 499 zuap = -zcoef0 * ( rhd (ji+1,jj ,jk) + rhd (ji,jj,jk) ) & 514 500 & * ( gde3w(ji+1,jj ,jk) - gde3w(ji,jj,jk) ) * r1_e1u(ji,jj) 515 zvap = -zcoef0 * ( rhd (ji ,jj+1,jk) + rhd (ji,jj,jk) + 2._wp * znad )&501 zvap = -zcoef0 * ( rhd (ji ,jj+1,jk) + rhd (ji,jj,jk) ) & 516 502 & * ( gde3w(ji ,jj+1,jk) - gde3w(ji,jj,jk) ) * r1_e2v(ji,jj) 517 503 ! … … 556 542 REAL(wp), DIMENSION(jpi,jpj,jpk,jpt), INTENT(inout) :: puu, pvv ! ocean velocities and RHS of momentum equation 557 543 !! 558 INTEGER :: ji, jj, jk, ikt, iktp1i, iktp1j ! dummy loop indices 559 REAL(wp) :: zcoef0, zuap, zvap, znad ! temporary scalars 544 INTEGER :: ji, jj, jk ! dummy loop indices 545 INTEGER :: ikt , ikti1, iktj1 ! local integer 546 REAL(wp) :: ze3w, ze3wi1, ze3wj1 ! local scalars 547 REAL(wp) :: zcoef0, zuap, zvap ! - - 560 548 REAL(wp), DIMENSION(jpi,jpj,jpk ) :: zhpi, zhpj 561 549 REAL(wp), DIMENSION(jpi,jpj,jpts) :: zts_top … … 564 552 ! 565 553 zcoef0 = - grav * 0.5_wp ! Local constant initialization 566 !567 znad=1._wp ! To use density and not density anomaly568 554 ! 569 555 ! ! iniitialised to 0. zhpi zhpi … … 581 567 CALL eos( zts_top, risfdep, zrhdtop_oce ) 582 568 583 !================================================================================== 584 !===== Compute surface value ===================================================== 585 !================================================================================== 569 ! !===========================! 570 ! !===== surface value =====! 571 ! !===========================! 586 572 DO_2D( 0, 0, 0, 0 ) 587 ikt = mikt(ji,jj) 588 iktp1i = mikt(ji+1,jj) 589 iktp1j = mikt(ji,jj+1) 590 ! hydrostatic pressure gradient along s-surfaces and ice shelf pressure 591 ! we assume ISF is in isostatic equilibrium 592 zhpi(ji,jj,1) = zcoef0 / e1u(ji,jj) * ( 0.5_wp * e3w(ji+1,jj,iktp1i,Kmm) & 593 & * ( 2._wp * znad + rhd(ji+1,jj,iktp1i) + zrhdtop_oce(ji+1,jj) ) & 594 & - 0.5_wp * e3w(ji,jj,ikt,Kmm) & 595 & * ( 2._wp * znad + rhd(ji,jj,ikt) + zrhdtop_oce(ji,jj) ) & 596 & + ( risfload(ji+1,jj) - risfload(ji,jj)) ) 597 zhpj(ji,jj,1) = zcoef0 / e2v(ji,jj) * ( 0.5_wp * e3w(ji,jj+1,iktp1j,Kmm) & 598 & * ( 2._wp * znad + rhd(ji,jj+1,iktp1j) + zrhdtop_oce(ji,jj+1) ) & 599 & - 0.5_wp * e3w(ji,jj,ikt,Kmm) & 600 & * ( 2._wp * znad + rhd(ji,jj,ikt) + zrhdtop_oce(ji,jj) ) & 601 & + ( risfload(ji,jj+1) - risfload(ji,jj)) ) 602 ! s-coordinate pressure gradient correction (=0 if z coordinate) 603 zuap = -zcoef0 * ( rhd (ji+1,jj,1) + rhd (ji,jj,1) + 2._wp * znad ) & 573 ikt = mikt(ji ,jj ) ; ze3w = e3w(ji ,jj ,ikt ,Kmm) 574 ikti1 = mikt(ji+1,jj ) ; ze3wi1 = e3w(ji+1,jj ,ikti1,Kmm) 575 iktj1 = mikt(ji ,jj+1) ; ze3wj1 = e3w(ji ,jj+1,iktj1,Kmm) 576 ! ! hydrostatic pressure gradient along s-surfaces and ice shelf pressure 577 ! ! we assume ISF is in isostatic equilibrium 578 zhpi(ji,jj,1) = zcoef0 * r1_e1u(ji,jj) * ( risfload(ji+1,jj) - risfload(ji,jj) & 579 & + 0.5_wp * ( ze3wi1 * ( rhd(ji+1,jj,ikti1) + zrhdtop_oce(ji+1,jj) ) & 580 & - ze3w * ( rhd(ji ,jj,ikt ) + zrhdtop_oce(ji ,jj) ) ) ) 581 zhpj(ji,jj,1) = zcoef0 * r1_e2v(ji,jj) * ( risfload(ji,jj+1) - risfload(ji,jj) & 582 & + 0.5_wp * ( ze3wj1 * ( rhd(ji,jj+1,iktj1) + zrhdtop_oce(ji,jj+1) ) & 583 & - ze3w * ( rhd(ji,jj ,ikt ) + zrhdtop_oce(ji,jj ) ) ) ) 584 ! ! s-coordinate pressure gradient correction (=0 if z coordinate) 585 zuap = -zcoef0 * ( rhd (ji+1,jj,1) + rhd (ji,jj,1) ) & 604 586 & * ( gde3w(ji+1,jj,1) - gde3w(ji,jj,1) ) * r1_e1u(ji,jj) 605 zvap = -zcoef0 * ( rhd (ji,jj+1,1) + rhd (ji,jj,1) + 2._wp * znad) &587 zvap = -zcoef0 * ( rhd (ji,jj+1,1) + rhd (ji,jj,1) ) & 606 588 & * ( gde3w(ji,jj+1,1) - gde3w(ji,jj,1) ) * r1_e2v(ji,jj) 607 ! add to the general momentum trend589 ! ! add to the general momentum trend 608 590 puu(ji,jj,1,Krhs) = puu(ji,jj,1,Krhs) + (zhpi(ji,jj,1) + zuap) * umask(ji,jj,1) 609 591 pvv(ji,jj,1,Krhs) = pvv(ji,jj,1,Krhs) + (zhpj(ji,jj,1) + zvap) * vmask(ji,jj,1) 610 592 END_2D 611 !==================================================================================612 !===== Compute interior value ===================================================== 613 !================================================================================== 614 ! interior value (2=<jk=<jpkm1)593 ! 594 ! !=============================! 595 ! !===== interior values =====! 596 ! !=============================! 615 597 DO_3D( 0, 0, 0, 0, 2, jpkm1 ) 616 ! hydrostatic pressure gradient along s-surfaces 598 ze3w = e3w(ji ,jj ,jk,Kmm) 599 ze3wi1 = e3w(ji+1,jj ,jk,Kmm) 600 ze3wj1 = e3w(ji ,jj+1,jk,Kmm) 601 ! ! hydrostatic pressure gradient along s-surfaces 617 602 zhpi(ji,jj,jk) = zhpi(ji,jj,jk-1) + zcoef0 / e1u(ji,jj) & 618 & * ( e3w(ji+1,jj,jk,Kmm) & 619 & * ( rhd(ji+1,jj,jk) + rhd(ji+1,jj,jk-1) + 2*znad ) * wmask(ji+1,jj,jk) & 620 & - e3w(ji ,jj,jk,Kmm) & 621 & * ( rhd(ji ,jj,jk) + rhd(ji ,jj,jk-1) + 2*znad ) * wmask(ji ,jj,jk) ) 603 & * ( ze3wi1 * ( rhd(ji+1,jj,jk) + rhd(ji+1,jj,jk-1) ) * wmask(ji+1,jj,jk) & 604 & - ze3w * ( rhd(ji ,jj,jk) + rhd(ji ,jj,jk-1) ) * wmask(ji ,jj,jk) ) 622 605 zhpj(ji,jj,jk) = zhpj(ji,jj,jk-1) + zcoef0 / e2v(ji,jj) & 623 & * ( e3w(ji,jj+1,jk,Kmm) & 624 & * ( rhd(ji,jj+1,jk) + rhd(ji,jj+1,jk-1) + 2*znad ) * wmask(ji,jj+1,jk) & 625 & - e3w(ji,jj ,jk,Kmm) & 626 & * ( rhd(ji,jj, jk) + rhd(ji,jj ,jk-1) + 2*znad ) * wmask(ji,jj ,jk) ) 627 ! s-coordinate pressure gradient correction 628 zuap = -zcoef0 * ( rhd (ji+1,jj ,jk) + rhd (ji,jj,jk) + 2._wp * znad ) & 606 & * ( ze3wj1 * ( rhd(ji,jj+1,jk) + rhd(ji,jj+1,jk-1) ) * wmask(ji,jj+1,jk) & 607 & - ze3w * ( rhd(ji,jj, jk) + rhd(ji,jj ,jk-1) ) * wmask(ji,jj ,jk) ) 608 ! ! s-coordinate pressure gradient correction 609 zuap = -zcoef0 * ( rhd (ji+1,jj ,jk) + rhd (ji,jj,jk) ) & 629 610 & * ( gde3w(ji+1,jj ,jk) - gde3w(ji,jj,jk) ) / e1u(ji,jj) 630 zvap = -zcoef0 * ( rhd (ji ,jj+1,jk) + rhd (ji,jj,jk) + 2._wp * znad) &611 zvap = -zcoef0 * ( rhd (ji ,jj+1,jk) + rhd (ji,jj,jk) ) & 631 612 & * ( gde3w(ji ,jj+1,jk) - gde3w(ji,jj,jk) ) / e2v(ji,jj) 632 ! add to the general momentum trend613 ! ! add to the general momentum trend 633 614 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + (zhpi(ji,jj,jk) + zuap) * umask(ji,jj,jk) 634 615 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + (zhpj(ji,jj,jk) + zvap) * vmask(ji,jj,jk) … … 658 639 LOGICAL :: ll_tmp1, ll_tmp2 ! local logical variables 659 640 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zhpi, zhpj 660 REAL(wp), DIMENSION(jpi,jpj,jpk) :: rhd_opt661 641 662 642 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zdzx, zdzy, zdzz ! Primitive grid differences ('delta_xyz') … … 721 701 z1_12 = 1.0_wp / 12._wp 722 702 723 724 !! 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(:,:,:)725 IF( .NOT. ln_linssh ) THEN726 rhd_opt(:,:,:) = rhd(:,:,:) + 1.0_wp ! for vvl option727 ELSE728 rhd_opt(:,:,:) = rhd(:,:,:)729 END IF730 731 703 !---------------------------------------------------------------------------------------- 732 704 ! 1. compute and store elementary vertical differences in provisional arrays … … 736 708 737 709 DO_3D( 1, 1, 1, 1, 2, jpkm1 ) 738 zdrhoz(ji,jj,jk) = rhd _opt (ji ,jj ,jk) - rhd_opt(ji,jj,jk-1)710 zdrhoz(ji,jj,jk) = rhd (ji ,jj ,jk) - rhd (ji,jj,jk-1) 739 711 zdzz (ji,jj,jk) = - gde3w(ji ,jj ,jk) + gde3w(ji,jj,jk-1) 740 712 END_3D … … 763 735 764 736 ! mb for sea-ice shelves we will need to re-write this upper boundary condition in the same form as the lower boundary condition 765 zdrho_k(:,:,1) = aco_bc_vrt * ( rhd _opt (:,:,2) - rhd_opt(:,:,1) ) - bco_bc_vrt * zdrho_k(:,:,2)737 zdrho_k(:,:,1) = aco_bc_vrt * ( rhd (:,:,2) - rhd (:,:,1) ) - bco_bc_vrt * zdrho_k(:,:,2) 766 738 zdz_k (:,:,1) = aco_bc_vrt * (-gde3w(:,:,2) + gde3w(:,:,1) ) - bco_bc_vrt * zdz_k (:,:,2) 767 739 … … 769 741 IF ( mbkt(ji,jj)>1 ) THEN 770 742 iktb = mbkt(ji,jj) 771 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)743 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) 772 744 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) 773 745 END IF … … 787 759 DO_2D( 0, 1, 0, 1) 788 760 z_rho_k(ji,jj,1) = grav * ( ssh(ji,jj,Kmm) + gde3w(ji,jj,1) ) & 789 & * ( rhd _opt(ji,jj,1) &790 & + 0.5_wp * ( rhd _opt (ji,jj,2) - rhd_opt(ji,jj,1) ) &761 & * ( rhd(ji,jj,1) & 762 & + 0.5_wp * ( rhd (ji,jj,2) - rhd (ji,jj,1) ) & 791 763 & * ( ssh (ji,jj,Kmm) + gde3w(ji,jj,1) ) & 792 764 & / ( - gde3w(ji,jj,2) + gde3w(ji,jj,1) ) ) … … 798 770 799 771 DO_3D( 0, 1, 0, 1, 2, jpkm1 ) 800 z_rho_k(ji,jj,jk) = zcoef0 * ( rhd _opt (ji,jj,jk) + rhd_opt(ji,jj,jk-1) ) &772 z_rho_k(ji,jj,jk) = zcoef0 * ( rhd (ji,jj,jk) + rhd (ji,jj,jk-1) ) & 801 773 & * ( - gde3w(ji,jj,jk) + gde3w(ji,jj,jk-1) ) & 802 774 & + z_grav_10 * ( & … … 804 776 & * ( - gde3w(ji,jj,jk) + gde3w(ji,jj,jk-1) - z1_12 * ( zdz_k (ji,jj,jk) + zdz_k (ji,jj,jk-1) ) ) & 805 777 & - ( zdz_k (ji,jj,jk) - zdz_k (ji,jj,jk-1) ) & 806 & * ( rhd _opt (ji,jj,jk) - rhd_opt(ji,jj,jk-1) - z1_12 * ( zdrho_k(ji,jj,jk) + zdrho_k(ji,jj,jk-1) ) ) &778 & * ( rhd (ji,jj,jk) - rhd (ji,jj,jk-1) - z1_12 * ( zdrho_k(ji,jj,jk) + zdrho_k(ji,jj,jk-1) ) ) & 807 779 & ) 808 780 END_3D … … 813 785 814 786 DO_3D( 1, 0, 1, 0, 1, jpkm1 ) 815 zdrhox(ji,jj,jk) = rhd _opt (ji+1,jj ,jk) - rhd_opt(ji,jj,jk )787 zdrhox(ji,jj,jk) = rhd (ji+1,jj ,jk) - rhd (ji,jj,jk ) 816 788 zdzx (ji,jj,jk) = - gde3w(ji+1,jj ,jk) + gde3w(ji,jj,jk ) 817 zdrhoy(ji,jj,jk) = rhd _opt (ji ,jj+1,jk) - rhd_opt(ji,jj,jk )789 zdrhoy(ji,jj,jk) = rhd (ji ,jj+1,jk) - rhd (ji,jj,jk ) 818 790 zdzy (ji,jj,jk) = - gde3w(ji ,jj+1,jk) + gde3w(ji,jj,jk ) 819 791 END_3D … … 870 842 IF (ji < jpi) THEN 871 843 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 872 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)844 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) 873 845 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) 874 846 END IF … … 877 849 IF (ji > 2) THEN 878 850 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 879 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)851 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) 880 852 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) 881 853 END IF … … 884 856 IF (jj < jpj) THEN 885 857 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 886 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)858 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) 887 859 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) 888 860 END IF … … 891 863 IF (jj > 2) THEN 892 864 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 893 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)865 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) 894 866 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) 895 867 END IF … … 908 880 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) 909 881 ! two -ve signs cancel in next two lines (within zcoef0 and because gde3w is a depth not a height) 910 z_rho_i(ji,jj,jk) = zcoef0 * ( rhd _opt (ji+1,jj,jk) + rhd_opt(ji,jj,jk) ) &882 z_rho_i(ji,jj,jk) = zcoef0 * ( rhd (ji+1,jj,jk) + rhd (ji,jj,jk) ) & 911 883 & * ( gde3w(ji+1,jj,jk) - gde3w(ji,jj,jk) ) 912 884 IF ( umask(ji-1, jj, jk) > 0.5 .OR. umask(ji+1, jj, jk) > 0.5 ) THEN … … 915 887 & * ( - gde3w(ji+1,jj,jk) + gde3w(ji,jj,jk) - z1_12 * ( zdz_i (ji+1,jj,jk) + zdz_i (ji,jj,jk) ) ) & 916 888 & - ( zdz_i (ji+1,jj,jk) - zdz_i (ji,jj,jk) ) & 917 & * ( rhd _opt (ji+1,jj,jk) - rhd_opt(ji,jj,jk) - z1_12 * ( zdrho_i(ji+1,jj,jk) + zdrho_i(ji,jj,jk) ) ) &889 & * ( rhd (ji+1,jj,jk) - rhd (ji,jj,jk) - z1_12 * ( zdrho_i(ji+1,jj,jk) + zdrho_i(ji,jj,jk) ) ) & 918 890 & ) 919 891 END IF 920 892 921 z_rho_j(ji,jj,jk) = zcoef0 * ( rhd _opt (ji,jj+1,jk) + rhd_opt(ji,jj,jk) ) &893 z_rho_j(ji,jj,jk) = zcoef0 * ( rhd (ji,jj+1,jk) + rhd (ji,jj,jk) ) & 922 894 & * ( gde3w(ji,jj+1,jk) - gde3w(ji,jj,jk) ) 923 895 IF ( vmask(ji, jj-1, jk) > 0.5 .OR. vmask(ji, jj+1, jk) > 0.5 ) THEN … … 926 898 & * ( - gde3w(ji,jj+1,jk) + gde3w(ji,jj,jk) - z1_12 * ( zdz_j (ji,jj+1,jk) + zdz_j (ji,jj,jk) ) ) & 927 899 & - ( zdz_j (ji,jj+1,jk) - zdz_j (ji,jj,jk) ) & 928 & * ( rhd _opt (ji,jj+1,jk) - rhd_opt(ji,jj,jk) - z1_12 * ( zdrho_j(ji,jj+1,jk) + zdrho_j(ji,jj,jk) ) ) &900 & * ( rhd (ji,jj+1,jk) - rhd (ji,jj,jk) - z1_12 * ( zdrho_j(ji,jj+1,jk) + zdrho_j(ji,jj,jk) ) ) & 929 901 & ) 930 902 END IF … … 934 906 ! 8. Integrate in the vertical 935 907 !------------------------------------------------------------- 936 908 ! 937 909 ! --------------- 938 910 ! Surface value … … 1000 972 REAL(wp) :: zrhdt1 1001 973 REAL(wp) :: zdpdx1, zdpdx2, zdpdy1, zdpdy2 974 REAL(wp), DIMENSION(jpi,jpj) :: zpgu, zpgv ! 2D workspace 975 REAL(wp), DIMENSION(jpi,jpj) :: zsshu_n, zsshv_n 1002 976 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zdept, zrhh 1003 977 REAL(wp), DIMENSION(jpi,jpj,jpk) :: zhpi, zu, zv, fsp, xsp, asp, bsp, csp, dsp 1004 REAL(wp), DIMENSION(jpi,jpj) :: zsshu_n, zsshv_n1005 978 REAL(wp), DIMENSION(:,:), ALLOCATABLE :: zcpx, zcpy !W/D pressure filter 1006 979 !!---------------------------------------------------------------------- … … 1015 988 zcoef0 = - grav 1016 989 znad = 1._wp 1017 IF( ln_linssh ) znad = 0._wp 1018 990 IF( ln_linssh ) znad = 1._wp 991 ! 992 ! --------------- 993 ! Surface pressure gradient to be removed 994 ! --------------- 995 DO_2D( 0, 0, 0, 0 ) 996 zpgu(ji,jj) = - grav * ( ssh(ji+1,jj,Kmm) - ssh(ji,jj,Kmm) ) * r1_e1u(ji,jj) 997 zpgv(ji,jj) = - grav * ( ssh(ji,jj+1,Kmm) - ssh(ji,jj,Kmm) ) * r1_e2v(ji,jj) 998 END_2D 999 ! 1019 1000 IF( ln_wd_il ) THEN 1020 1001 ALLOCATE( zcpx(jpi,jpj) , zcpy(jpi,jpj) ) … … 1082 1063 ! Transfer the depth of "T(:,:,:)" to vertical coordinate "zdept(:,:,:)" 1083 1064 DO_2D( 1, 1, 1, 1 ) 1084 zdept(ji,jj,1) = 0.5_wp * e3w(ji,jj,1,Kmm) - ssh(ji,jj,Kmm) * znad1065 zdept(ji,jj,1) = 0.5_wp * e3w(ji,jj,1,Kmm) - ssh(ji,jj,Kmm) 1085 1066 END_2D 1086 1067 … … 1133 1114 1134 1115 DO_2D( 0, 0, 0, 0 ) 1135 zu(ji,jj,1) = - ( e3u(ji,jj,1,Kmm) - zsshu_n(ji,jj) * znad)1136 zv(ji,jj,1) = - ( e3v(ji,jj,1,Kmm) - zsshv_n(ji,jj) * znad)1116 zu(ji,jj,1) = - ( e3u(ji,jj,1,Kmm) - zsshu_n(ji,jj) ) 1117 zv(ji,jj,1) = - ( e3v(ji,jj,1,Kmm) - zsshv_n(ji,jj) ) 1137 1118 END_2D 1138 1119 … … 1216 1197 zdpdx2 = zdpdx2 * zcpx(ji,jj) * wdrampu(ji,jj) 1217 1198 ENDIF 1218 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + (zdpdx1 + zdpdx2 ) * umask(ji,jj,jk)1199 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + (zdpdx1 + zdpdx2 - zpgu(ji,jj)) * umask(ji,jj,jk) 1219 1200 ENDIF 1220 1201 … … 1276 1257 ENDIF 1277 1258 1278 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + (zdpdy1 + zdpdy2 ) * vmask(ji,jj,jk)1259 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + (zdpdy1 + zdpdy2 - zpgv(ji,jj)) * vmask(ji,jj,jk) 1279 1260 ENDIF 1280 1261 ! -
NEMO/trunk/src/OCE/DYN/dynspg.F90
r14007 r14064 82 82 INTEGER :: ji, jj, jk ! dummy loop indices 83 83 REAL(wp) :: z2dt, zg_2, zintp, zgrho0r, zld ! local scalars 84 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zpice 85 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdu, ztrdv 84 REAL(wp) , DIMENSION(jpi,jpj) :: zpgu, zpgv ! 2D workspace 85 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zpice 86 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrdu, ztrdv 86 87 !!---------------------------------------------------------------------- 87 88 ! … … 99 100 ! 100 101 DO_2D( 0, 0, 0, 0 ) 101 spgu(ji,jj) = 0._wp102 spgv(ji,jj) = 0._wp102 zpgu(ji,jj) = 0._wp 103 zpgv(ji,jj) = 0._wp 103 104 END_2D 104 105 ! … … 106 107 zg_2 = grav * 0.5 107 108 DO_2D( 0, 0, 0, 0 ) ! gradient of Patm using inverse barometer ssh 108 spgu(ji,jj) = spgu(ji,jj) + zg_2 * ( ssh_ib (ji+1,jj) - ssh_ib (ji,jj) &109 zpgu(ji,jj) = zpgu(ji,jj) + zg_2 * ( ssh_ib (ji+1,jj) - ssh_ib (ji,jj) & 109 110 & + ssh_ibb(ji+1,jj) - ssh_ibb(ji,jj) ) * r1_e1u(ji,jj) 110 spgv(ji,jj) = spgv(ji,jj) + zg_2 * ( ssh_ib (ji,jj+1) - ssh_ib (ji,jj) &111 zpgv(ji,jj) = zpgv(ji,jj) + zg_2 * ( ssh_ib (ji,jj+1) - ssh_ib (ji,jj) & 111 112 & + ssh_ibb(ji,jj+1) - ssh_ibb(ji,jj) ) * r1_e2v(ji,jj) 112 113 END_2D … … 121 122 ! 122 123 DO_2D( 0, 0, 0, 0 ) ! add tide potential forcing 123 spgu(ji,jj) = spgu(ji,jj) + grav * ( pot_astro(ji+1,jj) - pot_astro(ji,jj) ) * r1_e1u(ji,jj)124 spgv(ji,jj) = spgv(ji,jj) + grav * ( pot_astro(ji,jj+1) - pot_astro(ji,jj) ) * r1_e2v(ji,jj)124 zpgu(ji,jj) = zpgu(ji,jj) + grav * ( pot_astro(ji+1,jj) - pot_astro(ji,jj) ) * r1_e1u(ji,jj) 125 zpgv(ji,jj) = zpgv(ji,jj) + grav * ( pot_astro(ji,jj+1) - pot_astro(ji,jj) ) * r1_e2v(ji,jj) 125 126 END_2D 126 127 ! … … 128 129 zld = rn_scal_load * grav 129 130 DO_2D( 0, 0, 0, 0 ) ! add scalar approximation for load potential 130 spgu(ji,jj) = spgu(ji,jj) + zld * ( pssh(ji+1,jj,Kmm) - pssh(ji,jj,Kmm) ) * r1_e1u(ji,jj)131 spgv(ji,jj) = spgv(ji,jj) + zld * ( pssh(ji,jj+1,Kmm) - pssh(ji,jj,Kmm) ) * r1_e2v(ji,jj)131 zpgu(ji,jj) = zpgu(ji,jj) + zld * ( pssh(ji+1,jj,Kmm) - pssh(ji,jj,Kmm) ) * r1_e1u(ji,jj) 132 zpgv(ji,jj) = zpgv(ji,jj) + zld * ( pssh(ji,jj+1,Kmm) - pssh(ji,jj,Kmm) ) * r1_e2v(ji,jj) 132 133 END_2D 133 134 ENDIF … … 140 141 zpice(:,:) = ( zintp * snwice_mass(:,:) + ( 1.- zintp ) * snwice_mass_b(:,:) ) * zgrho0r 141 142 DO_2D( 0, 0, 0, 0 ) 142 spgu(ji,jj) = spgu(ji,jj) + ( zpice(ji+1,jj) - zpice(ji,jj) ) * r1_e1u(ji,jj)143 spgv(ji,jj) = spgv(ji,jj) + ( zpice(ji,jj+1) - zpice(ji,jj) ) * r1_e2v(ji,jj)143 zpgu(ji,jj) = zpgu(ji,jj) + ( zpice(ji+1,jj) - zpice(ji,jj) ) * r1_e1u(ji,jj) 144 zpgv(ji,jj) = zpgv(ji,jj) + ( zpice(ji,jj+1) - zpice(ji,jj) ) * r1_e2v(ji,jj) 144 145 END_2D 145 146 DEALLOCATE( zpice ) … … 148 149 IF( ln_wave .and. ln_bern_srfc ) THEN !== Add J terms: depth-independent Bernoulli head 149 150 DO_2D( 0, 0, 0, 0 ) 150 spgu(ji,jj) = spgu(ji,jj) + ( bhd_wave(ji+1,jj) - bhd_wave(ji,jj) ) / e1u(ji,jj) !++ bhd_wave from wave model in m2/s2 [BHD parameters in WW3]151 spgv(ji,jj) = spgv(ji,jj) + ( bhd_wave(ji,jj+1) - bhd_wave(ji,jj) ) / e2v(ji,jj)151 zpgu(ji,jj) = zpgu(ji,jj) + ( bhd_wave(ji+1,jj) - bhd_wave(ji,jj) ) / e1u(ji,jj) !++ bhd_wave from wave model in m2/s2 [BHD parameters in WW3] 152 zpgv(ji,jj) = zpgv(ji,jj) + ( bhd_wave(ji,jj+1) - bhd_wave(ji,jj) ) / e2v(ji,jj) 152 153 END_2D 153 154 ENDIF 154 155 ! 155 156 DO_3D( 0, 0, 0, 0, 1, jpkm1 ) !== Add all terms to the general trend 156 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + spgu(ji,jj)157 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + spgv(ji,jj)157 puu(ji,jj,jk,Krhs) = puu(ji,jj,jk,Krhs) + zpgu(ji,jj) 158 pvv(ji,jj,jk,Krhs) = pvv(ji,jj,jk,Krhs) + zpgv(ji,jj) 158 159 END_3D 159 160 ! -
NEMO/trunk/src/OCE/DYN/dynspg_exp.F90
r13497 r14064 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/trunk/src/OCE/DYN/dynspg_ts.F90
r14053 r14064 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/trunk/src/OCE/ISF/isf_oce.F90
r13558 r14064 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/trunk/src/OCE/ISF/isfload.F90
r13295 r14064 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/trunk/src/OCE/ISF/isfstp.F90
r13237 r14064 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/trunk/src/OCE/oce.F90
r14053 r14064 50 50 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ub2_i_b, vb2_i_b !: Half step time integrated fluxes 51 51 #endif 52 !53 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: spgu, spgv !: horizontal surface pressure gradient54 52 55 53 !! interpolated gradient (only used in zps case) … … 91 89 ! 92 90 ierr(:) = 0 93 ALLOCATE( uu (jpi,jpj,jpk,jpt) , vv (jpi,jpj,jpk,jpt) ,&94 & ww (jpi,jpj,jpk) , hdiv(jpi,jpj,jpk) ,&95 & ts (jpi,jpj,jpk,jpts,jpt) ,&96 & rab_b(jpi,jpj,jpk,jpts) , rab_n(jpi,jpj,jpk,jpts) ,&97 & rn2b (jpi,jpj,jpk) , rn2 (jpi,jpj,jpk) ,&98 & rhd (jpi,jpj,jpk) , rhop (jpi,jpj,jpk) 91 ALLOCATE( uu (jpi,jpj,jpk,jpt) , vv (jpi,jpj,jpk,jpt) , & 92 & ww (jpi,jpj,jpk) , hdiv(jpi,jpj,jpk) , & 93 & ts (jpi,jpj,jpk,jpts,jpt) , & 94 & rab_b(jpi,jpj,jpk,jpts) , rab_n(jpi,jpj,jpk,jpts) , & 95 & rn2b (jpi,jpj,jpk) , rn2 (jpi,jpj,jpk) , & 96 & rhd (jpi,jpj,jpk) , rhop (jpi,jpj,jpk) , STAT=ierr(1) ) 99 97 ! 100 ALLOCATE( ssh(jpi,jpj,jpt) , uu_b(jpi,jpj,jpt) , vv_b(jpi,jpj,jpt) , & 101 & spgu (jpi,jpj) , spgv(jpi,jpj) , & 102 & gtsu(jpi,jpj,jpts), gtsv(jpi,jpj,jpts) , & 103 & gru(jpi,jpj) , grv(jpi,jpj) , & 104 & gtui(jpi,jpj,jpts), gtvi(jpi,jpj,jpts) , & 105 & grui(jpi,jpj) , grvi(jpi,jpj) , & 106 & riceload(jpi,jpj) , STAT=ierr(2) ) 98 ALLOCATE( ssh (jpi,jpj,jpt) , uu_b(jpi,jpj,jpt) , vv_b(jpi,jpj,jpt) , & 99 & gtsu(jpi,jpj,jpts) , gtsv(jpi,jpj,jpts) , & 100 & gru (jpi,jpj) , grv (jpi,jpj) , & 101 & gtui(jpi,jpj,jpts) , gtvi(jpi,jpj,jpts) , & 102 & grui(jpi,jpj) , grvi(jpi,jpj) , & 103 & riceload(jpi,jpj) , STAT=ierr(2) ) 107 104 ! 108 105 ALLOCATE( fraqsr_1lev(jpi,jpj) , STAT=ierr(3) )
Note: See TracChangeset
for help on using the changeset viewer.
