- Timestamp:
- 2016-11-28T18:21:42+01:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2016/dev_CNRS_2016/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk.F90
r7282 r7355 40 40 USE lib_fortran ! to use key_nosignedzero 41 41 #if defined key_lim3 42 USE ice , ONLY : u_ice, v_ice, jpl, pfrld, a_i_b 42 USE ice , ONLY : u_ice, v_ice, jpl, pfrld, a_i_b, at_i_b 43 43 USE limthd_dh ! for CALL lim_thd_snwblow 44 44 #elif defined key_lim2 … … 93 93 94 94 ! !!! Bulk parameters 95 REAL(wp), PARAMETER :: cpa = 1000.5 ! specific heat of air (only used for ice fluxes now...)96 REAL(wp), PARAMETER :: Ls = 2.839e6 ! latent heat of sublimation97 REAL(wp), PARAMETER :: Stef = 5.67e-8 ! Stefan Boltzmann constant98 REAL(wp), PARAMETER :: C ice = 1.4e-3 ! iovi 1.63e-3 ! transfer coefficient over ice99 REAL(wp), PARAMETER :: albo = 0.066 ! ocean albedo assumed to be constant95 REAL(wp), PARAMETER :: cpa = 1000.5 ! specific heat of air (only used for ice fluxes now...) 96 REAL(wp), PARAMETER :: Ls = 2.839e6 ! latent heat of sublimation 97 REAL(wp), PARAMETER :: Stef = 5.67e-8 ! Stefan Boltzmann constant 98 REAL(wp), PARAMETER :: Cd_ice = 1.4e-3 ! iovi 1.63e-3 ! transfer coefficient over ice 99 REAL(wp), PARAMETER :: albo = 0.066 ! ocean albedo assumed to be constant 100 100 ! 101 101 ! !!* Namelist namsbc_blk : bulk parameters … … 111 111 REAL(wp) :: rn_zqt ! z(q,t) : height of humidity and temperature measurements 112 112 REAL(wp) :: rn_zu ! z(u) : height of wind measurements 113 LOGICAL :: ln_Cd_L12 = .FALSE. ! Modify the drag ice-atm and oce-atm depending on ice concentration (from Lupkes et al. JGR2012) 114 ! 115 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: Cd_oce ! air-ocean drag (clem) 113 116 114 117 INTEGER :: nblk ! choice of the bulk algorithm … … 128 131 CONTAINS 129 132 133 INTEGER FUNCTION sbc_blk_alloc() 134 !!------------------------------------------------------------------- 135 !! *** ROUTINE sbc_blk_alloc *** 136 !!------------------------------------------------------------------- 137 ALLOCATE( Cd_oce(jpi,jpj) , STAT=sbc_blk_alloc ) 138 ! 139 IF( lk_mpp ) CALL mpp_sum ( sbc_blk_alloc ) 140 IF( sbc_blk_alloc /= 0 ) CALL ctl_warn('sbc_blk_alloc: failed to allocate arrays') 141 END FUNCTION sbc_blk_alloc 142 130 143 SUBROUTINE sbc_blk_init 131 144 !!--------------------------------------------------------------------- … … 153 166 & sn_tair, sn_prec, sn_snow, sn_slp, sn_tdif, & 154 167 & ln_NCAR, ln_COARE_3p0, ln_COARE_3p5, ln_ECMWF, & ! bulk algorithm 155 & cn_dir , ln_taudif, rn_zqt, rn_zu, rn_pfac, rn_efac, rn_vfac 156 !!--------------------------------------------------------------------- 168 & cn_dir , ln_taudif, rn_zqt, rn_zu, & 169 & rn_pfac, rn_efac, rn_vfac, ln_Cd_L12 170 !!--------------------------------------------------------------------- 171 ! 172 ! ! allocate sbc_blk_core array 173 IF( sbc_blk_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'sbc_blk : unable to allocate standard arrays' ) 157 174 ! 158 175 ! !** read bulk namelist … … 425 442 END IF 426 443 444 Cd_oce(:,:) = Cd(:,:) ! record value of pure ocean-atm. drag (clem) 445 427 446 DO jj = 1, jpj ! tau module, i and j component 428 447 DO ji = 1, jpi … … 549 568 INTEGER :: ji, jj ! dummy loop indices 550 569 ! 551 REAL(wp), DIMENSION(:,:) , POINTER :: z coef_wnorm570 REAL(wp), DIMENSION(:,:) , POINTER :: zrhoa 552 571 ! 553 572 REAL(wp) :: zwnorm_f, zwndi_f , zwndj_f ! relative wind module and components at F-point 554 573 REAL(wp) :: zwndi_t , zwndj_t ! relative wind components at T-point 574 REAL(wp), DIMENSION(:,:), POINTER :: Cd ! transfer coefficient for momentum (tau) 555 575 !!--------------------------------------------------------------------- 556 576 ! 557 577 IF( nn_timing == 1 ) CALL timing_start('blk_ice_tau') 558 578 ! 559 CALL wrk_alloc( jpi,jpj, zcoef_wnorm ) 579 CALL wrk_alloc( jpi,jpj, zrhoa ) 580 CALL wrk_alloc( jpi,jpj, Cd ) 581 582 Cd(:,:) = Cd_ice 583 584 ! Make ice-atm. drag dependent on ice concentration (see Lupkes et al. 2012) (clem) 585 #if defined key_lim3 586 IF( ln_Cd_L12 ) THEN 587 CALL Cdn10_Lupkes2012( Cd ) ! calculate new drag from Lupkes(2012) equations 588 ENDIF 589 #endif 560 590 561 591 ! local scalars ( place there for vector optimisation purposes) 562 592 ! Computing density of air! Way denser that 1.2 over sea-ice !!! 563 593 !! 564 zcoef_wnorm (:,:) = rho_air(sf(jp_tair)%fnow(:,:,1), sf(jp_humi)%fnow(:,:,1), sf(jp_slp)%fnow(:,:,1)) 565 zcoef_wnorm (:,:) = Cice * zcoef_wnorm (:,:) 594 zrhoa (:,:) = rho_air(sf(jp_tair)%fnow(:,:,1), sf(jp_humi)%fnow(:,:,1), sf(jp_slp)%fnow(:,:,1)) 566 595 567 596 !!gm brutal.... … … 584 613 zwndj_f = 0.25 * ( sf(jp_wndj)%fnow(ji-1,jj ,1) + sf(jp_wndj)%fnow(ji ,jj ,1) & 585 614 & + sf(jp_wndj)%fnow(ji-1,jj-1,1) + sf(jp_wndj)%fnow(ji ,jj-1,1) ) - rn_vfac * v_ice(ji,jj) 586 zwnorm_f = z coef_wnorm(ji,jj) * SQRT( zwndi_f * zwndi_f + zwndj_f * zwndj_f )615 zwnorm_f = zrhoa(ji,jj) * Cd(ji,jj) * SQRT( zwndi_f * zwndi_f + zwndj_f * zwndj_f ) 587 616 ! ... ice stress at I-point 588 617 utau_ice(ji,jj) = zwnorm_f * zwndi_f … … 610 639 DO jj = 2, jpjm1 611 640 DO ji = fs_2, fs_jpim1 ! vect. opt. 612 utau_ice(ji,jj) = 0.5 *zcoef_wnorm(ji,jj) * ( wndm_ice(ji+1,jj ) + wndm_ice(ji,jj) ) &641 utau_ice(ji,jj) = 0.5 * zrhoa(ji,jj) * Cd(ji,jj) * ( wndm_ice(ji+1,jj ) + wndm_ice(ji,jj) ) & 613 642 & * ( 0.5 * (sf(jp_wndi)%fnow(ji+1,jj,1) + sf(jp_wndi)%fnow(ji,jj,1) ) - rn_vfac * u_ice(ji,jj) ) 614 vtau_ice(ji,jj) = 0.5 *zcoef_wnorm(ji,jj) * ( wndm_ice(ji,jj+1 ) + wndm_ice(ji,jj) ) &643 vtau_ice(ji,jj) = 0.5 * zrhoa(ji,jj) * Cd(ji,jj) * ( wndm_ice(ji,jj+1 ) + wndm_ice(ji,jj) ) & 615 644 & * ( 0.5 * (sf(jp_wndj)%fnow(ji,jj+1,1) + sf(jp_wndj)%fnow(ji,jj,1) ) - rn_vfac * v_ice(ji,jj) ) 616 645 END DO … … 626 655 CALL prt_ctl(tab2d_1=wndm_ice , clinfo1=' blk_ice: wndm_ice : ') 627 656 ENDIF 628 629 CALL wrk_dealloc( jpi,jpj, zcoef_wnorm )630 657 631 658 IF( nn_timing == 1 ) CALL timing_stop('blk_ice_tau') … … 659 686 REAL(wp), DIMENSION(:,:) , POINTER :: zevap, zsnw ! evaporation and snw distribution after wind blowing (LIM3) 660 687 REAL(wp), DIMENSION(:,:) , POINTER :: zrhoa 688 REAL(wp), DIMENSION(:,:) , POINTER :: Cd ! transfer coefficient for momentum (tau) 661 689 !!--------------------------------------------------------------------- 662 690 ! … … 665 693 CALL wrk_alloc( jpi,jpj,jpl, z_qlw, z_qsb, z_dqlw, z_dqsb ) 666 694 CALL wrk_alloc( jpi,jpj, zrhoa) 695 CALL wrk_alloc( jpi,jpj, Cd ) 696 697 Cd(:,:) = Cd_ice 698 699 ! Make ice-atm. drag dependent on ice concentration (see Lupkes et al. 2012) (clem) 700 #if defined key_lim3 701 IF( ln_Cd_L12 ) THEN 702 CALL Cdn10_Lupkes2012( Cd ) ! calculate new drag from Lupkes(2012) equations 703 ENDIF 704 #endif 705 667 706 ! 668 707 ! local scalars ( place there for vector optimisation purposes) 669 708 zcoef_dqlw = 4.0 * 0.95 * Stef 670 zcoef_dqla = -Ls * Cice *11637800. * (-5897.8)709 zcoef_dqla = -Ls * 11637800. * (-5897.8) 671 710 ! 672 711 zrhoa(:,:) = rho_air( sf(jp_tair)%fnow(:,:,1), sf(jp_humi)%fnow(:,:,1), sf(jp_slp)%fnow(:,:,1) ) … … 696 735 ! ... turbulent heat fluxes 697 736 ! Sensible Heat 698 z_qsb(ji,jj,jl) = zrhoa(ji,jj) * cpa * C ice* wndm_ice(ji,jj) * ( ptsu(ji,jj,jl) - sf(jp_tair)%fnow(ji,jj,1) )737 z_qsb(ji,jj,jl) = zrhoa(ji,jj) * cpa * Cd(ji,jj) * wndm_ice(ji,jj) * ( ptsu(ji,jj,jl) - sf(jp_tair)%fnow(ji,jj,1) ) 699 738 ! Latent Heat 700 qla_ice(ji,jj,jl) = rn_efac * MAX( 0.e0, zrhoa(ji,jj) * Ls * C ice* wndm_ice(ji,jj) &739 qla_ice(ji,jj,jl) = rn_efac * MAX( 0.e0, zrhoa(ji,jj) * Ls * Cd(ji,jj) * wndm_ice(ji,jj) & 701 740 & * ( 11637800. * EXP( -5897.8 / ptsu(ji,jj,jl) ) / zrhoa(ji,jj) - sf(jp_humi)%fnow(ji,jj,1) ) ) 702 741 ! Latent heat sensitivity for ice (Dqla/Dt) … … 708 747 709 748 ! Sensible heat sensitivity (Dqsb_ice/Dtn_ice) 710 z_dqsb(ji,jj,jl) = zrhoa(ji,jj) *cpa*Cice* wndm_ice(ji,jj)749 z_dqsb(ji,jj,jl) = zrhoa(ji,jj) * cpa * Cd(ji,jj) * wndm_ice(ji,jj) 711 750 712 751 ! ----------------------------! … … 782 821 CALL wrk_dealloc( jpi,jpj,jpl, z_qlw, z_qsb, z_dqlw, z_dqsb ) 783 822 CALL wrk_dealloc( jpi,jpj, zrhoa ) 823 CALL wrk_dealloc( jpi,jpj, Cd ) 784 824 ! 785 825 IF( nn_timing == 1 ) CALL timing_stop('blk_ice_flx') … … 906 946 END FUNCTION L_vap 907 947 948 949 #if defined key_lim3 950 SUBROUTINE Cdn10_Lupkes2012( Cd ) 951 !!---------------------------------------------------------------------- 952 !! *** ROUTINE Cdn10_Lupkes2012 *** 953 !! 954 !! ** Purpose : Recompute the ice-atm drag at 10m height to make 955 !! it dependent on edges at leads, melt ponds and flows. 956 !! After some approximations, this can be resumed to a dependency 957 !! on ice concentration. 958 !! 959 !! ** Method : The parameterization is taken from Lupkes et al. (2012) eq.(50) 960 !! with the highest level of approximation: level4, eq.(59) 961 !! The generic drag over a cell partly covered by ice can be re-written as follows: 962 !! 963 !! Cd = Cdw * (1-A) + Cdi * A + Ce * (1-A)**(nu+1/(10*beta)) * A**mu 964 !! 965 !! Ce = 2.23e-3 , as suggested by Lupkes (eq. 59) 966 !! nu = mu = beta = 1 , as suggested by Lupkes (eq. 59) 967 !! A is the concentration of ice minus melt ponds (if any) 968 !! 969 !! This new drag has a parabolic shape (as a function of A) starting at 970 !! Cdw(say 1.5e-3) for A=0, reaching 1.97e-3 for A~0.5 971 !! and going down to Cdi(say 1.4e-3) for A=1 972 !! 973 !! It is theoretically applicable to all ice conditions !(not only MIZ) 974 !! => see Lupkes et al (2013) 975 !! 976 !! ** References : Lupkes et al. JGR 2012 (theory) 977 !! Lupkes et al. GRL 2013 (application to GCM) 978 !! 979 !!---------------------------------------------------------------------- 980 REAL(wp), DIMENSION(:,:), INTENT(inout) :: Cd 981 REAL(wp), PARAMETER :: zCe = 2.23e-03_wp 982 REAL(wp), PARAMETER :: znu = 1._wp 983 REAL(wp), PARAMETER :: zmu = 1._wp 984 REAL(wp), PARAMETER :: zbeta = 1._wp 985 REAL(wp) :: zcoef 986 !!---------------------------------------------------------------------- 987 zcoef = znu + 1._wp / ( 10._wp * zbeta ) 988 989 ! generic drag over a cell partly covered by ice 990 !!Cd(:,:) = Cd_oce(:,:) * ( 1._wp - at_i_b(:,:) ) + & !! pure ocean drag 991 !! & Cd_ice * at_i_b(:,:) + & !! pure ice drag 992 !! & zCe * ( 1._wp - at_i_b(:,:) )**zcoef * !at_i_b(:,:)**zmu ! change due to sea-ice morphology 993 994 ! ice-atm drag 995 Cd(:,:) = Cd_ice + & ! pure ice drag 996 & zCe * ( 1._wp - at_i_b(:,:) )**zcoef * at_i_b(:,:)**(zmu-1._wp) ! change due to sea-ice morphology 997 998 END SUBROUTINE Cdn10_Lupkes2012 999 #endif 1000 1001 908 1002 !!====================================================================== 909 1003 END MODULE sbcblk
Note: See TracChangeset
for help on using the changeset viewer.