Changeset 1705 for trunk/NEMO
- Timestamp:
- 2009-11-03T17:37:00+01:00 (14 years ago)
- Location:
- trunk/NEMO/OPA_SRC
- Files:
-
- 5 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/OPA_SRC/SBC/sbc_oce.F90
r1695 r1705 29 29 ! !: = 1 global mean of e-p-r set to zero at each nn_fsbc time step 30 30 ! !: = 2 annual global mean of e-p-r set to zero 31 INTEGER , PUBLIC :: nn_ico_cpl = 0 !: ice-ocean coupling indicator31 INTEGER , PUBLIC :: nn_ico_cpl = 0 !: ice-ocean coupling indicator 32 32 ! !: = 0 LIM-3 old case 33 33 ! !: = 1 stresses computed using now ocean velocity … … 37 37 !! Ocean Surface Boundary Condition fields 38 38 !!---------------------------------------------------------------------- 39 LOGICAL , PUBLIC :: lhftau = .FALSE. !: HF tau contribution: mean of stress module - module of the mean stress 39 40 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: utau !: sea surface i-stress (ocean referential) [N/m2] 40 41 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: vtau !: sea surface j-stress (ocean referential) [N/m2] 41 42 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: taum !: module of sea surface stress (at T-point) [N/m2] 42 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: wndm !: wind speed module at T-point (=|U10m-Uoce|) [m/s] Used only in PISCES 43 !! wndm is used only in PISCES to compute gases exchanges at the surface of the free ocean or in the leads in sea-ice parts 44 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: wndm !: wind speed module at T-point (=|U10m-Uoce|) [m/s] 43 45 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: qsr !: sea heat flux: solar [W/m2] 44 46 REAL(wp), PUBLIC, DIMENSION(jpi,jpj) :: qns !: sea heat flux: non solar [W/m2] -
trunk/NEMO/OPA_SRC/SBC/sbcblk_core.F90
r1695 r1705 43 43 PUBLIC blk_ice_core ! routine called in sbc_ice_lim module 44 44 45 INTEGER , PARAMETER :: jpfld = 8! maximum number of files to read45 INTEGER , PARAMETER :: jpfld = 9 ! maximum number of files to read 46 46 INTEGER , PARAMETER :: jp_wndi = 1 ! index of 10m wind velocity (i-component) (m/s) at T-point 47 47 INTEGER , PARAMETER :: jp_wndj = 2 ! index of 10m wind velocity (j-component) (m/s) at T-point … … 52 52 INTEGER , PARAMETER :: jp_prec = 7 ! index of total precipitation (rain+snow) (Kg/m2/s) 53 53 INTEGER , PARAMETER :: jp_snow = 8 ! index of snow (solid prcipitation) (kg/m2/s) 54 INTEGER , PARAMETER :: jp_tdif = 9 ! index of tau diff associated to HF tau (N/m2) at T-point 54 55 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf ! structure of input fields (file informations, fields read) 55 56 … … 63 64 64 65 ! !!* Namelist namsbc_core : CORE bulk parameters 65 LOGICAL :: ln_2m = .FALSE. ! logical flag for height of air temp. and hum 66 REAL(wp) :: rn_pfac = 1. ! multiplication factor for precipitation 66 LOGICAL :: ln_2m = .FALSE. ! logical flag for height of air temp. and hum 67 LOGICAL :: ln_taudif = .FALSE. ! logical flag to use the "mean of stress module - module of mean stress" data 68 REAL(wp) :: rn_pfac = 1. ! multiplication factor for precipitation 67 69 68 70 !! * Substitutions … … 93 95 !! the total precipitation (rain+snow) (Kg/m2/s) 94 96 !! the snow (solid prcipitation) (kg/m2/s) 97 !! OPTIONAL parameter (see ln_taudif namelist flag): 98 !! the tau diff associated to HF tau (N/m2) at T-point 95 99 !! (2) CALL blk_oce_core 96 100 !! … … 110 114 INTEGER :: ierror ! return error code 111 115 INTEGER :: ifpr ! dummy loop indice 116 INTEGER :: jfld ! dummy loop arguments 112 117 !! 113 118 CHARACTER(len=100) :: cn_dir ! Root directory for location of core files … … 115 120 TYPE(FLD_N) :: sn_wndi, sn_wndj, sn_humi, sn_qsr ! informations about the fields to be read 116 121 TYPE(FLD_N) :: sn_qlw , sn_tair, sn_prec, sn_snow ! " " 117 NAMELIST/namsbc_core/ cn_dir, ln_2m, rn_pfac, sn_wndi, sn_wndj, sn_humi, sn_qsr, & 118 & sn_qlw , sn_tair, sn_prec, sn_snow 122 TYPE(FLD_N) :: sn_tdif ! " " 123 NAMELIST/namsbc_core/ cn_dir , ln_2m , ln_taudif, rn_pfac, & 124 & sn_wndi, sn_wndj, sn_humi , sn_qsr , & 125 & sn_qlw , sn_tair, sn_prec , sn_snow, sn_tdif 119 126 !!--------------------------------------------------------------------- 120 127 … … 136 143 sn_prec = FLD_N( 'precip' , -1. , 'precip' , .true. , .false. , 'yearly' , '' , '' ) 137 144 sn_snow = FLD_N( 'snow' , -1. , 'snow' , .true. , .false. , 'yearly' , '' , '' ) 145 sn_tdif = FLD_N( 'taudif' , 24. , 'taudif' , .true. , .false. , 'yearly' , '' , '' ) 138 146 ! 139 147 REWIND( numnam ) ! ... read in namlist namsbc_core … … 145 153 slf_i(jp_tair) = sn_tair ; slf_i(jp_humi) = sn_humi 146 154 slf_i(jp_prec) = sn_prec ; slf_i(jp_snow) = sn_snow 155 slf_i(jp_tdif) = sn_tdif 156 ! 157 ! do we use HF tau information? 158 lhftau = ln_taudif 159 jfld = jpfld - COUNT( (/.NOT. lhftau/) ) 147 160 ! 148 161 ! set sf structure 149 ALLOCATE( sf(j pfld), STAT=ierror )162 ALLOCATE( sf(jfld), STAT=ierror ) 150 163 IF( ierror > 0 ) THEN 151 164 CALL ctl_stop( 'sbc_blk_core: unable to allocate sf structure' ) ; RETURN 152 165 ENDIF 153 DO ifpr= 1, j pfld166 DO ifpr= 1, jfld 154 167 ALLOCATE( sf(ifpr)%fnow(jpi,jpj) ) 155 168 ALLOCATE( sf(ifpr)%fdta(jpi,jpj,2) ) … … 291 304 END DO 292 305 END DO 306 307 ! ... add the HF tau contribution to the wind stress module? 308 IF( lhftau ) THEN 309 !CDIR COLLAPSE 310 taum(:,:) = taum(:,:) + sf(jp_tdif)%fnow(:,:) 311 ENDIF 312 CALL iom_put( "taum_oce", taum ) ! output wind stress module 313 293 314 ! ... utau, vtau at U- and V_points, resp. 294 315 ! Note the use of 0.5*(2-umask) in order to unmask the stress along coastlines -
trunk/NEMO/OPA_SRC/SBC/sbccpl.F90
r1698 r1705 419 419 ! ! ------------------------- ! 420 420 srcv(jpr_taum)%clname = 'O_TauMod' ; IF( TRIM(cn_rcv_taumod) == 'coupled' ) srcv(jpr_taum)%laction = .TRUE. 421 lhftau = srcv(jpr_taum)%laction 421 422 422 423 #if defined key_cpl_carbon_cycle … … 693 694 taum(:,:) = frcv(:,:,jpr_taum) 694 695 wndm(:,:) = frcv(:,:,jpr_w10m) 696 CALL iom_put( "taum_oce", taum ) ! output wind stress module 695 697 ! 696 698 ENDIF -
trunk/NEMO/OPA_SRC/SBC/sbcmod.F90
r1649 r1705 253 253 CALL iom_put( "utau", utau ) ! i-wind stress (stress can be updated at 254 254 CALL iom_put( "vtau", vtau ) ! j-wind stress each time step in sea-ice) 255 CALL iom_put( "wspd", wndm ) ! wind speed module 255 CALL iom_put( "taum", taum ) ! wind stress module 256 CALL iom_put( "wspd", wndm ) ! wind speed module 256 257 ! 257 258 IF(ln_ctl) THEN ! print mean trends (used for debugging) -
trunk/NEMO/OPA_SRC/ZDF/zdftke.F90
r1695 r1705 65 65 66 66 ! !!! ** Namelist namzdf_tke ** 67 LOGICAL :: ln_mxl0 = .FALSE. ! mixing length scale surface value as function of wind stress or not 68 INTEGER :: nn_mxl = 2 ! type of mixing length (=0/1/2/3) 69 REAL(wp) :: rn_lmin0 = 0.4_wp ! surface min value of mixing length [m] 70 REAL(wp) :: rn_lmin = 0.1_wp ! interior min value of mixing length [m] 71 INTEGER :: nn_pdl = 1 ! Prandtl number or not (ratio avt/avm) (=0/1) 72 REAL(wp) :: rn_ediff = 0.1_wp ! coefficient for avt: avt=rn_ediff*mxl*sqrt(e) 73 REAL(wp) :: rn_ediss = 0.7_wp ! coefficient of the Kolmogoroff dissipation 74 REAL(wp) :: rn_ebb = 3.75_wp ! coefficient of the surface input of tke 75 REAL(wp) :: rn_emin = 0.7071e-6_wp ! minimum value of tke [m2/s2] 76 REAL(wp) :: rn_emin0 = 1.e-4_wp ! surface minimum value of tke [m2/s2] 77 REAL(wp) :: rn_bshear= 1.e-20 ! background shear (>0) 78 INTEGER :: nn_etau = 0 ! type of depth penetration of surface tke (=0/1/2) 79 INTEGER :: nn_htau = 0 ! type of tke profile of penetration (=0/1) 80 REAL(wp) :: rn_efr = 1.0_wp ! fraction of TKE surface value which penetrates in the ocean 81 LOGICAL :: ln_lc = .FALSE. ! Langmuir cells (LC) as a source term of TKE or not 82 REAL(wp) :: rn_lc = 0.15_wp ! coef to compute vertical velocity of Langmuir cells 67 LOGICAL :: ln_mxl0 = .FALSE. ! mixing length scale surface value as function of wind stress or not 68 INTEGER :: nn_mxl = 2 ! type of mixing length (=0/1/2/3) 69 REAL(wp) :: rn_lmin0 = 0.4_wp ! surface min value of mixing length [m] 70 REAL(wp) :: rn_lmin = 0.1_wp ! interior min value of mixing length [m] 71 INTEGER :: nn_pdl = 1 ! Prandtl number or not (ratio avt/avm) (=0/1) 72 REAL(wp) :: rn_ediff = 0.1_wp ! coefficient for avt: avt=rn_ediff*mxl*sqrt(e) 73 REAL(wp) :: rn_ediss = 0.7_wp ! coefficient of the Kolmogoroff dissipation 74 REAL(wp) :: rn_ebb = 3.75_wp ! coefficient of the surface input of tke 75 REAL(wp) :: rn_emin = 0.7071e-6_wp ! minimum value of tke [m2/s2] 76 REAL(wp) :: rn_emin0 = 1.e-4_wp ! surface minimum value of tke [m2/s2] 77 REAL(wp) :: rn_bshear = 1.e-20 ! background shear (>0) 78 INTEGER :: nn_etau = 0 ! type of depth penetration of surface tke (=0/1/2) 79 INTEGER :: nn_htau = 0 ! type of tke profile of penetration (=0/1) 80 REAL(wp) :: rn_efr = 1.0_wp ! fraction of TKE surface value which penetrates in the ocean 81 REAL(wp) :: rn_addhft = 0.0_wp ! add offset applied to HF tau 82 REAL(wp) :: rn_sclhft = 1.0_wp ! scale factor applied to HF tau 83 LOGICAL :: ln_lc = .FALSE. ! Langmuir cells (LC) as a source term of TKE or not 84 REAL(wp) :: rn_lc = 0.15_wp ! coef to compute vertical velocity of Langmuir cells 83 85 84 86 REAL(wp) :: ri_cri ! critic Richardson number (deduced from rn_ediff and rn_ediss values) … … 177 179 USE oce, zd_lw => ta ! use ta as workspace 178 180 !! 179 INTEGER :: ji, jj, jk ! dummy loop arguments 180 REAL(wp) :: zbbrau, zesh2 ! temporary scalars 181 REAL(wp) :: zfact1, zfact2, zfact3 ! - - 182 REAL(wp) :: ztx2 , zty2 , zcof ! - - 183 REAL(wp) :: zus , zwlc , zind ! - - 184 REAL(wp) :: zzd_up, zzd_lw ! - - 181 INTEGER :: ji, jj, jk ! dummy loop arguments 185 182 INTEGER :: ikbu, ikbv, ikbum1, ikbvm1 ! temporary scalar 186 183 INTEGER :: ikbt, ikbumm1, ikbvmm1 ! temporary scalar 187 REAL(wp) :: zebot ! temporary scalars 184 REAL(wp) :: zrhoa = 1.22 ! Air density kg/m3 185 REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient 186 REAL(wp) :: zbbrau, zesh2 ! temporary scalars 187 REAL(wp) :: zfact1, zfact2, zfact3 ! - - 188 REAL(wp) :: ztx2 , zty2 , zcof ! - - 189 REAL(wp) :: ztau , zdif ! - - 190 REAL(wp) :: zus , zwlc , zind ! - - 191 REAL(wp) :: zzd_up, zzd_lw ! - - 192 REAL(wp) :: zebot ! - - 188 193 INTEGER , DIMENSION(jpi,jpj) :: imlc ! 2D workspace 189 194 REAL(wp), DIMENSION(jpi,jpj) :: zhlc ! - - … … 242 247 END DO 243 248 ! !* finite Langmuir Circulation depth 249 zcof = 0.5 * 0.016 * 0.016 / ( zrhoa * zcdrag ) 244 250 imlc(:,:) = mbathy(:,:) ! Initialization to the number of w ocean point mbathy 245 251 DO jk = jpkm1, 2, -1 246 252 DO jj = 1, jpj ! Last w-level at which zpelc>=0.5*us*us 247 253 DO ji = 1, jpi ! with us=0.016*wind(starting from jpk-1) 248 zus = 0.000128 * wndm(ji,jj) * wndm(ji,jj)254 zus = zcof * taum(ji,jj) 249 255 IF( zpelc(ji,jj,jk) > zus ) imlc(ji,jj) = jk 250 256 END DO … … 265 271 hlc(:,:) = zhlc(:,:) * tmask(:,:,1) ! c1d configuration: save finite Langmuir Circulation depth 266 272 # endif 273 zcof = 0.016 / SQRT( zrhoa * zcdrag ) 267 274 !CDIR NOVERRCHK 268 275 DO jk = 2, jpkm1 !* TKE Langmuir circulation source term added to en … … 271 278 !CDIR NOVERRCHK 272 279 DO ji = fs_2, fs_jpim1 ! vector opt. 273 !!gm replace here wndn by a formulation with the stress module 274 zus = 0.016 * wndm(ji,jj) ! Stokes drift 280 zus = zcof * SQRT( taum(ji,jj) ) ! Stokes drift 275 281 ! ! vertical velocity due to LC 276 282 zind = 0.5 - SIGN( 0.5, fsdepw(ji,jj,jk) - zhlc(ji,jj) ) … … 371 377 ! ! TKE due to surface and internal wave breaking 372 378 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 373 IF( nn_etau == 1 ) THEN !* 379 IF( nn_etau == 1 ) THEN !* penetration throughout the water column 374 380 DO jk = 2, jpkm1 375 381 DO jj = 2, jpjm1 … … 386 392 en(ji,jj,jk) = en(ji,jj,jk) + rn_efr * en(ji,jj,1) * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) ) & 387 393 & * ( 1.e0 - fr_i(ji,jj) ) * tmask(ji,jj,jk) 394 END DO 395 END DO 396 ELSEIF( nn_etau == 3 ) THEN !* penetration throughout the water column 397 !CDIR NOVERRCHK 398 DO jk = 2, jpkm1 399 !CDIR NOVERRCHK 400 DO jj = 2, jpjm1 401 !CDIR NOVERRCHK 402 DO ji = fs_2, fs_jpim1 ! vector opt. 403 ztx2 = utau(ji-1,jj ) + utau(ji,jj) 404 zty2 = vtau(ji ,jj-1) + vtau(ji,jj) 405 ztau = 0.5 * SQRT( ztx2 * ztx2 + zty2 * zty2 ) ! module of the mean stress 406 zdif = taum(ji,jj) - ztau ! mean of the module - module of the mean 407 zdif = rn_sclhft * MAX( 0.e0, zdif + rn_addhft ) ! apply some modifications... 408 en(ji,jj,jk) = en(ji,jj,jk) + zbbrau * zdif * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) ) & 409 & * ( 1.e0 - fr_i(ji,jj) ) * tmask(ji,jj,jk) 410 END DO 388 411 END DO 389 412 END DO … … 642 665 INTEGER :: ji, jj, jk ! dummy loop indices 643 666 !! 644 NAMELIST/namzdf_tke/ rn_ediff, rn_ediss , rn_ebb, rn_emin, & 645 & rn_emin0, rn_bshear, nn_mxl, ln_mxl0, & 646 & rn_lmin , rn_lmin0 , nn_pdl, nn_etau, & 647 & nn_htau , rn_efr , ln_lc , rn_lc 667 NAMELIST/namzdf_tke/ rn_ediff, rn_ediss , rn_ebb , rn_emin , & 668 & rn_emin0, rn_bshear, nn_mxl , ln_mxl0 , & 669 & rn_lmin , rn_lmin0 , nn_pdl , nn_etau , & 670 & nn_htau , rn_efr , rn_addhft, rn_sclhft, & 671 & ln_lc , rn_lc 648 672 !!---------------------------------------------------------------------- 649 673 … … 658 682 WRITE(numout,*) '~~~~~~~~' 659 683 WRITE(numout,*) ' Namelist namzdf_tke : set tke mixing parameters' 660 WRITE(numout,*) ' coef. to compute avt rn_ediff = ', rn_ediff 661 WRITE(numout,*) ' Kolmogoroff dissipation coef. rn_ediss = ', rn_ediss 662 WRITE(numout,*) ' tke surface input coef. rn_ebb = ', rn_ebb 663 WRITE(numout,*) ' minimum value of tke rn_emin = ', rn_emin 664 WRITE(numout,*) ' surface minimum value of tke rn_emin0 = ', rn_emin0 665 WRITE(numout,*) ' background shear (>0) rn_bshear= ', rn_bshear 666 WRITE(numout,*) ' mixing length type nn_mxl = ', nn_mxl 667 WRITE(numout,*) ' prandl number flag nn_pdl = ', nn_pdl 668 WRITE(numout,*) ' surface mixing length = F(stress) or not ln_mxl0 = ', ln_mxl0 669 WRITE(numout,*) ' surface mixing length minimum value rn_lmin0 = ', rn_lmin0 670 WRITE(numout,*) ' interior mixing length minimum value rn_lmin0 = ', rn_lmin0 671 WRITE(numout,*) ' test param. to add tke induced by wind nn_etau = ', nn_etau 672 WRITE(numout,*) ' flag for computation of exp. tke profile nn_htau = ', nn_htau 673 WRITE(numout,*) ' % of rn_emin0 which pene. the thermocline rn_efr = ', rn_efr 674 WRITE(numout,*) ' flag to take into acc. Langmuir circ. ln_lc = ', ln_lc 675 WRITE(numout,*) ' coef to compute verticla velocity of LC rn_lc = ', rn_lc 684 WRITE(numout,*) ' coef. to compute avt rn_ediff = ', rn_ediff 685 WRITE(numout,*) ' Kolmogoroff dissipation coef. rn_ediss = ', rn_ediss 686 WRITE(numout,*) ' tke surface input coef. rn_ebb = ', rn_ebb 687 WRITE(numout,*) ' minimum value of tke rn_emin = ', rn_emin 688 WRITE(numout,*) ' surface minimum value of tke rn_emin0 = ', rn_emin0 689 WRITE(numout,*) ' background shear (>0) rn_bshear = ', rn_bshear 690 WRITE(numout,*) ' mixing length type nn_mxl = ', nn_mxl 691 WRITE(numout,*) ' prandl number flag nn_pdl = ', nn_pdl 692 WRITE(numout,*) ' surface mixing length = F(stress) or not ln_mxl0 = ', ln_mxl0 693 WRITE(numout,*) ' surface mixing length minimum value rn_lmin0 = ', rn_lmin0 694 WRITE(numout,*) ' interior mixing length minimum value rn_lmin0 = ', rn_lmin0 695 WRITE(numout,*) ' test param. to add tke induced by wind nn_etau = ', nn_etau 696 WRITE(numout,*) ' flag for computation of exp. tke profile nn_htau = ', nn_htau 697 WRITE(numout,*) ' fraction of en which pene. the thermocline rn_efr = ', rn_efr 698 WRITE(numout,*) ' add offset applied to HF tau rn_addhft = ', rn_addhft 699 WRITE(numout,*) ' scale factor applied to HF tau rn_sclhft = ', rn_sclhft 700 WRITE(numout,*) ' flag to take into acc. Langmuir circ. ln_lc = ', ln_lc 701 WRITE(numout,*) ' coef to compute verticla velocity of LC rn_lc = ', rn_lc 676 702 WRITE(numout,*) 677 703 WRITE(numout,*) ' critical Richardson nb with your parameters ri_cri = ', ri_cri … … 679 705 680 706 ! !* Check of some namelist values 681 IF( nn_mxl < 0 .OR. nn_mxl > 3 ) CALL ctl_stop( 'bad flag: nn_mxl is 0, 1 or 2 ' ) 682 IF( nn_pdl < 0 .OR. nn_pdl > 1 ) CALL ctl_stop( 'bad flag: nn_pdl is 0 or 1 ' ) 683 IF( nn_htau < 0 .OR. nn_htau > 1 ) CALL ctl_stop( 'bad flag: nn_htau is 0 or 1 ' ) 684 IF( rn_lc < 0.15 .OR. rn_lc > 0.2 ) CALL ctl_stop( 'bad value: rn_lc must be between 0.15 and 0.2 ' ) 707 IF( nn_mxl < 0 .OR. nn_mxl > 3 ) CALL ctl_stop( 'bad flag: nn_mxl is 0, 1 or 2 ' ) 708 IF( nn_pdl < 0 .OR. nn_pdl > 1 ) CALL ctl_stop( 'bad flag: nn_pdl is 0 or 1 ' ) 709 IF( nn_htau < 0 .OR. nn_htau > 2 ) CALL ctl_stop( 'bad flag: nn_htau is 0, 1 or 2 ' ) 710 IF( rn_lc < 0.15 .OR. rn_lc > 0.2 ) CALL ctl_stop( 'bad value: rn_lc must be between 0.15 and 0.2 ' ) 711 IF( rn_sclhft == 0. .AND. nn_etau == 3 ) CALL ctl_stop( 'force null HF tau to penetrate the thermocline...' ) 712 IF( .NOT. lhftau .AND. nn_etau == 3 ) CALL ctl_stop( 'bad flag: nn_etau == 3 must be used with HF tau' ) 685 713 686 714 IF( nn_etau == 2 ) CALL zdf_mxl( nit000 ) ! Initialization of nmln … … 697 725 END DO 698 726 END DO 727 CASE( 2 ) ! constant depth penetration (here 30 meters) 728 htau(:,:) = 30.e0 699 729 END SELECT 700 730 ENDIF
Note: See TracChangeset
for help on using the changeset viewer.