- Timestamp:
- 2015-12-07T16:11:45+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_MetOffice_merge_2015/NEMOGCM/NEMO/OPA_SRC/DYN/dynhpg.F90
r5930 r6012 45 45 USE wrk_nemo ! Memory Allocation 46 46 USE timing ! Timing 47 USE iom 47 48 48 49 IMPLICIT NONE … … 130 131 INTEGER :: ioptio = 0 ! temporary integer 131 132 INTEGER :: ios ! Local integer output status for namelist read 133 !! 134 INTEGER :: ji, jj, jk, ikt ! dummy loop indices ISF 135 REAL(wp) :: znad 136 REAL(wp), POINTER, DIMENSION(:,:,:) :: ztstop, zrhd ! hypothesys on isf density 137 REAL(wp), POINTER, DIMENSION(:,:) :: zrhdtop_isf ! density at bottom of ISF 138 REAL(wp), POINTER, DIMENSION(:,:) :: ziceload ! density at bottom of ISF 132 139 !! 133 140 NAMELIST/namdyn_hpg/ ln_hpg_zco, ln_hpg_zps, ln_hpg_sco, & … … 190 197 IF( ioptio /= 1 ) CALL ctl_stop( 'NO or several hydrostatic pressure gradient options used' ) 191 198 ! 192 ! initialisation of ice load 193 riceload(:,:)=0.0 199 ! initialisation of ice shelf load 200 IF ( .NOT. ln_isfcav ) riceload(:,:)=0.0 201 IF ( ln_isfcav ) THEN 202 CALL wrk_alloc( jpi,jpj, 2, ztstop) 203 CALL wrk_alloc( jpi,jpj,jpk, zrhd ) 204 CALL wrk_alloc( jpi,jpj, zrhdtop_isf, ziceload) 205 ! 206 IF(lwp) WRITE(numout,*) 207 IF(lwp) WRITE(numout,*) 'dyn:hpg_isf : hydrostatic pressure gradient trend for ice shelf' 208 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~' 209 210 ! To use density and not density anomaly 211 znad=1._wp 212 213 ! assume water displaced by the ice shelf is at T=-1.9 and S=34.4 (rude) 214 ztstop(:,:,1)=-1.9_wp ; ztstop(:,:,2)=34.4_wp 215 216 ! compute density of the water displaced by the ice shelf 217 DO jk = 1, jpk 218 CALL eos(ztstop(:,:,:),fsdept_n(:,:,jk),zrhd(:,:,jk)) 219 END DO 220 221 ! compute rhd at the ice/oce interface (ice shelf side) 222 CALL eos(ztstop,risfdep,zrhdtop_isf) 223 224 ! Surface value + ice shelf gradient 225 ! compute pressure due to ice shelf load (used to compute hpgi/j for all the level from 1 to miku/v) 226 ! divided by 2 later 227 ziceload = 0._wp 228 DO jj = 1, jpj 229 DO ji = 1, jpi 230 ikt=mikt(ji,jj) 231 ziceload(ji,jj) = ziceload(ji,jj) + (znad + zrhd(ji,jj,1) ) * fse3w(ji,jj,1) * (1._wp - tmask(ji,jj,1)) 232 DO jk=2,ikt-1 233 ziceload(ji,jj) = ziceload(ji,jj) + (2._wp * znad + zrhd(ji,jj,jk-1) + zrhd(ji,jj,jk)) * fse3w(ji,jj,jk) & 234 & * (1._wp - tmask(ji,jj,jk)) 235 END DO 236 IF (ikt >= 2) ziceload(ji,jj) = ziceload(ji,jj) + (2._wp * znad + zrhdtop_isf(ji,jj) + zrhd(ji,jj,ikt-1)) & 237 & * ( risfdep(ji,jj) - gdept_1d(ikt-1) ) 238 END DO 239 END DO 240 riceload(:,:)=ziceload(:,:) ! need to be saved for diaar5 241 242 CALL wrk_dealloc( jpi,jpj, 2, ztstop) 243 CALL wrk_dealloc( jpi,jpj,jpk, zrhd ) 244 CALL wrk_dealloc( jpi,jpj, zrhdtop_isf, ziceload) 245 END IF 194 246 ! 195 247 END SUBROUTINE dyn_hpg_init … … 447 499 SUBROUTINE hpg_isf( kt ) 448 500 !!--------------------------------------------------------------------- 449 !! *** ROUTINE hpg_ sco***501 !! *** ROUTINE hpg_isf *** 450 502 !! 451 503 !! ** Method : s-coordinate case. Jacobian scheme. … … 466 518 INTEGER, INTENT(in) :: kt ! ocean time-step index 467 519 !! 468 INTEGER :: ji, jj, jk, ik u, ikv, ikt, iktp1i, iktp1j! dummy loop indices469 REAL(wp) :: zcoef0, zuap, zvap, znad , ze3wu, ze3wv, zuapint, zvapint, zhpjint, zhpiint, zdzwt, zdzwtjp1, zdzwtip1! temporary scalars470 REAL(wp), POINTER, DIMENSION(:,:,:) :: zhpi, zhpj , zrhd520 INTEGER :: ji, jj, jk, ikt, iktp1i, iktp1j ! dummy loop indices 521 REAL(wp) :: zcoef0, zuap, zvap, znad ! temporary scalars 522 REAL(wp), POINTER, DIMENSION(:,:,:) :: zhpi, zhpj 471 523 REAL(wp), POINTER, DIMENSION(:,:,:) :: ztstop 472 REAL(wp), POINTER, DIMENSION(:,:) :: ze3w, zp, zrhdtop_isf, zrhdtop_oce, ziceload, zdept, zpshpi, zpshpj 473 !!---------------------------------------------------------------------- 474 ! 475 CALL wrk_alloc( jpi,jpj, 2, ztstop) 476 CALL wrk_alloc( jpi,jpj,jpk, zhpi, zhpj, zrhd) 477 CALL wrk_alloc( jpi,jpj, ze3w, zp, zrhdtop_isf, zrhdtop_oce, ziceload, zdept, zpshpi, zpshpj) 478 ! 479 IF( kt == nit000 ) THEN 480 IF(lwp) WRITE(numout,*) 481 IF(lwp) WRITE(numout,*) 'dyn:hpg_isf : hydrostatic pressure gradient trend for ice shelf' 482 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ s-coordinate case, OPA original scheme used' 483 ENDIF 484 524 REAL(wp), POINTER, DIMENSION(:,:) :: zrhdtop_oce 525 !!---------------------------------------------------------------------- 526 ! 527 CALL wrk_alloc( jpi,jpj, 2, ztstop) 528 CALL wrk_alloc( jpi,jpj,jpk, zhpi, zhpj) 529 CALL wrk_alloc( jpi,jpj, zrhdtop_oce ) 530 ! 485 531 ! Local constant initialization 486 532 zcoef0 = - grav * 0.5_wp 533 487 534 ! To use density and not density anomaly 488 ! IF ( lk_vvl ) THEN ; znad = 1._wp ! Variable volume489 ! ELSE ; znad = 0._wp ! Fixed volume490 ! ENDIF491 535 znad=1._wp 536 492 537 ! iniitialised to 0. zhpi zhpi 493 538 zhpi(:,:,:)=0._wp ; zhpj(:,:,:)=0._wp 494 539 495 ! Partial steps: top & bottom before horizontal gradient496 !jc CALL zps_hde_isf( kt, jpts, tsn, gtsu, gtsv, gtui, gtvi, &497 !jc & rhd, gru , grv , aru , arv , gzu , gzv , ge3ru , ge3rv , &498 !jc & grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi )499 500 !==================================================================================501 !=====Compute iceload and contribution of the half first wet layer =================502 !===================================================================================503 504 ! assume water displaced by the ice shelf is at T=-1.9 and S=34.4 (rude)505 ztstop(:,:,1)=-1.9_wp ; ztstop(:,:,2)=34.4_wp506 507 ! compute density of the water displaced by the ice shelf508 zrhd = rhd ! save rhd509 DO jk = 1, jpk510 zdept(:,:)=gdept_1d(jk)511 CALL eos(ztstop(:,:,:),zdept(:,:),rhd(:,:,jk))512 END DO513 WHERE ( tmask(:,:,:) == 1._wp)514 rhd(:,:,:) = zrhd(:,:,:) ! replace wet cell by the saved rhd515 END WHERE516 517 ! compute rhd at the ice/oce interface (ice shelf side)518 CALL eos(ztstop,risfdep,zrhdtop_isf)519 520 540 ! compute rhd at the ice/oce interface (ocean side) 541 ! usefull to reduce residual current in the test case ISOMIP with no melting 521 542 DO ji=1,jpi 522 543 DO jj=1,jpj … … 527 548 END DO 528 549 CALL eos(ztstop,risfdep,zrhdtop_oce) 529 ! 530 ! Surface value + ice shelf gradient 531 ! compute pressure due to ice shelf load (used to compute hpgi/j for all the level from 1 to miku/v) 532 ziceload = 0._wp 533 DO jj = 1, jpj 534 DO ji = 1, jpi ! vector opt. 535 ikt=mikt(ji,jj) 536 ziceload(ji,jj) = ziceload(ji,jj) + (znad + rhd(ji,jj,1) ) * fse3w(ji,jj,1) * (1._wp - tmask(ji,jj,1)) 537 DO jk=2,ikt-1 538 ziceload(ji,jj) = ziceload(ji,jj) + (2._wp * znad + rhd(ji,jj,jk-1) + rhd(ji,jj,jk)) * fse3w(ji,jj,jk) & 539 & * (1._wp - tmask(ji,jj,jk)) 540 END DO 541 IF (ikt .GE. 2) ziceload(ji,jj) = ziceload(ji,jj) + (2._wp * znad + zrhdtop_isf(ji,jj) + rhd(ji,jj,ikt-1)) & 542 & * ( risfdep(ji,jj) - gdept_1d(ikt-1) ) 543 END DO 544 END DO 545 riceload(:,:) = 0.0_wp ; riceload(:,:)=ziceload(:,:) ! need to be saved for diaar5 546 ! compute zp from z=0 to first T wet point (correction due to zps not yet applied) 550 551 !================================================================================== 552 !===== Compute surface value ===================================================== 553 !================================================================================== 547 554 DO jj = 2, jpjm1 548 555 DO ji = fs_2, fs_jpim1 ! vector opt. … … 550 557 ! hydrostatic pressure gradient along s-surfaces and ice shelf pressure 551 558 ! we assume ISF is in isostatic equilibrium 552 zhpi(ji,jj,1) = zcoef0 / e1u(ji,jj) * ( 0.5_wp * fse3w(ji+1,jj 553 & * ( 2._wp * znad + rhd(ji+1,jj ,iktp1i) + zrhdtop_oce(ji+1,jj) ) &554 & - 0.5_wp * fse3w(ji ,jj ,ikt )&555 & * ( 2._wp * znad + rhd(ji ,jj ,ikt ) + zrhdtop_oce(ji ,jj ) )&556 & + ( ziceload(ji+1,jj) - ziceload(ji,jj)))557 zhpj(ji,jj,1) = zcoef0 / e2v(ji,jj) * ( 0.5_wp * fse3w(ji 558 & * ( 2._wp * znad + rhd(ji ,jj+1,iktp1j) + zrhdtop_oce(ji,jj+1) ) &559 & - 0.5_wp * fse3w(ji ,jj ,ikt )&560 & * ( 2._wp * znad + rhd(ji ,jj ,ikt ) + zrhdtop_oce(ji ,jj ) )&561 & + ( ziceload(ji,jj+1) - ziceload(ji,jj) ))559 zhpi(ji,jj,1) = zcoef0 / e1u(ji,jj) * ( 0.5_wp * fse3w(ji+1,jj,iktp1i) & 560 & * ( 2._wp * znad + rhd(ji+1,jj,iktp1i) + zrhdtop_oce(ji+1,jj) ) & 561 & - 0.5_wp * fse3w(ji,jj,ikt) & 562 & * ( 2._wp * znad + rhd(ji,jj,ikt) + zrhdtop_oce(ji,jj) ) & 563 & + ( riceload(ji+1,jj) - riceload(ji,jj)) ) 564 zhpj(ji,jj,1) = zcoef0 / e2v(ji,jj) * ( 0.5_wp * fse3w(ji,jj+1,iktp1j) & 565 & * ( 2._wp * znad + rhd(ji,jj+1,iktp1j) + zrhdtop_oce(ji,jj+1) ) & 566 & - 0.5_wp * fse3w(ji,jj,ikt) & 567 & * ( 2._wp * znad + rhd(ji,jj,ikt) + zrhdtop_oce(ji,jj) ) & 568 & + ( riceload(ji,jj+1) - riceload(ji,jj)) ) 562 569 ! s-coordinate pressure gradient correction (=0 if z coordinate) 563 570 zuap = -zcoef0 * ( rhd (ji+1,jj,1) + rhd (ji,jj,1) + 2._wp * znad ) & … … 570 577 END DO 571 578 END DO 572 !==================================================================================573 !===== Compute partial cell contribution for the top cell =========================574 !==================================================================================575 DO jj = 2, jpjm1576 DO ji = fs_2, fs_jpim1 ! vector opt.577 iku = miku(ji,jj) ;578 zpshpi(ji,jj)=0.0_wp ; zpshpj(ji,jj)=0.0_wp579 ze3wu = (gdepw_0(ji+1,jj,iku+1) - gdept_0(ji+1,jj,iku)) - (gdepw_0(ji,jj,iku+1) - gdept_0(ji,jj,iku))580 ! u direction581 IF ( iku .GT. 1 ) THEN582 ! case iku583 zhpi(ji,jj,iku) = zcoef0 / e1u(ji,jj) * ze3wu &584 & * ( rhd (ji+1,jj,iku) + rhd (ji,jj,iku) &585 & + SIGN(1._wp,ze3wu) * grui(ji,jj) + 2._wp * znad )586 ! corrective term ( = 0 if z coordinate )587 zuap = -zcoef0 * ( arui(ji,jj) + 2._wp * znad ) * gzui(ji,jj) / e1u(ji,jj)588 ! zhpi will be added in interior loop589 ua(ji,jj,iku) = ua(ji,jj,iku) + zuap590 ! in case of 2 cell water column, need to save the pressure gradient to compute the bottom pressure591 IF (mbku(ji,jj) == iku + 1) zpshpi(ji,jj) = zhpi(ji,jj,iku)592 593 ! case iku + 1 (remove the zphi term added in the interior loop and compute the one corrected for zps)594 zhpiint = zcoef0 / e1u(ji,jj) &595 & * ( fse3w(ji+1,jj ,iku+1) * ( (rhd(ji+1,jj,iku+1) + znad) &596 & + (rhd(ji+1,jj,iku ) + znad) ) * tmask(ji+1,jj,iku) &597 & - fse3w(ji ,jj ,iku+1) * ( (rhd(ji ,jj,iku+1) + znad) &598 & + (rhd(ji ,jj,iku ) + znad) ) * tmask(ji ,jj,iku) )599 zhpi(ji,jj,iku+1) = zcoef0 / e1u(ji,jj) * ge3rui(ji,jj) - zhpiint600 END IF601 602 ! v direction603 ikv = mikv(ji,jj)604 ze3wv = (gdepw_0(ji,jj+1,ikv+1) - gdept_0(ji,jj+1,ikv)) - (gdepw_0(ji,jj,ikv+1) - gdept_0(ji,jj,ikv))605 IF ( ikv .GT. 1 ) THEN606 ! case ikv607 zhpj(ji,jj,ikv) = zcoef0 / e2v(ji,jj) * ze3wv &608 & * ( rhd(ji,jj+1,ikv) + rhd (ji,jj,ikv) &609 & + SIGN(1._wp,ze3wv) * grvi(ji,jj) + 2._wp * znad )610 ! corrective term ( = 0 if z coordinate )611 zvap = -zcoef0 * ( arvi(ji,jj) + 2._wp * znad ) * gzvi(ji,jj) / e2v(ji,jj)612 ! zhpi will be added in interior loop613 va(ji,jj,ikv) = va(ji,jj,ikv) + zvap614 ! in case of 2 cell water column, need to save the pressure gradient to compute the bottom pressure615 IF (mbkv(ji,jj) == ikv + 1) zpshpj(ji,jj) = zhpj(ji,jj,ikv)616 617 ! case ikv + 1 (remove the zphj term added in the interior loop and compute the one corrected for zps)618 zhpjint = zcoef0 / e2v(ji,jj) &619 & * ( fse3w(ji ,jj+1,ikv+1) * ( (rhd(ji,jj+1,ikv+1) + znad) &620 & + (rhd(ji,jj+1,ikv ) + znad) ) * tmask(ji,jj+1,ikv) &621 & - fse3w(ji ,jj ,ikv+1) * ( (rhd(ji,jj ,ikv+1) + znad) &622 & + (rhd(ji,jj ,ikv ) + znad) ) * tmask(ji,jj ,ikv) )623 zhpj(ji,jj,ikv+1) = zcoef0 / e2v(ji,jj) * ge3rvi(ji,jj) - zhpjint624 END IF625 END DO626 END DO627 579 628 580 !================================================================================== 629 581 !===== Compute interior value ===================================================== 630 582 !================================================================================== 631 632 DO jj = 2, jpjm1 633 DO ji = fs_2, fs_jpim1 ! vector opt. 634 iku=miku(ji,jj); ikv=mikv(ji,jj) 635 DO jk = 2, jpkm1 583 ! interior value (2=<jk=<jpkm1) 584 DO jk = 2, jpkm1 585 DO jj = 2, jpjm1 586 DO ji = fs_2, fs_jpim1 ! vector opt. 636 587 ! hydrostatic pressure gradient along s-surfaces 637 ! zhpi is masked for the first wet cell (contribution already done in the upper bloc) 638 zhpi(ji,jj,jk) = zhpi(ji,jj,jk) + zhpi(ji,jj,jk-1) & 639 & + zcoef0 / e1u(ji,jj) & 640 & * ( fse3w(ji+1,jj ,jk) * ( (rhd(ji+1,jj,jk ) + znad) & 641 & + (rhd(ji+1,jj,jk-1) + znad) ) * tmask(ji+1,jj,jk-1) & 642 & - fse3w(ji ,jj ,jk) * ( (rhd(ji ,jj,jk ) + znad) & 643 & + (rhd(ji ,jj,jk-1) + znad) ) * tmask(ji ,jj,jk-1) ) 588 zhpi(ji,jj,jk) = zhpi(ji,jj,jk-1) + zcoef0 / e1u(ji,jj) & 589 & * ( fse3w(ji+1,jj,jk) * ( rhd(ji+1,jj,jk) + rhd(ji+1,jj,jk-1) + 2*znad ) * wmask(ji+1,jj,jk) & 590 & - fse3w(ji ,jj,jk) * ( rhd(ji ,jj,jk) + rhd(ji ,jj,jk-1) + 2*znad ) * wmask(ji ,jj,jk) ) 591 zhpj(ji,jj,jk) = zhpj(ji,jj,jk-1) + zcoef0 / e2v(ji,jj) & 592 & * ( fse3w(ji,jj+1,jk) * ( rhd(ji,jj+1,jk) + rhd(ji,jj+1,jk-1) + 2*znad ) * wmask(ji,jj+1,jk) & 593 & - fse3w(ji,jj ,jk) * ( rhd(ji,jj, jk) + rhd(ji,jj ,jk-1) + 2*znad ) * wmask(ji,jj ,jk) ) 644 594 ! s-coordinate pressure gradient correction 645 ! corrective term, we mask this term for the first wet level beneath the ice shelf (contribution done in the upper bloc) 646 zuap = - zcoef0 * ( rhd (ji+1,jj ,jk) + rhd (ji,jj,jk) + 2._wp * znad ) & 647 & * ( fsde3w(ji+1,jj ,jk) - fsde3w(ji,jj,jk) ) / e1u(ji,jj) * umask(ji,jj,jk-1) 648 ua(ji,jj,jk) = ua(ji,jj,jk) + ( zhpi(ji,jj,jk) + zuap) * umask(ji,jj,jk) 649 650 ! hydrostatic pressure gradient along s-surfaces 651 ! zhpi is masked for the first wet cell (contribution already done in the upper bloc) 652 zhpj(ji,jj,jk) = zhpj(ji,jj,jk) + zhpj(ji,jj,jk-1) & 653 & + zcoef0 / e2v(ji,jj) & 654 & * ( fse3w(ji ,jj+1,jk) * ( (rhd(ji,jj+1,jk ) + znad) & 655 & + (rhd(ji,jj+1,jk-1) + znad) ) * tmask(ji,jj+1,jk-1) & 656 & - fse3w(ji ,jj ,jk) * ( (rhd(ji,jj ,jk ) + znad) & 657 & + (rhd(ji,jj ,jk-1) + znad) ) * tmask(ji,jj ,jk-1) ) 658 ! s-coordinate pressure gradient correction 659 ! corrective term, we mask this term for the first wet level beneath the ice shelf (contribution done in the upper bloc) 660 zvap = - zcoef0 * ( rhd (ji ,jj+1,jk) + rhd (ji,jj,jk) + 2._wp * znad ) & 661 & * ( fsde3w(ji ,jj+1,jk) - fsde3w(ji,jj,jk) ) / e2v(ji,jj) * vmask(ji,jj,jk-1) 595 zuap = -zcoef0 * ( rhd (ji+1,jj ,jk) + rhd (ji,jj,jk) + 2._wp * znad ) & 596 & * ( fsde3w(ji+1,jj ,jk) - fsde3w(ji,jj,jk) ) / e1u(ji,jj) 597 zvap = -zcoef0 * ( rhd (ji ,jj+1,jk) + rhd (ji,jj,jk) + 2._wp * znad ) & 598 & * ( fsde3w(ji ,jj+1,jk) - fsde3w(ji,jj,jk) ) / e2v(ji,jj) 662 599 ! add to the general momentum trend 663 va(ji,jj,jk) = va(ji,jj,jk) + ( zhpj(ji,jj,jk) + zvap ) * vmask(ji,jj,jk) 600 ua(ji,jj,jk) = ua(ji,jj,jk) + (zhpi(ji,jj,jk) + zuap) * umask(ji,jj,jk) 601 va(ji,jj,jk) = va(ji,jj,jk) + (zhpj(ji,jj,jk) + zvap) * vmask(ji,jj,jk) 664 602 END DO 665 603 END DO 666 604 END DO 667 668 !================================================================================== 669 !===== Compute bottom cell contribution (partial cell) ============================ 670 !================================================================================== 671 672 DO jj = 2, jpjm1 673 DO ji = 2, jpim1 674 iku = mbku(ji,jj) 675 ikv = mbkv(ji,jj) 676 677 IF (iku .GT. 1) THEN 678 ! remove old value (interior case) 679 zuap = -zcoef0 * ( rhd (ji+1,jj ,iku) + rhd (ji,jj,iku) + 2._wp * znad ) & 680 & * ( fsde3w(ji+1,jj ,iku) - fsde3w(ji,jj,iku) ) / e1u(ji,jj) 681 ua(ji,jj,iku) = ua(ji,jj,iku) - zhpi(ji,jj,iku) - zuap 682 ! put new value 683 ! -zpshpi to avoid double contribution of the partial step in the top layer 684 zuap = -zcoef0 * ( aru(ji,jj) + 2._wp * znad ) * gzu(ji,jj) / e1u(ji,jj) 685 zhpi(ji,jj,iku) = zhpi(ji,jj,iku-1) + zcoef0 / e1u(ji,jj) * ge3ru(ji,jj) - zpshpi(ji,jj) 686 ua(ji,jj,iku) = ua(ji,jj,iku) + zhpi(ji,jj,iku) + zuap 687 END IF 688 ! v direction 689 IF (ikv .GT. 1) THEN 690 ! remove old value (interior case) 691 zvap = -zcoef0 * ( rhd (ji ,jj+1,ikv) + rhd (ji,jj,ikv) + 2._wp * znad ) & 692 & * ( fsde3w(ji ,jj+1,ikv) - fsde3w(ji,jj,ikv) ) / e2v(ji,jj) 693 va(ji,jj,ikv) = va(ji,jj,ikv) - zhpj(ji,jj,ikv) - zvap 694 ! put new value 695 ! -zpshpj to avoid double contribution of the partial step in the top layer 696 zvap = -zcoef0 * ( arv(ji,jj) + 2._wp * znad ) * gzv(ji,jj) / e2v(ji,jj) 697 zhpj(ji,jj,ikv) = zhpj(ji,jj,ikv-1) + zcoef0 / e2v(ji,jj) * ge3rv(ji,jj) - zpshpj(ji,jj) 698 va(ji,jj,ikv) = va(ji,jj,ikv) + zhpj(ji,jj,ikv) + zvap 699 END IF 700 END DO 701 END DO 702 703 ! set back to original density value into the ice shelf cell (maybe useless because it is masked) 704 rhd = zrhd 705 ! 706 CALL wrk_dealloc( jpi,jpj,2, ztstop) 707 CALL wrk_dealloc( jpi,jpj,jpk, zhpi, zhpj, zrhd) 708 CALL wrk_dealloc( jpi,jpj, ze3w, zp, zrhdtop_isf, zrhdtop_oce, ziceload, zdept, zpshpi, zpshpj) 605 ! 606 CALL wrk_dealloc( jpi,jpj,2 , ztstop) 607 CALL wrk_dealloc( jpi,jpj,jpk, zhpi, zhpj) 608 CALL wrk_dealloc( jpi,jpj , zrhdtop_oce ) 709 609 ! 710 610 END SUBROUTINE hpg_isf
Note: See TracChangeset
for help on using the changeset viewer.
