- Timestamp:
- 2014-12-02T10:38:20+01:00 (9 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_MERGE_2014/NEMOGCM/NEMO/OPA_SRC/SBC/sbcisf.F90
r4938 r4946 54 54 REAL(wp) , PUBLIC, ALLOCATABLE, SAVE, DIMENSION (:,:) :: risfLeff !:effective length (Leff) BG03 nn_isf==2 55 55 REAL(wp) , PUBLIC, ALLOCATABLE, SAVE, DIMENSION (:,:) :: ttbl, stbl, utbl, vtbl !:top boundary layer variable at T point 56 INTEGER(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION (:,:) :: misfkt, misfkb !:Level of ice shelf base 56 #ifdef key_agrif 57 ! AGRIF can not handle these arrays as integers. The reason is a mystery but problems avoided by declaring them as reals 58 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION (:,:) :: misfkt, misfkb !:Level of ice shelf base 57 59 !: (first wet level and last level include in the tbl) 60 #else 61 INTEGER, PUBLIC, ALLOCATABLE, SAVE, DIMENSION (:,:) :: misfkt, misfkb !:Level of ice shelf base 62 #endif 63 58 64 59 65 REAL(wp), PUBLIC, SAVE :: rcpi = 2000.0_wp ! phycst ? … … 303 309 sbc_isf_alloc = 0 ! set to zero if no array to be allocated 304 310 IF( .NOT. ALLOCATED( qisf ) ) THEN 305 ALLOCATE( risf_tsc(jpi,jpj,jpts), risf_tsc_b(jpi,jpj,jpts), qisf(jpi,jpj), fwfisf(jpi,jpj), & 306 & fwfisf_b(jpi,jpj), misfkt(jpi,jpj), rhisf_tbl(jpi,jpj), r1_hisf_tbl(jpi,jpj), & 307 & rzisf_tbl(jpi,jpj), misfkb(jpi,jpj), ttbl(jpi,jpj), stbl(jpi,jpj), utbl(jpi,jpj), & 308 & vtbl(jpi, jpj), risfLeff(jpi,jpj), rhisf_tbl_0(jpi,jpj), ralpha(jpi,jpj), STAT= sbc_isf_alloc ) 311 ALLOCATE( risf_tsc(jpi,jpj,jpts), risf_tsc_b(jpi,jpj,jpts) , & 312 & qisf(jpi,jpj) , fwfisf(jpi,jpj) , fwfisf_b(jpi,jpj) , & 313 & rhisf_tbl(jpi,jpj), r1_hisf_tbl(jpi,jpj), rzisf_tbl(jpi,jpj) , & 314 & ttbl(jpi,jpj) , stbl(jpi,jpj) , utbl(jpi,jpj) , & 315 & vtbl(jpi, jpj) , risfLeff(jpi,jpj) , rhisf_tbl_0(jpi,jpj), & 316 & ralpha(jpi,jpj) , misfkt(jpi,jpj) , misfkb(jpi,jpj) , & 317 & STAT= sbc_isf_alloc ) 309 318 ! 310 319 IF( lk_mpp ) CALL mpp_sum ( sbc_isf_alloc ) … … 363 372 ! Calculate freezing temperature 364 373 zpress = grav*rau0*fsdept(ji,jj,ik)*1.e-04 365 zt_frz = tfreez1D(tsb(ji,jj,ik,jp_sal), zpress)374 zt_frz = eos_fzp(tsb(ji,jj,ik,jp_sal), zpress) 366 375 zt_sum = zt_sum + (tsn(ji,jj,ik,jp_tem)-zt_frz) * fse3t(ji,jj,ik) * tmask(ji,jj,ik) ! sum temp 367 376 ENDDO … … 445 454 zti(:,:)=tinsitu( ttbl, stbl, zpress ) 446 455 ! Calculate freezing temperature 447 zfrz(:,:)= tfreez( sss_m(:,:), zpress )456 zfrz(:,:)=eos_fzp( sss_m(:,:), zpress ) 448 457 449 458 … … 526 535 nit = nit + 1 527 536 IF (nit .GE. 51) THEN 528 WRITE(numout,*) "sbcisf : too many iteration ... ", zhtflx, zhtflx_b, zgammat, zgammas, nn_gammablk, ji, jj, mikt(ji,jj), narea 537 WRITE(numout,*) "sbcisf : too many iteration ... ", & 538 & zhtflx, zhtflx_b, zgammat, zgammas, nn_gammablk, ji, jj, mikt(ji,jj), narea 529 539 CALL ctl_stop( 'STOP', 'sbc_isf_hol99 : too many iteration ...' ) 530 540 END IF … … 584 594 REAL(wp) :: zgmolet, zgmoles, zgturb ! contribution of modelecular sublayer and turbulence 585 595 REAL(wp) :: zcoef ! temporary coef 586 REAL(wp) :: zrhos, zalbet, zbeta, zthermal, zhalin 587 REAL(wp) :: zt, zs, zh 596 REAL(wp) :: zdep 588 597 REAL(wp), PARAMETER :: zxsiN = 0.052 ! dimensionless constant 589 598 REAL(wp), PARAMETER :: epsln = 1.0e-20 ! a small positive number 590 599 REAL(wp), PARAMETER :: znu = 1.95e-6 ! kinamatic viscosity of sea water (m2.s-1) 591 600 REAL(wp) :: rcs = 1.0e-3_wp ! conversion: mm/s ==> m/s 601 REAL(wp), DIMENSION(2) :: zts, zab 592 602 !!--------------------------------------------------------------------- 593 603 ! … … 656 666 657 667 !! compute bouyancy 658 IF( nn_eos < 1) THEN 659 zt = ttbl(ji,jj) 660 zs = stbl(ji,jj) - 35.0 661 zh = fsdepw(ji,jj,ikt) 662 ! potential volumic mass 663 zrhos = rhop(ji,jj,ikt) 664 zalbet = ( ( ( - 0.255019e-07 * zt + 0.298357e-05 ) * zt & ! ratio alpha/beta 665 & - 0.203814e-03 ) * zt & 666 & + 0.170907e-01 ) * zt & 667 & + 0.665157e-01 & 668 & + ( - 0.678662e-05 * zs & 669 & - 0.846960e-04 * zt + 0.378110e-02 ) * zs & 670 & + ( ( - 0.302285e-13 * zh & 671 & - 0.251520e-11 * zs & 672 & + 0.512857e-12 * zt * zt ) * zh & 673 & - 0.164759e-06 * zs & 674 & +( 0.791325e-08 * zt - 0.933746e-06 ) * zt & 675 & + 0.380374e-04 ) * zh 676 677 zbeta = ( ( -0.415613e-09 * zt + 0.555579e-07 ) * zt & ! beta 678 & - 0.301985e-05 ) * zt & 679 & + 0.785567e-03 & 680 & + ( 0.515032e-08 * zs & 681 & + 0.788212e-08 * zt - 0.356603e-06 ) * zs & 682 & +( ( 0.121551e-17 * zh & 683 & - 0.602281e-15 * zs & 684 & - 0.175379e-14 * zt + 0.176621e-12 ) * zh & 685 & + 0.408195e-10 * zs & 686 & + ( - 0.213127e-11 * zt + 0.192867e-09 ) * zt & 687 & - 0.121555e-07 ) * zh 688 689 zthermal = zbeta * zalbet / ( rcp * zrhos + epsln ) 690 zhalin = zbeta * stbl(ji,jj) * rcs 691 ELSE 692 zrhos = rhop(ji,jj,ikt) + rau0 * ( 1. - tmask(ji,jj,ikt) ) 693 zthermal = rn_alpha / ( rcp * zrhos + epsln ) 694 zhalin = rn_beta * stbl(ji,jj) * rcs 695 ENDIF 668 zts(jp_tem) = ttbl(ji,jj) 669 zts(jp_sal) = stbl(ji,jj) 670 zdep = fsdepw(ji,jj,ikt) 671 ! 672 CALL eos_rab( zts, zdep, zab ) 673 ! 696 674 !! compute length scale 697 zbuofdep = grav * ( z thermal * zqhisf - zhalin* zqwisf ) !!!!!!!!!!!!!!!!!!!!!!!!!!!!675 zbuofdep = grav * ( zab(jp_tem) * zqhisf - zab(jp_sal) * zqwisf ) !!!!!!!!!!!!!!!!!!!!!!!!!!!! 698 676 699 677 !! compute Monin Obukov Length … … 766 744 ! level partially include in ice shelf boundary layer 767 745 zhk = SUM( fse3t_n(ji, jj, ikt:ikb - 1)) * r1_hisf_tbl(ji,jj) 768 IF (cptin == 'T') varout(ji,jj) = varout(ji,jj) + varin(ji,jj,ikb) * (1._wp - zhk) 769 IF (cptin == 'U') varout(ji,jj) = varout(ji,jj) + 0.5_wp * (varin(ji,jj,ikb) + varin(ji-1,jj,ikb)) * (1._wp - zhk) 770 IF (cptin == 'V') varout(ji,jj) = varout(ji,jj) + 0.5_wp * (varin(ji,jj,ikb) + varin(ji,jj-1,ikb)) * (1._wp - zhk) 746 IF (cptin == 'T') & 747 & varout(ji,jj) = varout(ji,jj) + varin(ji,jj,ikb) * (1._wp - zhk) 748 IF (cptin == 'U') & 749 & varout(ji,jj) = varout(ji,jj) + 0.5_wp * (varin(ji,jj,ikb) + varin(ji-1,jj,ikb)) * (1._wp - zhk) 750 IF (cptin == 'V') & 751 & varout(ji,jj) = varout(ji,jj) + 0.5_wp * (varin(ji,jj,ikb) + varin(ji,jj-1,ikb)) * (1._wp - zhk) 771 752 END IF 772 753 END DO … … 796 777 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: phdivn ! horizontal divergence 797 778 !! 798 INTEGER (wp):: ji, jj, jk ! dummy loop indices799 INTEGER (wp):: ikt, ikb800 INTEGER (wp):: nk_isf779 INTEGER :: ji, jj, jk ! dummy loop indices 780 INTEGER :: ikt, ikb 781 INTEGER :: nk_isf 801 782 REAL(wp) :: zhk, z1_hisf_tbl, zhisf_tbl 802 783 REAL(wp) :: zfact ! local scalar … … 834 815 ! level fully include in the ice shelf boundary layer 835 816 DO jk = ikt, ikb - 1 836 phdivn(ji,jj,jk) = phdivn(ji,jj,jk) + ( fwfisf(ji,jj) + fwfisf_b(ji,jj) ) * r1_hisf_tbl(ji,jj) * r1_rau0 * zfact 817 phdivn(ji,jj,jk) = phdivn(ji,jj,jk) + ( fwfisf(ji,jj) + fwfisf_b(ji,jj) ) & 818 & * r1_hisf_tbl(ji,jj) * r1_rau0 * zfact 837 819 END DO 838 820 ! level partially include in ice shelf boundary layer 839 phdivn(ji,jj,ikb) = phdivn(ji,jj,ikb) + ( fwfisf(ji,jj) + fwfisf_b(ji,jj) ) * r1_hisf_tbl(ji,jj) * r1_rau0 * zfact * ralpha(ji,jj) 821 phdivn(ji,jj,ikb) = phdivn(ji,jj,ikb) + ( fwfisf(ji,jj) & 822 & + fwfisf_b(ji,jj) ) * r1_hisf_tbl(ji,jj) * r1_rau0 * zfact * ralpha(ji,jj) 840 823 !== ice shelf melting mass distributed over several levels ==! 841 824 END DO
Note: See TracChangeset
for help on using the changeset viewer.