- Timestamp:
- 2017-11-24T17:56:51+01:00 (6 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk.F90
r8637 r8813 15 15 !! 3.7 ! 2014-06 (L. Brodeau) simplification and optimization of CORE bulk 16 16 !! 4.0 ! 2016-06 (L. Brodeau) sbcblk_core becomes sbcblk and is not restricted to the CORE algorithm anymore 17 !! 17 !! ! ==> based on AeroBulk (http://aerobulk.sourceforge.net/) 18 18 !! 4.0 ! 2016-10 (G. Madec) introduce a sbc_blk_init routine 19 !! 4.0 ! 2016-10 (M. Vancoppenolle) Introduce Jules emulator (M. Vancoppenolle) 19 20 !!---------------------------------------------------------------------- 20 21 … … 23 24 !! sbc_blk : bulk formulation as ocean surface boundary condition 24 25 !! blk_oce : computes momentum, heat and freshwater fluxes over ocean 25 !! blk_ice : computes momentum, heat and freshwater fluxes over sea ice26 26 !! rho_air : density of (moist) air (depends on T_air, q_air and SLP 27 27 !! cp_air : specific heat of (moist) air (depends spec. hum. q_air) 28 28 !! q_sat : saturation humidity as a function of SLP and temperature 29 29 !! L_vap : latent heat of vaporization of water as a function of temperature 30 !! sea-ice case only : 31 !! blk_ice_tau : provide the air-ice stress 32 !! blk_ice_flx : provide the heat and mass fluxes at air-ice interface 33 !! blk_ice_qcn : provide ice surface temperature and snow/ice conduction flux (emulating JULES coupler) 34 !! Cdn10_Lupkes2012 : Lupkes et al. (2012) air-ice drag 35 !! Cdn10_Lupkes2015 : Lupkes et al. (2015) air-ice drag 30 36 !!---------------------------------------------------------------------- 31 37 USE oce ! ocean dynamics and tracers … … 42 48 USE ice , ONLY : u_ice, v_ice, jpl, a_i_b, at_i_b, tm_su 43 49 USE icethd_dh ! for CALL ice_thd_snwblow 50 USE icethd_zdf,ONLY: rn_cnd_s ! for blk_ice_qcn 44 51 #endif 45 52 USE sbcblk_algo_ncar ! => turb_ncar : NCAR - CORE (Large & Yeager, 2009) … … 63 70 PUBLIC blk_ice_tau ! routine called in icestp module 64 71 PUBLIC blk_ice_flx ! routine called in icestp module 65 #endif 72 PUBLIC blk_ice_qcn ! routine called in icestp module 73 #endif 66 74 67 75 !!Lolo: should ultimately be moved in the module with all physical constants ? … … 111 119 LOGICAL :: ln_Cd_L15 = .FALSE. ! Modify the drag ice-atm depending on ice concentration (from Lupkes et al. JGR2015) 112 120 ! 113 !!cr REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: Cd_oce ! air-ocean drag (clem)114 121 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: Cd_atm ! transfer coefficient for momentum (tau) 115 122 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: Ch_atm ! transfer coefficient for sensible heat (Q_sens) … … 139 146 !! *** ROUTINE sbc_blk_alloc *** 140 147 !!------------------------------------------------------------------- 141 !!cr ALLOCATE( Cd_oce(jpi,jpj) , STAT=sbc_blk_alloc )142 148 ALLOCATE( Cd_atm (jpi,jpj), Ch_atm (jpi,jpj), Ce_atm (jpi,jpj), t_zu(jpi,jpj), q_zu(jpi,jpj), & 143 149 & cdn_oce(jpi,jpj), chn_oce(jpi,jpj), cen_oce(jpi,jpj), STAT=sbc_blk_alloc ) … … 147 153 END FUNCTION sbc_blk_alloc 148 154 155 149 156 SUBROUTINE sbc_blk_init 150 157 !!--------------------------------------------------------------------- … … 154 161 !! 155 162 !! ** Method : 156 !!157 !! C A U T I O N : never mask the surface stress fields158 !! the stress is assumed to be in the (i,j) mesh referential159 !!160 !! ** Action :161 163 !! 162 164 !!---------------------------------------------------------------------- … … 285 287 !! 286 288 !! ** Purpose : provide at each time step the surface ocean fluxes 287 !! (momentum, heat, freshwater and runoff)289 !! (momentum, heat, freshwater and runoff) 288 290 !! 289 291 !! ** Method : (1) READ each fluxes in NetCDF files: … … 566 568 END SUBROUTINE blk_oce 567 569 570 571 572 FUNCTION rho_air( ptak, pqa, pslp ) 573 !!------------------------------------------------------------------------------- 574 !! *** FUNCTION rho_air *** 575 !! 576 !! ** Purpose : compute density of (moist) air using the eq. of state of the atmosphere 577 !! 578 !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://sourceforge.net/p/aerobulk) 579 !!------------------------------------------------------------------------------- 580 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: ptak ! air temperature [K] 581 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pqa ! air specific humidity [kg/kg] 582 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pslp ! pressure in [Pa] 583 REAL(wp), DIMENSION(jpi,jpj) :: rho_air ! density of moist air [kg/m^3] 584 !!------------------------------------------------------------------------------- 585 ! 586 rho_air = pslp / ( R_dry*ptak * ( 1._wp + rctv0*pqa ) ) 587 ! 588 END FUNCTION rho_air 589 590 591 FUNCTION cp_air( pqa ) 592 !!------------------------------------------------------------------------------- 593 !! *** FUNCTION cp_air *** 594 !! 595 !! ** Purpose : Compute specific heat (Cp) of moist air 596 !! 597 !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://sourceforge.net/p/aerobulk) 598 !!------------------------------------------------------------------------------- 599 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pqa ! air specific humidity [kg/kg] 600 REAL(wp), DIMENSION(jpi,jpj) :: cp_air ! specific heat of moist air [J/K/kg] 601 !!------------------------------------------------------------------------------- 602 ! 603 Cp_air = Cp_dry + Cp_vap * pqa 604 ! 605 END FUNCTION cp_air 606 607 608 FUNCTION q_sat( ptak, pslp ) 609 !!---------------------------------------------------------------------------------- 610 !! *** FUNCTION q_sat *** 611 !! 612 !! ** Purpose : Specific humidity at saturation in [kg/kg] 613 !! Based on accurate estimate of "e_sat" 614 !! aka saturation water vapor (Goff, 1957) 615 !! 616 !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://sourceforge.net/p/aerobulk) 617 !!---------------------------------------------------------------------------------- 618 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: ptak ! air temperature [K] 619 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pslp ! sea level atmospheric pressure [Pa] 620 REAL(wp), DIMENSION(jpi,jpj) :: q_sat ! Specific humidity at saturation [kg/kg] 621 ! 622 INTEGER :: ji, jj ! dummy loop indices 623 REAL(wp) :: ze_sat, ztmp ! local scalar 624 !!---------------------------------------------------------------------------------- 625 ! 626 DO jj = 1, jpj 627 DO ji = 1, jpi 628 ! 629 ztmp = rt0 / ptak(ji,jj) 630 ! 631 ! Vapour pressure at saturation [hPa] : WMO, (Goff, 1957) 632 ze_sat = 10.**( 10.79574*(1. - ztmp) - 5.028*LOG10(ptak(ji,jj)/rt0) & 633 & + 1.50475*10.**(-4)*(1. - 10.**(-8.2969*(ptak(ji,jj)/rt0 - 1.)) ) & 634 & + 0.42873*10.**(-3)*(10.**(4.76955*(1. - ztmp)) - 1.) + 0.78614 ) 635 ! 636 q_sat(ji,jj) = reps0 * ze_sat/( 0.01_wp*pslp(ji,jj) - (1._wp - reps0)*ze_sat ) ! 0.01 because SLP is in [Pa] 637 ! 638 END DO 639 END DO 640 ! 641 END FUNCTION q_sat 642 643 644 FUNCTION gamma_moist( ptak, pqa ) 645 !!---------------------------------------------------------------------------------- 646 !! *** FUNCTION gamma_moist *** 647 !! 648 !! ** Purpose : Compute the moist adiabatic lapse-rate. 649 !! => http://glossary.ametsoc.org/wiki/Moist-adiabatic_lapse_rate 650 !! => http://www.geog.ucsb.edu/~joel/g266_s10/lecture_notes/chapt03/oh10_3_01/oh10_3_01.html 651 !! 652 !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://sourceforge.net/p/aerobulk) 653 !!---------------------------------------------------------------------------------- 654 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: ptak ! air temperature [K] 655 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pqa ! specific humidity [kg/kg] 656 REAL(wp), DIMENSION(jpi,jpj) :: gamma_moist ! moist adiabatic lapse-rate 657 ! 658 INTEGER :: ji, jj ! dummy loop indices 659 REAL(wp) :: zrv, ziRT ! local scalar 660 !!---------------------------------------------------------------------------------- 661 ! 662 DO jj = 1, jpj 663 DO ji = 1, jpi 664 zrv = pqa(ji,jj) / (1. - pqa(ji,jj)) 665 ziRT = 1. / (R_dry*ptak(ji,jj)) ! 1/RT 666 gamma_moist(ji,jj) = grav * ( 1. + cevap*zrv*ziRT ) / ( Cp_dry + cevap*cevap*zrv*reps0*ziRT/ptak(ji,jj) ) 667 END DO 668 END DO 669 ! 670 END FUNCTION gamma_moist 671 672 673 FUNCTION L_vap( psst ) 674 !!--------------------------------------------------------------------------------- 675 !! *** FUNCTION L_vap *** 676 !! 677 !! ** Purpose : Compute the latent heat of vaporization of water from temperature 678 !! 679 !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://sourceforge.net/p/aerobulk) 680 !!---------------------------------------------------------------------------------- 681 REAL(wp), DIMENSION(jpi,jpj) :: L_vap ! latent heat of vaporization [J/kg] 682 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: psst ! water temperature [K] 683 !!---------------------------------------------------------------------------------- 684 ! 685 L_vap = ( 2.501 - 0.00237 * ( psst(:,:) - rt0) ) * 1.e6 686 ! 687 END FUNCTION L_vap 688 568 689 #if defined key_lim3 690 !!---------------------------------------------------------------------- 691 !! 'key_lim3' ESIM sea-ice model 692 !!---------------------------------------------------------------------- 693 !! blk_ice_tau : provide the air-ice stress 694 !! blk_ice_flx : provide the heat and mass fluxes at air-ice interface 695 !! blk_ice_qcn : provide ice surface temperature and snow/ice conduction flux (emulating JULES coupler) 696 !! Cdn10_Lupkes2012 : Lupkes et al. (2012) air-ice drag 697 !! Cdn10_Lupkes2015 : Lupkes et al. (2015) air-ice drag 698 !!---------------------------------------------------------------------- 569 699 570 700 SUBROUTINE blk_ice_tau … … 688 818 689 819 690 SUBROUTINE blk_ice_flx( ptsu, p alb )820 SUBROUTINE blk_ice_flx( ptsu, phs, phi, palb ) 691 821 !!--------------------------------------------------------------------- 692 822 !! *** ROUTINE blk_ice_flx *** 693 823 !! 694 !! ** Purpose : provide the surface boundary condition over sea-ice824 !! ** Purpose : provide the heat and mass fluxes at air-ice interface 695 825 !! 696 826 !! ** Method : compute heat and freshwater exchanged … … 701 831 !!--------------------------------------------------------------------- 702 832 REAL(wp), DIMENSION(:,:,:), INTENT(in) :: ptsu ! sea ice surface temperature 833 REAL(wp), DIMENSION(:,:,:), INTENT(in) :: phs ! snow thickness 834 REAL(wp), DIMENSION(:,:,:), INTENT(in) :: phi ! ice thickness 703 835 REAL(wp), DIMENSION(:,:,:), INTENT(in) :: palb ! ice albedo (all skies) 704 836 !! … … 707 839 REAL(wp) :: zcoef_dqlw, zcoef_dqla ! - - 708 840 REAL(wp) :: zztmp, z1_lsub ! - - 841 REAL(wp) :: zfrqsr_tr_i ! sea ice shortwave fraction transmitted below through the ice 842 REAL(wp) :: zfr1, zfr2, zfac ! local variables 709 843 REAL(wp), DIMENSION(jpi,jpj,jpl) :: z_qlw ! long wave heat flux over ice 710 844 REAL(wp), DIMENSION(jpi,jpj,jpl) :: z_qsb ! sensible heat flux over ice … … 717 851 IF( ln_timing ) CALL timing_start('blk_ice_flx') 718 852 ! 719 ! 720 ! local scalars 721 zcoef_dqlw = 4.0 * 0.95 * Stef 722 zcoef_dqla = -Ls * 11637800. * (-5897.8) 853 zcoef_dqlw = 4.0 * 0.95 * Stef ! local scalars 854 zcoef_dqla = -Ls * 11637800. * (-5897.8) 723 855 ! 724 856 zrhoa(:,:) = rho_air( sf(jp_tair)%fnow(:,:,1), sf(jp_humi)%fnow(:,:,1), sf(jp_slp)%fnow(:,:,1) ) … … 815 947 END DO 816 948 817 !-------------------------------------------------------------------- 818 ! FRACTIONs of net shortwave radiation which is not absorbed in the 819 ! thin surface layer and penetrates inside the ice cover 820 ! ( Maykut and Untersteiner, 1971 ; Ebert and Curry, 1993 ) 821 ! 822 fr1_i0(:,:) = ( 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ) 823 fr2_i0(:,:) = ( 0.82 * ( 1.0 - cldf_ice ) + 0.65 * cldf_ice ) 949 ! --- absorbed and transmitted shortwave radiation (W/m2) --- ! 950 ! 951 ! former coding was 952 ! fr1_i0(:,:) = ( 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ) 953 ! fr2_i0(:,:) = ( 0.82 * ( 1.0 - cldf_ice ) + 0.65 * cldf_ice ) 954 955 ! --- surface transmission parameter (i0, Grenfell Maykut 77) --- ! 956 zfr1 = ( 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ) ! standard value 957 zfr2 = ( 0.82 * ( 1.0 - cldf_ice ) + 0.65 * cldf_ice ) ! zfr2 such that zfr1 + zfr2 to equal 1 958 959 qsr_ice_tr(:,:,:) = 0._wp 960 961 DO jl = 1, jpl 962 DO jj = 1, jpj 963 DO ji = 1, jpi 964 ! 965 zfac = MAX( 0._wp , 1._wp - phi(ji,jj,jl) * 10._wp ) ! linear weighting factor: =0 for phi=0, =1 for phi = 0.1 966 zfrqsr_tr_i = zfr1 + zfac * zfr2 ! below 10 cm, linearly increase zfrqsr_tr_i until 1 at zero thickness 967 ! 968 IF ( phs(ji,jj,jl) <= 0._wp ) THEN ; zfrqsr_tr_i = zfr1 + zfac * zfr2 969 ELSE ; zfrqsr_tr_i = 0._wp ! snow fully opaque 970 ENDIF 971 ! 972 qsr_ice_tr(ji,jj,jl) = zfrqsr_tr_i * qsr_ice(ji,jj,jl) ! transmitted solar radiation 973 ! 974 END DO 975 END DO 976 END DO 824 977 ! 825 978 IF(ln_ctl) THEN … … 836 989 END SUBROUTINE blk_ice_flx 837 990 838 #endif 839 840 FUNCTION rho_air( ptak, pqa, pslp ) 841 !!------------------------------------------------------------------------------- 842 !! *** FUNCTION rho_air *** 843 !! 844 !! ** Purpose : compute density of (moist) air using the eq. of state of the atmosphere 845 !! 846 !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://sourceforge.net/p/aerobulk) 847 !!------------------------------------------------------------------------------- 848 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: ptak ! air temperature [K] 849 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pqa ! air specific humidity [kg/kg] 850 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pslp ! pressure in [Pa] 851 REAL(wp), DIMENSION(jpi,jpj) :: rho_air ! density of moist air [kg/m^3] 852 !!------------------------------------------------------------------------------- 853 ! 854 rho_air = pslp / ( R_dry*ptak * ( 1._wp + rctv0*pqa ) ) 855 ! 856 END FUNCTION rho_air 857 858 859 FUNCTION cp_air( pqa ) 860 !!------------------------------------------------------------------------------- 861 !! *** FUNCTION cp_air *** 862 !! 863 !! ** Purpose : Compute specific heat (Cp) of moist air 864 !! 865 !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://sourceforge.net/p/aerobulk) 866 !!------------------------------------------------------------------------------- 867 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pqa ! air specific humidity [kg/kg] 868 REAL(wp), DIMENSION(jpi,jpj) :: cp_air ! specific heat of moist air [J/K/kg] 869 !!------------------------------------------------------------------------------- 870 ! 871 Cp_air = Cp_dry + Cp_vap * pqa 872 ! 873 END FUNCTION cp_air 874 875 876 FUNCTION q_sat( ptak, pslp ) 877 !!---------------------------------------------------------------------------------- 878 !! *** FUNCTION q_sat *** 879 !! 880 !! ** Purpose : Specific humidity at saturation in [kg/kg] 881 !! Based on accurate estimate of "e_sat" 882 !! aka saturation water vapor (Goff, 1957) 883 !! 884 !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://sourceforge.net/p/aerobulk) 885 !!---------------------------------------------------------------------------------- 886 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: ptak ! air temperature [K] 887 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pslp ! sea level atmospheric pressure [Pa] 888 REAL(wp), DIMENSION(jpi,jpj) :: q_sat ! Specific humidity at saturation [kg/kg] 889 ! 890 INTEGER :: ji, jj ! dummy loop indices 891 REAL(wp) :: ze_sat, ztmp ! local scalar 892 !!---------------------------------------------------------------------------------- 893 ! 894 DO jj = 1, jpj 895 DO ji = 1, jpi 896 ! 897 ztmp = rt0 / ptak(ji,jj) 898 ! 899 ! Vapour pressure at saturation [hPa] : WMO, (Goff, 1957) 900 ze_sat = 10.**( 10.79574*(1. - ztmp) - 5.028*LOG10(ptak(ji,jj)/rt0) & 901 & + 1.50475*10.**(-4)*(1. - 10.**(-8.2969*(ptak(ji,jj)/rt0 - 1.)) ) & 902 & + 0.42873*10.**(-3)*(10.**(4.76955*(1. - ztmp)) - 1.) + 0.78614 ) 991 992 SUBROUTINE blk_ice_qcn( k_monocat, ptsu, ptb, phs, phi ) 993 !!--------------------------------------------------------------------- 994 !! *** ROUTINE blk_ice_qcn *** 995 !! 996 !! ** Purpose : Compute surface temperature and snow/ice conduction flux 997 !! to force sea ice / snow thermodynamics 998 !! in the case JULES coupler is emulated 999 !! 1000 !! ** Method : compute surface energy balance assuming neglecting heat storage 1001 !! following the 0-layer Semtner (1976) approach 1002 !! 1003 !! ** Outputs : - ptsu : sea-ice / snow surface temperature (K) 1004 !! - qcn_ice : surface inner conduction flux (W/m2) 1005 !! 1006 !!--------------------------------------------------------------------- 1007 INTEGER , INTENT(in ) :: k_monocat ! single-category option 1008 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: ptsu ! sea ice / snow surface temperature 1009 REAL(wp), DIMENSION(:,:) , INTENT(in ) :: ptb ! sea ice base temperature 1010 REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: phs ! snow thickness 1011 REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: phi ! sea ice thickness 1012 ! 1013 INTEGER , PARAMETER :: nit = 10 ! number of iterations 1014 REAL(wp), PARAMETER :: zepsilon = 0.1_wp ! characteristic thickness for enhanced conduction 1015 ! 1016 INTEGER :: ji, jj, jl ! dummy loop indices 1017 INTEGER :: iter ! local integer 1018 REAL(wp) :: zfac, zfac2, zfac3 ! local scalars 1019 REAL(wp) :: zkeff_h, ztsu ! 1020 REAL(wp) :: zqc, zqnet ! 1021 REAL(wp) :: zhe, zqa0 ! 1022 REAL(wp), DIMENSION(jpi,jpj,jpl) :: zgfac ! enhanced conduction factor 1023 !!--------------------------------------------------------------------- 1024 1025 IF( nn_timing == 1 ) CALL timing_start('blk_ice_qcn') 1026 1027 ! -------------------------------------! 1028 ! I Enhanced conduction factor ! 1029 ! -------------------------------------! 1030 ! Emulates the enhancement of conduction by unresolved thin ice (k_monocat = 1/3) 1031 ! Fichefet and Morales Maqueda, JGR 1997 1032 ! 1033 zgfac(:,:,:) = 1._wp 1034 1035 SELECT CASE ( k_monocat ) 1036 ! 1037 CASE ( 1 , 3 ) 1038 ! 1039 zfac = 1._wp / ( rn_cnd_s + rcdic ) 1040 zfac2 = EXP(1._wp) * 0.5_wp * zepsilon 1041 zfac3 = 2._wp / zepsilon 1042 ! 1043 DO jl = 1, jpl 1044 DO jj = 1 , jpj 1045 DO ji = 1, jpi 1046 ! ! Effective thickness 1047 zhe = ( rn_cnd_s * phi(ji,jj,jl) + rcdic * phs(ji,jj,jl) ) * zfac 1048 ! ! Enhanced conduction factor 1049 IF( zhe >= zfac2 ) & 1050 zgfac(ji,jj,jl) = MIN( 2._wp, ( 0.5_wp + 0.5 * LOG( zhe * zfac3 ) ) ) 1051 END DO 1052 END DO 1053 END DO 1054 ! 1055 END SELECT 1056 1057 ! -------------------------------------------------------------! 1058 ! II Surface temperature and conduction flux ! 1059 ! -------------------------------------------------------------! 1060 ! 1061 zfac = rcdic * rn_cnd_s 1062 ! ! ========================== ! 1063 DO jl = 1, jpl ! Loop over ice categories ! 1064 ! ! ========================== ! 1065 DO jj = 1 , jpj 1066 DO ji = 1, jpi 1067 ! ! Effective conductivity of the snow-ice system divided by thickness 1068 zkeff_h = zfac * zgfac(ji,jj,jl) / ( rcdic * phs(ji,jj,jl) + rn_cnd_s * phi(ji,jj,jl) ) 1069 ! ! Store initial surface temperature 1070 ztsu = ptsu(ji,jj,jl) 1071 ! ! Net initial atmospheric heat flux 1072 zqa0 = qsr_ice(ji,jj,jl) - qsr_ice_tr(ji,jj,jl) + qns_ice(ji,jj,jl) 903 1073 ! 904 q_sat(ji,jj) = reps0 * ze_sat/( 0.01_wp*pslp(ji,jj) - (1._wp - reps0)*ze_sat ) ! 0.01 because SLP is in [Pa] 905 ! 906 END DO 907 END DO 908 ! 909 END FUNCTION q_sat 910 911 912 FUNCTION gamma_moist( ptak, pqa ) 913 !!---------------------------------------------------------------------------------- 914 !! *** FUNCTION gamma_moist *** 915 !! 916 !! ** Purpose : Compute the moist adiabatic lapse-rate. 917 !! => http://glossary.ametsoc.org/wiki/Moist-adiabatic_lapse_rate 918 !! => http://www.geog.ucsb.edu/~joel/g266_s10/lecture_notes/chapt03/oh10_3_01/oh10_3_01.html 919 !! 920 !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://sourceforge.net/p/aerobulk) 921 !!---------------------------------------------------------------------------------- 922 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: ptak ! air temperature [K] 923 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pqa ! specific humidity [kg/kg] 924 REAL(wp), DIMENSION(jpi,jpj) :: gamma_moist ! moist adiabatic lapse-rate 925 ! 926 INTEGER :: ji, jj ! dummy loop indices 927 REAL(wp) :: zrv, ziRT ! local scalar 928 !!---------------------------------------------------------------------------------- 929 ! 930 DO jj = 1, jpj 931 DO ji = 1, jpi 932 zrv = pqa(ji,jj) / (1. - pqa(ji,jj)) 933 ziRT = 1. / (R_dry*ptak(ji,jj)) ! 1/RT 934 gamma_moist(ji,jj) = grav * ( 1. + cevap*zrv*ziRT ) / ( Cp_dry + cevap*cevap*zrv*reps0*ziRT/ptak(ji,jj) ) 935 END DO 936 END DO 937 ! 938 END FUNCTION gamma_moist 939 940 941 FUNCTION L_vap( psst ) 942 !!--------------------------------------------------------------------------------- 943 !! *** FUNCTION L_vap *** 944 !! 945 !! ** Purpose : Compute the latent heat of vaporization of water from temperature 946 !! 947 !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://sourceforge.net/p/aerobulk) 948 !!---------------------------------------------------------------------------------- 949 REAL(wp), DIMENSION(jpi,jpj) :: L_vap ! latent heat of vaporization [J/kg] 950 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: psst ! water temperature [K] 951 !!---------------------------------------------------------------------------------- 952 ! 953 L_vap = ( 2.501 - 0.00237 * ( psst(:,:) - rt0) ) * 1.e6 954 ! 955 END FUNCTION L_vap 956 957 #if defined key_lim3 1074 DO iter = 1, nit ! --- Iteration loop 1075 ! ! Conduction heat flux through snow-ice system (>0 downwards) 1076 zqc = zkeff_h * ( ztsu - ptb(ji,jj) ) 1077 ! ! Surface energy budget 1078 zqnet = zqa0 + dqns_ice(ji,jj,jl) * ( ztsu - ptsu(ji,jj,jl) ) - zqc 1079 ! ! Temperature update 1080 ztsu = ztsu - zqnet / ( dqns_ice(ji,jj,jl) - zkeff_h ) 1081 END DO 1082 ! 1083 ptsu (ji,jj,jl) = MIN( rt0, ztsu ) 1084 qcn_ice(ji,jj,jl) = zkeff_h * ( ptsu(ji,jj,jl) - ptb(ji,jj) ) 1085 ! 1086 END DO 1087 END DO 1088 ! 1089 END DO 1090 ! 1091 IF( nn_timing == 1 ) CALL timing_stop('blk_ice_qcn') 1092 ! 1093 END SUBROUTINE blk_ice_qcn 1094 958 1095 959 1096 SUBROUTINE Cdn10_Lupkes2012( Cd ) … … 961 1098 !! *** ROUTINE Cdn10_Lupkes2012 *** 962 1099 !! 963 !! ** Purpose : Recompute the ice-atm drag at 10m height to make964 !! it dependent on edges at leads, melt ponds and flows.1100 !! ** Purpose : Recompute the neutral air-ice drag referenced at 10m 1101 !! to make it dependent on edges at leads, melt ponds and flows. 965 1102 !! After some approximations, this can be resumed to a dependency 966 1103 !! on ice concentration. … … 1012 1149 !! *** ROUTINE Cdn10_Lupkes2015 *** 1013 1150 !! 1014 !! ** pUrpose : 1lternative turbulent transfert coefficients formulation1151 !! ** pUrpose : Alternative turbulent transfert coefficients formulation 1015 1152 !! between sea-ice and atmosphere with distinct momentum 1016 1153 !! and heat coefficients depending on sea-ice concentration
Note: See TracChangeset
for help on using the changeset viewer.