- Timestamp:
- 2016-06-30T17:17:35+02:00 (8 years ago)
- Location:
- branches/2016/dev_v3_6_STABLE_r6506_AGRIF_LIM3/NEMOGCM/NEMO/OPA_SRC/SBC
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2016/dev_v3_6_STABLE_r6506_AGRIF_LIM3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90
r6399 r6763 16 16 !! 3.4 ! 2011-11 (C. Harris) Fill arrays required by CICE 17 17 !! 3.7 ! 2014-06 (L. Brodeau) simplification and optimization of CORE bulk 18 !! 4.0 ! 2016-06 (C. Rousset) Add new param of drags with sea-ice (Lupkes at al 2012) 18 19 !!---------------------------------------------------------------------- 19 20 … … 45 46 USE lib_fortran ! to use key_nosignedzero 46 47 #if defined key_lim3 47 USE ice, ONLY : u_ice, v_ice, jpl, pfrld, a_i_b 48 USE ice, ONLY : u_ice, v_ice, jpl, pfrld, a_i_b, at_i_b 48 49 USE limthd_dh ! for CALL lim_thd_snwblow 49 50 #elif defined key_lim2 … … 60 61 PUBLIC blk_ice_core_flx ! routine called in sbc_ice_lim module 61 62 #endif 62 PUBLIC turb_core_2z ! routine calle sin sbcblk_mfs module63 PUBLIC turb_core_2z ! routine called in sbcblk_mfs module 63 64 64 65 INTEGER , PARAMETER :: jpfld = 9 ! maximum number of files to read … … 76 77 77 78 ! !!! CORE bulk parameters 78 REAL(wp), PARAMETER :: rhoa = 1.22 ! air density79 REAL(wp), PARAMETER :: cpa = 1000.5 ! specific heat of air80 REAL(wp), PARAMETER :: Lv = 2.5e6 ! latent heat of vaporization81 REAL(wp), PARAMETER :: Ls = 2.839e6 ! latent heat of sublimation82 REAL(wp), PARAMETER :: Stef = 5.67e-8 ! Stefan Boltzmann constant83 REAL(wp), PARAMETER :: C ice = 1.4e-3 ! iovi 1.63e-3! transfer coefficient over ice84 REAL(wp), PARAMETER :: albo = 0.066 ! ocean albedo assumed to be constant79 REAL(wp), PARAMETER :: rhoa = 1.22 ! air density 80 REAL(wp), PARAMETER :: cpa = 1000.5 ! specific heat of air 81 REAL(wp), PARAMETER :: Lv = 2.5e6 ! latent heat of vaporization 82 REAL(wp), PARAMETER :: Ls = 2.839e6 ! latent heat of sublimation 83 REAL(wp), PARAMETER :: Stef = 5.67e-8 ! Stefan Boltzmann constant 84 REAL(wp), PARAMETER :: Cd_ice = 1.4e-3 ! transfer coefficient over ice 85 REAL(wp), PARAMETER :: albo = 0.066 ! ocean albedo assumed to be constant 85 86 86 87 ! !!* Namelist namsbc_core : CORE bulk parameters … … 91 92 REAL(wp) :: rn_zqt ! z(q,t) : height of humidity and temperature measurements 92 93 REAL(wp) :: rn_zu ! z(u) : height of wind measurements 93 94 LOGICAL :: ln_Cd_L12 = .FALSE. ! Modify the drag ice-atm and oce-atm depending on ice concentration (from Lupkes et al. JGR2012) 95 96 ! 97 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: Cd_oce ! air-ocean drag (clem) 98 94 99 !! * Substitutions 95 100 # include "domzgr_substitute.h90" … … 102 107 CONTAINS 103 108 109 INTEGER FUNCTION sbc_blk_core_alloc() 110 !!------------------------------------------------------------------- 111 !! *** ROUTINE sbc_blk_core_alloc (clem) *** 112 !!------------------------------------------------------------------- 113 ALLOCATE( Cd_oce(jpi,jpj) , STAT=sbc_blk_core_alloc ) 114 ! 115 IF( lk_mpp ) CALL mpp_sum( sbc_blk_core_alloc ) 116 IF( sbc_blk_core_alloc /= 0 ) CALL ctl_warn('sbc_blk_core_alloc: failed to allocate arrays') 117 END FUNCTION sbc_blk_core_alloc 118 119 104 120 SUBROUTINE sbc_blk_core( kt ) 105 121 !!--------------------------------------------------------------------- … … 151 167 & sn_wndi, sn_wndj, sn_humi , sn_qsr , & 152 168 & sn_qlw , sn_tair, sn_prec , sn_snow, & 153 & sn_tdif, rn_zqt, rn_zu 169 & sn_tdif, rn_zqt, rn_zu , ln_Cd_L12 154 170 !!--------------------------------------------------------------------- 155 171 ! … … 157 173 IF( kt == nit000 ) THEN ! First call kt=nit000 ! 158 174 ! ! ====================== ! 175 ! 176 ! ! allocate sbc_blk_core array (clem) 177 IF( sbc_blk_core_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'sbc_blk_core : unable to allocate standard arrays' ) 159 178 ! 160 179 REWIND( numnam_ref ) ! Namelist namsbc_core in reference namelist : CORE bulk parameters … … 319 338 & Cd, Ch, Ce, zt_zu, zq_zu ) 320 339 340 ! Make ocean-atm. drag dependent on ice concentration (see Lupkes et al. 2012) (clem) 341 #if defined key_lim3 342 IF( ln_Cd_L12 ) THEN 343 344 Cd_oce(:,:) = Cd(:,:) ! record value of pure ocean-atm. drag 345 346 CALL Cdn10_Lupkes2012( Cd ) ! calculate new drag from Lupkes(2012) equations 347 348 ENDIF 349 #endif 350 321 351 ! ... tau module, i and j component 322 352 DO jj = 1, jpj … … 437 467 !!--------------------------------------------------------------------- 438 468 INTEGER :: ji, jj ! dummy loop indices 439 REAL(wp) :: zcoef_wnorm, zcoef_wnorm2440 469 REAL(wp) :: zwnorm_f, zwndi_f , zwndj_f ! relative wind module and components at F-point 441 470 REAL(wp) :: zwndi_t , zwndj_t ! relative wind components at T-point 471 REAL(wp), DIMENSION(:,:), POINTER :: Cd ! transfer coefficient for momentum (tau) 442 472 !!--------------------------------------------------------------------- 443 473 ! 444 474 IF( nn_timing == 1 ) CALL timing_start('blk_ice_core_tau') 445 475 ! 446 ! local scalars ( place there for vector optimisation purposes) 447 zcoef_wnorm = rhoa * Cice 448 zcoef_wnorm2 = rhoa * Cice * 0.5 476 CALL wrk_alloc( jpi,jpj, Cd ) 477 478 Cd(:,:) = Cd_ice 479 480 ! Make ice-atm. drag dependent on ice concentration (see Lupkes et al. 2012) (clem) 481 #if defined key_lim3 482 IF( ln_Cd_L12 ) THEN 483 484 CALL Cdn10_Lupkes2012( Cd ) ! calculate new drag from Lupkes(2012) equations 485 486 ENDIF 487 #endif 449 488 450 489 !!gm brutal.... … … 467 506 zwndj_f = 0.25 * ( sf(jp_wndj)%fnow(ji-1,jj ,1) + sf(jp_wndj)%fnow(ji ,jj ,1) & 468 507 & + sf(jp_wndj)%fnow(ji-1,jj-1,1) + sf(jp_wndj)%fnow(ji ,jj-1,1) ) - rn_vfac * v_ice(ji,jj) 469 zwnorm_f = zcoef_wnorm* SQRT( zwndi_f * zwndi_f + zwndj_f * zwndj_f )508 zwnorm_f = rhoa * Cd(ji,jj) * SQRT( zwndi_f * zwndi_f + zwndj_f * zwndj_f ) 470 509 ! ... ice stress at I-point 471 510 utau_ice(ji,jj) = zwnorm_f * zwndi_f … … 476 515 zwndj_t = sf(jp_wndj)%fnow(ji,jj,1) - rn_vfac * 0.25 * ( v_ice(ji,jj+1) + v_ice(ji+1,jj+1) & 477 516 & + v_ice(ji,jj ) + v_ice(ji+1,jj ) ) 478 wndm_ice(ji,jj) 517 wndm_ice(ji,jj) = SQRT( zwndi_t * zwndi_t + zwndj_t * zwndj_t ) * tmask(ji,jj,1) 479 518 END DO 480 519 END DO … … 493 532 DO jj = 2, jpjm1 494 533 DO ji = fs_2, fs_jpim1 ! vect. opt. 495 utau_ice(ji,jj) = zcoef_wnorm2* ( wndm_ice(ji+1,jj ) + wndm_ice(ji,jj) ) &534 utau_ice(ji,jj) = rhoa * Cd(ji,jj) * 0.5_wp * ( wndm_ice(ji+1,jj ) + wndm_ice(ji,jj) ) & 496 535 & * ( 0.5 * (sf(jp_wndi)%fnow(ji+1,jj,1) + sf(jp_wndi)%fnow(ji,jj,1) ) - rn_vfac * u_ice(ji,jj) ) 497 vtau_ice(ji,jj) = zcoef_wnorm2* ( wndm_ice(ji,jj+1 ) + wndm_ice(ji,jj) ) &536 vtau_ice(ji,jj) = rhoa * Cd(ji,jj) * 0.5_wp * ( wndm_ice(ji,jj+1 ) + wndm_ice(ji,jj) ) & 498 537 & * ( 0.5 * (sf(jp_wndj)%fnow(ji,jj+1,1) + sf(jp_wndj)%fnow(ji,jj,1) ) - rn_vfac * v_ice(ji,jj) ) 499 538 END DO … … 510 549 ENDIF 511 550 551 CALL wrk_dealloc( jpi,jpj, Cd ) 552 512 553 IF( nn_timing == 1 ) CALL timing_stop('blk_ice_core_tau') 513 554 … … 548 589 ! local scalars ( place there for vector optimisation purposes) 549 590 zcoef_dqlw = 4.0 * 0.95 * Stef 550 zcoef_dqla = -Ls * C ice * 11637800. * (-5897.8)551 zcoef_dqsb = rhoa * cpa * C ice591 zcoef_dqla = -Ls * Cd_ice * 11637800. * (-5897.8) 592 zcoef_dqsb = rhoa * cpa * Cd_ice 552 593 553 594 zztmp = 1. / ( 1. - albo ) … … 575 616 ! ... turbulent heat fluxes 576 617 ! Sensible Heat 577 z_qsb(ji,jj,jl) = rhoa * cpa * C ice * wndm_ice(ji,jj) * ( ptsu(ji,jj,jl) - sf(jp_tair)%fnow(ji,jj,1) )618 z_qsb(ji,jj,jl) = rhoa * cpa * Cd_ice * wndm_ice(ji,jj) * ( ptsu(ji,jj,jl) - sf(jp_tair)%fnow(ji,jj,1) ) 578 619 ! Latent Heat 579 qla_ice(ji,jj,jl) = rn_efac * MAX( 0.e0, rhoa * Ls * C ice * wndm_ice(ji,jj) &620 qla_ice(ji,jj,jl) = rn_efac * MAX( 0.e0, rhoa * Ls * Cd_ice * wndm_ice(ji,jj) & 580 621 & * ( 11637800. * EXP( -5897.8 / ptsu(ji,jj,jl) ) / rhoa - sf(jp_humi)%fnow(ji,jj,1) ) ) 581 622 ! Latent heat sensitivity for ice (Dqla/Dt) … … 820 861 ! 821 862 END DO 822 863 823 864 CALL wrk_dealloc( jpi,jpj, U_zu, Ce_n10, Ch_n10, sqrt_Cd_n10, sqrt_Cd ) 824 865 CALL wrk_dealloc( jpi,jpj, zeta_u, stab ) … … 903 944 END FUNCTION psi_h 904 945 946 947 SUBROUTINE Cdn10_Lupkes2012( Cd ) 948 !!---------------------------------------------------------------------- 949 !! *** ROUTINE Cdn10_Lupkes2012 *** 950 !! 951 !! ** Purpose : Recompute the ice-atm and ocean-atm drags at 10m height to make 952 !! them dependent on edges at leads, melt ponds and flows. 953 !! After some approximations, this can be resumed to a dependency 954 !! on ice concentration. 955 !! 956 !! ** Method : The parameterization is taken from Lupkes et al. (2012) eq.(50) 957 !! with the highest level of approximation: level4, eq.(59) 958 !! The drag can be re-written as follows: 959 !! 960 !! Cd = Cdw * (1-A) + Cdi * A + Ce * (1-A)**(nu+1/(10*beta)) * A**mu 961 !! 962 !! Ce = 2.23e-3 , as suggested by Lupkes (eq. 59) 963 !! nu = mu = beta = 1 , as suggested by Lupkes (eq. 59) 964 !! A is the concentration of ice minus melt ponds (if any) 965 !! 966 !! This new drag has a parabolic shape (as a function of A) starting at 967 !! Cdw(say 1.5e-3) for A=0, reaching 1.97e-3 for A~0.5 968 !! and going down to Cdi(say 1.4e-3) for A=1 969 !! 970 !! It is theoretically applicable to all ice conditions (not only MIZ) 971 !! => see Lupkes et al (2013) 972 !! 973 !! ** References : Lupkes et al. JGR 2012 (theory) 974 !! Lupkes et al. GRL 2013 (application to GCM) 975 !! 976 !!---------------------------------------------------------------------- 977 REAL(wp), DIMENSION(:,:), INTENT(inout) :: Cd 978 REAL(wp), PARAMETER :: zCe = 2.23e-03_wp 979 REAL(wp), PARAMETER :: znu = 1._wp 980 REAL(wp), PARAMETER :: zmu = 1._wp 981 REAL(wp), PARAMETER :: zbeta = 1._wp 982 REAL(wp) :: zcoef 983 !!---------------------------------------------------------------------- 984 zcoef = znu + 1._wp / ( 10._wp * zbeta ) 985 986 Cd(:,:) = Cd_oce(:,:) * ( 1._wp - at_i_b(:,:) ) + & ! pure ocean drag 987 & Cd_ice * at_i_b(:,:) + & ! pure ice drag 988 & zCe * ( 1._wp - at_i_b(:,:) )**zcoef * at_i_b(:,:)**zmu ! change due to sea-ice morphology 989 990 END SUBROUTINE Cdn10_Lupkes2012 991 905 992 !!====================================================================== 906 993 END MODULE sbcblk_core -
branches/2016/dev_v3_6_STABLE_r6506_AGRIF_LIM3/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90
r6746 r6763 617 617 u_ice_b(:,:) = u_ice(:,:) 618 618 v_ice_b(:,:) = v_ice(:,:) 619 at_i_b (:,:) = SUM( a_i_b(:,:,:), dim=3 ) 619 620 620 621 END SUBROUTINE sbc_lim_bef
Note: See TracChangeset
for help on using the changeset viewer.