Changeset 2375 for branches/nemo_v3_3_beta
- Timestamp:
- 2010-11-11T15:32:35+01:00 (13 years ago)
- Location:
- branches/nemo_v3_3_beta/NEMOGCM
- Files:
-
- 7 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/nemo_v3_3_beta/NEMOGCM/CONFIG/GYRE/EXP00/namelist
r2371 r2375 544 544 rn_ediff = 0.1 ! coef. for vertical eddy coef. (avt=rn_ediff*mxl*sqrt(e) ) 545 545 rn_ediss = 0.7 ! coef. of the Kolmogoroff dissipation 546 rn_ebb = 3.75 ! coef. of the surface input of tke547 rn_emin = 1.e- 5! minimum value of tke [m2/s2]546 rn_ebb = 67.83 ! coef. of the surface input of tke (=67.83 suggested when ln_mxl0=T) 547 rn_emin = 1.e-6 ! minimum value of tke [m2/s2] 548 548 rn_emin0 = 1.e-4 ! surface minimum value of tke [m2/s2] 549 rn_bshear = 1.e-20 ! background shear (>0)550 549 nn_mxl = 2 ! mixing length: = 0 bounded by the distance to surface and bottom 551 550 ! = 1 bounded by the local vertical scale factor 552 551 ! = 2 first vertical derivative of mixing length bounded by 1 553 ! = 3 same criteria as case 2 but applied in a different way552 ! = 3 as =2 with distinct disspipative an mixing length scale 554 553 nn_pdl = 1 ! Prandtl number function of richarson number (=1, avt=pdl(Ri)*avm) or not (=0, avt=avm) 555 ln_mxl0 = .false. ! mixing length scale surface value as function of wind stress (T) or not (F) 556 rn_lmin = 0.4 ! interior buoyancy lenght scale minimum value 557 rn_lmin0 = 0.4 ! surface buoyancy lenght scale minimum value 558 nn_etau = 0 ! exponentially deceasing penetration of tke due to internal & intertial waves 559 ! = 0 no penetration ( O(2 km) resolution) 560 ! = 1 additional tke source (rn_efr * en) 561 ! = 2 additional tke source (rn_efr * en) applied only at the base of the mixed layer 562 ! = 3 additional tke source (HF contribution: mean of stress module - module of mean stress) 563 nn_htau = 1 ! type of exponential decrease of tke penetration 554 ln_mxl0 = .true. ! surface mixing length scale = F(wind stress) (T) or not (F) 555 rn_mxl0 = 0.04 ! surface buoyancy lenght scale minimum value 556 ln_lc = .true. ! Langmuir cell parameterisation (Axell 2002) 557 rn_lc = 0.15 ! coef. associated to Langmuir cells 558 nn_etau = 0 ! penetration of tke below the mixed layer (ML) due to internal & intertial waves 559 ! = 0 no penetration 560 ! = 1 add a tke source below the ML 561 ! = 2 add a tke source just at the base of the ML 562 ! = 3 as = 1 applied on HF part of the stress ("key_coupled") 563 rn_efr = 0.05 ! fraction of surface tke value which penetrates below the ML (nn_etau=1 or 2) 564 nn_htau = 1 ! type of exponential decrease of tke penetration below the ML 564 565 ! = 0 constant 10 m length scale 565 ! = 1 0.5m at the equator to 30m at high latitudes 566 ! = 2 30 meters constant depth penetration 567 ! otion used only if nn_etau = 1 or 2: 568 rn_efr = 0.05 ! fraction of surface tke value which penetrates inside the ocean 569 ! otion used only if nn_etau = 3: 570 rn_addhft = -1.e-3 ! add offset applied to the "mean of stress module - module of mean stress" (always kept > 0) 571 rn_sclhft = 1. ! scale factor applied to the "mean of stress module - module of mean stress" 572 ln_lc = .false. ! Langmuir cell parameterisation 573 rn_lc = 0.15 ! coef. associated to Langmuir cells 566 ! = 1 0.5m at the equator to 30m poleward of 40 degrees 574 567 / 575 568 !------------------------------------------------------------------------ -
branches/nemo_v3_3_beta/NEMOGCM/CONFIG/GYRE_LOBSTER/EXP00/namelist
r2371 r2375 544 544 rn_ediff = 0.1 ! coef. for vertical eddy coef. (avt=rn_ediff*mxl*sqrt(e) ) 545 545 rn_ediss = 0.7 ! coef. of the Kolmogoroff dissipation 546 rn_ebb = 3.75 ! coef. of the surface input of tke547 rn_emin = 1.e- 5! minimum value of tke [m2/s2]546 rn_ebb = 67.83 ! coef. of the surface input of tke (=67.83 suggested when ln_mxl0=T) 547 rn_emin = 1.e-6 ! minimum value of tke [m2/s2] 548 548 rn_emin0 = 1.e-4 ! surface minimum value of tke [m2/s2] 549 rn_bshear = 1.e-20 ! background shear (>0)550 549 nn_mxl = 2 ! mixing length: = 0 bounded by the distance to surface and bottom 551 550 ! = 1 bounded by the local vertical scale factor 552 551 ! = 2 first vertical derivative of mixing length bounded by 1 553 ! = 3 same criteria as case 2 but applied in a different way552 ! = 3 as =2 with distinct disspipative an mixing length scale 554 553 nn_pdl = 1 ! Prandtl number function of richarson number (=1, avt=pdl(Ri)*avm) or not (=0, avt=avm) 555 ln_mxl0 = .false. ! mixing length scale surface value as function of wind stress (T) or not (F) 556 rn_lmin = 0.4 ! interior buoyancy lenght scale minimum value 557 rn_lmin0 = 0.4 ! surface buoyancy lenght scale minimum value 558 nn_etau = 0 ! exponentially deceasing penetration of tke due to internal & intertial waves 559 ! = 0 no penetration ( O(2 km) resolution) 560 ! = 1 additional tke source (rn_efr * en) 561 ! = 2 additional tke source (rn_efr * en) applied only at the base of the mixed layer 562 ! = 3 additional tke source (HF contribution: mean of stress module - module of mean stress) 563 nn_htau = 1 ! type of exponential decrease of tke penetration 554 ln_mxl0 = .true. ! surface mixing length scale = F(wind stress) (T) or not (F) 555 rn_mxl0 = 0.04 ! surface buoyancy lenght scale minimum value 556 ln_lc = .true. ! Langmuir cell parameterisation (Axell 2002) 557 rn_lc = 0.15 ! coef. associated to Langmuir cells 558 nn_etau = 0 ! penetration of tke below the mixed layer (ML) due to internal & intertial waves 559 ! = 0 no penetration 560 ! = 1 add a tke source below the ML 561 ! = 2 add a tke source just at the base of the ML 562 ! = 3 as = 1 applied on HF part of the stress ("key_coupled") 563 rn_efr = 0.05 ! fraction of surface tke value which penetrates below the ML (nn_etau=1 or 2) 564 nn_htau = 1 ! type of exponential decrease of tke penetration below the ML 564 565 ! = 0 constant 10 m length scale 565 ! = 1 0.5m at the equator to 30m at high latitudes 566 ! = 2 30 meters constant depth penetration 567 ! otion used only if nn_etau = 1 or 2: 568 rn_efr = 0.05 ! fraction of surface tke value which penetrates inside the ocean 569 ! otion used only if nn_etau = 3: 570 rn_addhft = -1.e-3 ! add offset applied to the "mean of stress module - module of mean stress" (always kept > 0) 571 rn_sclhft = 1. ! scale factor applied to the "mean of stress module - module of mean stress" 572 ln_lc = .false. ! Langmuir cell parameterisation 573 rn_lc = 0.15 ! coef. associated to Langmuir cells 566 ! = 1 0.5m at the equator to 30m poleward of 40 degrees 574 567 / 575 568 !------------------------------------------------------------------------ -
branches/nemo_v3_3_beta/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/1_namelist
r2371 r2375 552 552 rn_ediff = 0.1 ! coef. for vertical eddy coef. (avt=rn_ediff*mxl*sqrt(e) ) 553 553 rn_ediss = 0.7 ! coef. of the Kolmogoroff dissipation 554 rn_ebb = 6 0. ! coef. of the surface input of tke554 rn_ebb = 67.83 ! coef. of the surface input of tke (=67.83 suggested when ln_mxl0=T) 555 555 rn_emin = 1.e-6 ! minimum value of tke [m2/s2] 556 556 rn_emin0 = 1.e-4 ! surface minimum value of tke [m2/s2] 557 rn_bshear = 1.e-20 ! background shear (>0)558 557 nn_mxl = 2 ! mixing length: = 0 bounded by the distance to surface and bottom 559 558 ! = 1 bounded by the local vertical scale factor 560 559 ! = 2 first vertical derivative of mixing length bounded by 1 561 ! = 3 same criteria as case 2 but applied in a different way560 ! = 3 as =2 with distinct disspipative an mixing length scale 562 561 nn_pdl = 1 ! Prandtl number function of richarson number (=1, avt=pdl(Ri)*avm) or not (=0, avt=avm) 563 ln_mxl0 = .false. ! mixing length scale surface value as function of wind stress (T) or not (F) 564 rn_lmin = 0.001 ! interior buoyancy lenght scale minimum value 565 rn_lmin0 = 0.01 ! surface buoyancy lenght scale minimum value 566 nn_etau = 0 ! exponentially deceasing penetration of tke due to internal & intertial waves 567 ! = 0 no penetration ( O(2 km) resolution) 568 ! = 1 additional tke source (rn_efr * en) 569 ! = 2 additional tke source (rn_efr * en) applied only at the base of the mixed layer 570 ! = 3 additional tke source (HF contribution: mean of stress module - module of mean stress) 571 nn_htau = 1 ! type of exponential decrease of tke penetration 562 ln_mxl0 = .true. ! surface mixing length scale = F(wind stress) (T) or not (F) 563 rn_mxl0 = 0.04 ! surface buoyancy lenght scale minimum value 564 ln_lc = .true. ! Langmuir cell parameterisation (Axell 2002) 565 rn_lc = 0.15 ! coef. associated to Langmuir cells 566 nn_etau = 1 ! penetration of tke below the mixed layer (ML) due to internal & intertial waves 567 ! = 0 no penetration 568 ! = 1 add a tke source below the ML 569 ! = 2 add a tke source just at the base of the ML 570 ! = 3 as = 1 applied on HF part of the stress ("key_coupled") 571 rn_efr = 0.05 ! fraction of surface tke value which penetrates below the ML (nn_etau=1 or 2) 572 nn_htau = 1 ! type of exponential decrease of tke penetration below the ML 572 573 ! = 0 constant 10 m length scale 573 ! = 1 0.5m at the equator to 30m at high latitudes 574 ! = 2 30 meters constant depth penetration 575 ! otion used only if nn_etau = 1 or 2: 576 rn_efr = 0.05 ! fraction of surface tke value which penetrates inside the ocean 577 ! otion used only if nn_etau = 3: 578 rn_addhft = -1.e-3 ! add offset applied to the "mean of stress module - module of mean stress" (always kept > 0) 579 rn_sclhft = 1. ! scale factor applied to the "mean of stress module - module of mean stress" 580 ln_lc = .false. ! Langmuir cell parameterisation 581 rn_lc = 0.15 ! coef. associated to Langmuir cells 574 ! = 1 0.5m at the equator to 30m poleward of 40 degrees 582 575 / 583 576 !------------------------------------------------------------------------ -
branches/nemo_v3_3_beta/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/namelist
r2371 r2375 589 589 rn_ediff = 0.1 ! coef. for vertical eddy coef. (avt=rn_ediff*mxl*sqrt(e) ) 590 590 rn_ediss = 0.7 ! coef. of the Kolmogoroff dissipation 591 rn_ebb = 6 0. ! coef. of the surface input of tke591 rn_ebb = 67.83 ! coef. of the surface input of tke (=67.83 suggested when ln_mxl0=T) 592 592 rn_emin = 1.e-6 ! minimum value of tke [m2/s2] 593 593 rn_emin0 = 1.e-4 ! surface minimum value of tke [m2/s2] 594 rn_bshear = 1.e-20 ! background shear (>0)595 594 nn_mxl = 2 ! mixing length: = 0 bounded by the distance to surface and bottom 596 595 ! = 1 bounded by the local vertical scale factor 597 596 ! = 2 first vertical derivative of mixing length bounded by 1 598 ! = 3 same criteria as case 2 but applied in a different way597 ! = 3 as =2 with distinct disspipative an mixing length scale 599 598 nn_pdl = 1 ! Prandtl number function of richarson number (=1, avt=pdl(Ri)*avm) or not (=0, avt=avm) 600 ln_mxl0 = .false. ! mixing length scale surface value as function of wind stress (T) or not (F) 601 rn_lmin = 0.001 ! interior buoyancy lenght scale minimum value 602 rn_lmin0 = 0.01 ! surface buoyancy lenght scale minimum value 603 nn_etau = 0 ! exponentially deceasing penetration of tke due to internal & intertial waves 604 ! = 0 no penetration ( O(2 km) resolution) 605 ! = 1 additional tke source (rn_efr * en) 606 ! = 2 additional tke source (rn_efr * en) applied only at the base of the mixed layer 607 ! = 3 additional tke source (HF contribution: mean of stress module - module of mean stress) 608 nn_htau = 1 ! type of exponential decrease of tke penetration 599 ln_mxl0 = .true. ! surface mixing length scale = F(wind stress) (T) or not (F) 600 rn_mxl0 = 0.04 ! surface buoyancy lenght scale minimum value 601 ln_lc = .true. ! Langmuir cell parameterisation (Axell 2002) 602 rn_lc = 0.15 ! coef. associated to Langmuir cells 603 nn_etau = 1 ! penetration of tke below the mixed layer (ML) due to internal & intertial waves 604 ! = 0 no penetration 605 ! = 1 add a tke source below the ML 606 ! = 2 add a tke source just at the base of the ML 607 ! = 3 as = 1 applied on HF part of the stress ("key_coupled") 608 rn_efr = 0.05 ! fraction of surface tke value which penetrates below the ML (nn_etau=1 or 2) 609 nn_htau = 1 ! type of exponential decrease of tke penetration below the ML 609 610 ! = 0 constant 10 m length scale 610 ! = 1 0.5m at the equator to 30m at high latitudes 611 ! = 2 30 meters constant depth penetration 612 ! otion used only if nn_etau = 1 or 2: 613 rn_efr = 0.05 ! fraction of surface tke value which penetrates inside the ocean 614 ! otion used only if nn_etau = 3: 615 rn_addhft = -1.e-3 ! add offset applied to the "mean of stress module - module of mean stress" (always kept > 0) 616 rn_sclhft = 1. ! scale factor applied to the "mean of stress module - module of mean stress" 617 ln_lc = .false. ! Langmuir cell parameterisation 618 rn_lc = 0.15 ! coef. associated to Langmuir cells 611 ! = 1 0.5m at the equator to 30m poleward of 40 degrees 619 612 / 620 613 !------------------------------------------------------------------------ -
branches/nemo_v3_3_beta/NEMOGCM/CONFIG/ORCA2_LIM_PISCES/EXP00/namelist
r2371 r2375 589 589 rn_ediff = 0.1 ! coef. for vertical eddy coef. (avt=rn_ediff*mxl*sqrt(e) ) 590 590 rn_ediss = 0.7 ! coef. of the Kolmogoroff dissipation 591 rn_ebb = 6 0. ! coef. of the surface input of tke591 rn_ebb = 67.83 ! coef. of the surface input of tke (=67.83 suggested when ln_mxl0=T) 592 592 rn_emin = 1.e-6 ! minimum value of tke [m2/s2] 593 593 rn_emin0 = 1.e-4 ! surface minimum value of tke [m2/s2] 594 rn_bshear = 1.e-20 ! background shear (>0)595 594 nn_mxl = 2 ! mixing length: = 0 bounded by the distance to surface and bottom 596 595 ! = 1 bounded by the local vertical scale factor 597 596 ! = 2 first vertical derivative of mixing length bounded by 1 598 ! = 3 same criteria as case 2 but applied in a different way597 ! = 3 as =2 with distinct disspipative an mixing length scale 599 598 nn_pdl = 1 ! Prandtl number function of richarson number (=1, avt=pdl(Ri)*avm) or not (=0, avt=avm) 600 ln_mxl0 = .false. ! mixing length scale surface value as function of wind stress (T) or not (F) 601 rn_lmin = 0.001 ! interior buoyancy lenght scale minimum value 602 rn_lmin0 = 0.01 ! surface buoyancy lenght scale minimum value 603 nn_etau = 0 ! exponentially deceasing penetration of tke due to internal & intertial waves 604 ! = 0 no penetration ( O(2 km) resolution) 605 ! = 1 additional tke source (rn_efr * en) 606 ! = 2 additional tke source (rn_efr * en) applied only at the base of the mixed layer 607 ! = 3 additional tke source (HF contribution: mean of stress module - module of mean stress) 608 nn_htau = 1 ! type of exponential decrease of tke penetration 599 ln_mxl0 = .true. ! surface mixing length scale = F(wind stress) (T) or not (F) 600 rn_mxl0 = 0.04 ! surface buoyancy lenght scale minimum value 601 ln_lc = .true. ! Langmuir cell parameterisation (Axell 2002) 602 rn_lc = 0.15 ! coef. associated to Langmuir cells 603 nn_etau = 1 ! penetration of tke below the mixed layer (ML) due to internal & intertial waves 604 ! = 0 no penetration 605 ! = 1 add a tke source below the ML 606 ! = 2 add a tke source just at the base of the ML 607 ! = 3 as = 1 applied on HF part of the stress ("key_coupled") 608 rn_efr = 0.05 ! fraction of surface tke value which penetrates below the ML (nn_etau=1 or 2) 609 nn_htau = 1 ! type of exponential decrease of tke penetration below the ML 609 610 ! = 0 constant 10 m length scale 610 ! = 1 0.5m at the equator to 30m at high latitudes 611 ! = 2 30 meters constant depth penetration 612 ! otion used only if nn_etau = 1 or 2: 613 rn_efr = 0.05 ! fraction of surface tke value which penetrates inside the ocean 614 ! otion used only if nn_etau = 3: 615 rn_addhft = -1.e-3 ! add offset applied to the "mean of stress module - module of mean stress" (always kept > 0) 616 rn_sclhft = 1. ! scale factor applied to the "mean of stress module - module of mean stress" 617 ln_lc = .false. ! Langmuir cell parameterisation 618 rn_lc = 0.15 ! coef. associated to Langmuir cells 611 ! = 1 0.5m at the equator to 30m poleward of 40 degrees 619 612 / 620 613 !------------------------------------------------------------------------ -
branches/nemo_v3_3_beta/NEMOGCM/CONFIG/ORCA2_OFF_PISCES/EXP00/namelist
r2371 r2375 589 589 rn_ediff = 0.1 ! coef. for vertical eddy coef. (avt=rn_ediff*mxl*sqrt(e) ) 590 590 rn_ediss = 0.7 ! coef. of the Kolmogoroff dissipation 591 rn_ebb = 6 0. ! coef. of the surface input of tke591 rn_ebb = 67.83 ! coef. of the surface input of tke (=67.83 suggested when ln_mxl0=T) 592 592 rn_emin = 1.e-6 ! minimum value of tke [m2/s2] 593 593 rn_emin0 = 1.e-4 ! surface minimum value of tke [m2/s2] 594 rn_bshear = 1.e-20 ! background shear (>0)595 594 nn_mxl = 2 ! mixing length: = 0 bounded by the distance to surface and bottom 596 595 ! = 1 bounded by the local vertical scale factor 597 596 ! = 2 first vertical derivative of mixing length bounded by 1 598 ! = 3 same criteria as case 2 but applied in a different way597 ! = 3 as =2 with distinct disspipative an mixing length scale 599 598 nn_pdl = 1 ! Prandtl number function of richarson number (=1, avt=pdl(Ri)*avm) or not (=0, avt=avm) 600 ln_mxl0 = .false. ! mixing length scale surface value as function of wind stress (T) or not (F) 601 rn_lmin = 0.001 ! interior buoyancy lenght scale minimum value 602 rn_lmin0 = 0.01 ! surface buoyancy lenght scale minimum value 603 nn_etau = 0 ! exponentially deceasing penetration of tke due to internal & intertial waves 604 ! = 0 no penetration ( O(2 km) resolution) 605 ! = 1 additional tke source (rn_efr * en) 606 ! = 2 additional tke source (rn_efr * en) applied only at the base of the mixed layer 607 ! = 3 additional tke source (HF contribution: mean of stress module - module of mean stress) 608 nn_htau = 1 ! type of exponential decrease of tke penetration 599 ln_mxl0 = .true. ! surface mixing length scale = F(wind stress) (T) or not (F) 600 rn_mxl0 = 0.04 ! surface buoyancy lenght scale minimum value 601 ln_lc = .true. ! Langmuir cell parameterisation (Axell 2002) 602 rn_lc = 0.15 ! coef. associated to Langmuir cells 603 nn_etau = 1 ! penetration of tke below the mixed layer (ML) due to internal & intertial waves 604 ! = 0 no penetration 605 ! = 1 add a tke source below the ML 606 ! = 2 add a tke source just at the base of the ML 607 ! = 3 as = 1 applied on HF part of the stress ("key_coupled") 608 rn_efr = 0.05 ! fraction of surface tke value which penetrates below the ML (nn_etau=1 or 2) 609 nn_htau = 1 ! type of exponential decrease of tke penetration below the ML 609 610 ! = 0 constant 10 m length scale 610 ! = 1 0.5m at the equator to 30m at high latitudes 611 ! = 2 30 meters constant depth penetration 612 ! otion used only if nn_etau = 1 or 2: 613 rn_efr = 0.05 ! fraction of surface tke value which penetrates inside the ocean 614 ! otion used only if nn_etau = 3: 615 rn_addhft = -1.e-3 ! add offset applied to the "mean of stress module - module of mean stress" (always kept > 0) 616 rn_sclhft = 1. ! scale factor applied to the "mean of stress module - module of mean stress" 617 ln_lc = .false. ! Langmuir cell parameterisation 618 rn_lc = 0.15 ! coef. associated to Langmuir cells 611 ! = 1 0.5m at the equator to 30m poleward of 40 degrees 619 612 / 620 613 !------------------------------------------------------------------------ -
branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/ZDF/zdftke.F90
r2287 r2375 61 61 62 62 #if defined key_c1d 63 ! !!* 1D cfg only *64 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: e_dis, e_mix 65 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: e_pdl, e_ric 63 ! !!** 1D cfg only ** ('key_c1d') 64 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: e_dis, e_mix !: dissipation and mixing turbulent lengh scales 65 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: e_pdl, e_ric !: prandl and local Richardson numbers 66 66 #endif 67 67 68 ! !!! ** Namelist namzdf_tke ** 69 LOGICAL :: ln_mxl0 = .FALSE. ! mixing length scale surface value as function of wind stress or not 70 INTEGER :: nn_mxl = 2 ! type of mixing length (=0/1/2/3) 71 REAL(wp) :: rn_lmin0 = 0.4_wp ! surface min value of mixing length [m] 72 REAL(wp) :: rn_lmin = 0.1_wp ! interior min value of mixing length [m] 73 INTEGER :: nn_pdl = 1 ! Prandtl number or not (ratio avt/avm) (=0/1) 74 REAL(wp) :: rn_ediff = 0.1_wp ! coefficient for avt: avt=rn_ediff*mxl*sqrt(e) 75 REAL(wp) :: rn_ediss = 0.7_wp ! coefficient of the Kolmogoroff dissipation 76 REAL(wp) :: rn_ebb = 3.75_wp ! coefficient of the surface input of tke 77 REAL(wp) :: rn_emin = 0.7071e-6_wp ! minimum value of tke [m2/s2] 78 REAL(wp) :: rn_emin0 = 1.e-4_wp ! surface minimum value of tke [m2/s2] 79 REAL(wp) :: rn_bshear = 1.e-20 ! background shear (>0) 80 INTEGER :: nn_etau = 0 ! type of depth penetration of surface tke (=0/1/2) 81 INTEGER :: nn_htau = 0 ! type of tke profile of penetration (=0/1) 82 REAL(wp) :: rn_efr = 1.0_wp ! fraction of TKE surface value which penetrates in the ocean 83 REAL(wp) :: rn_addhft = 0.0_wp ! add offset applied to HF tau 84 REAL(wp) :: rn_sclhft = 1.0_wp ! scale factor applied to HF tau 85 LOGICAL :: ln_lc = .FALSE. ! Langmuir cells (LC) as a source term of TKE or not 86 REAL(wp) :: rn_lc = 0.15_wp ! coef to compute vertical velocity of Langmuir cells 87 88 REAL(wp) :: ri_cri ! critic Richardson number (deduced from rn_ediff and rn_ediss values) 89 90 REAL(wp), DIMENSION(jpi,jpj) :: htau ! depth of tke penetration (nn_htau) 91 REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) :: en ! now turbulent kinetic energy [m2/s2] 92 REAL(wp), DIMENSION(jpi,jpj,jpk) :: dissl ! now mixing lenght of dissipation 68 ! !!** Namelist namzdf_tke ** 69 LOGICAL :: ln_mxl0 = .FALSE. ! mixing length scale surface value as function of wind stress or not 70 INTEGER :: nn_mxl = 2 ! type of mixing length (=0/1/2/3) 71 REAL(wp) :: rn_mxl0 = 0.04_wp ! surface min value of mixing length (kappa*z_o=0.4*0.1 m) [m] 72 INTEGER :: nn_pdl = 1 ! Prandtl number or not (ratio avt/avm) (=0/1) 73 REAL(wp) :: rn_ediff = 0.1_wp ! coefficient for avt: avt=rn_ediff*mxl*sqrt(e) 74 REAL(wp) :: rn_ediss = 0.7_wp ! coefficient of the Kolmogoroff dissipation 75 REAL(wp) :: rn_ebb = 3.75_wp ! coefficient of the surface input of tke 76 REAL(wp) :: rn_emin = 0.7071e-6_wp ! minimum value of tke [m2/s2] 77 REAL(wp) :: rn_emin0 = 1.e-4_wp ! surface minimum value of tke [m2/s2] 78 REAL(wp) :: rn_bshear = 1.e-20_wp ! background shear (>0) currently a numerical threshold (do not change it) 79 INTEGER :: nn_etau = 0 ! type of depth penetration of surface tke (=0/1/2/3) 80 INTEGER :: nn_htau = 0 ! type of tke profile of penetration (=0/1) 81 REAL(wp) :: rn_efr = 1.0_wp ! fraction of TKE surface value which penetrates in the ocean 82 LOGICAL :: ln_lc = .FALSE. ! Langmuir cells (LC) as a source term of TKE or not 83 REAL(wp) :: rn_lc = 0.15_wp ! coef to compute vertical velocity of Langmuir cells 84 85 REAL(wp) :: ri_cri ! critic Richardson number (deduced from rn_ediff and rn_ediss values) 86 REAL(wp) :: rmxl_min ! minimum mixing length value (deduced from rn_ediff and rn_emin values) [m] 87 REAL(wp) :: rhftau_add = 1.e-3_wp ! add offset applied to HF part of taum (nn_etau=3) 88 REAL(wp) :: rhftau_scl = 1.0_wp ! scale factor applied to HF part of taum (nn_etau=3) 89 90 REAL(wp), DIMENSION(jpi,jpj,jpk), PUBLIC :: en ! now turbulent kinetic energy [m2/s2] 91 92 REAL(wp), DIMENSION(jpi,jpj) :: htau ! depth of tke penetration (nn_htau) 93 REAL(wp), DIMENSION(jpi,jpj,jpk) :: dissl ! now mixing lenght of dissipation 93 94 94 95 !! * Substitutions … … 98 99 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 99 100 !! $Id$ 100 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)101 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 101 102 !!---------------------------------------------------------------------- 102 103 103 CONTAINS 104 104 … … 165 165 !! 166 166 !! ** Method : - TKE surface boundary condition 167 !! - source term due to Langmuir cells ( ln_lc=T)167 !! - source term due to Langmuir cells (Axell JGR 2002) (ln_lc=T) 168 168 !! - source term due to shear (saved in avmu, avmv arrays) 169 169 !! - Now TKE : resolution of the TKE equation by inverting … … 197 197 ! 198 198 zbbrau = rn_ebb / rau0 ! Local constant initialisation 199 zfact1 = -.5 * rdt200 zfact2 = 1.5 * rdt * rn_ediss201 zfact3 = 0.5 * rn_ediss199 zfact1 = -.5_wp * rdt 200 zfact2 = 1.5_wp * rdt * rn_ediss 201 zfact3 = 0.5_wp * rn_ediss 202 202 ! 203 203 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< … … 245 245 ! 246 246 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 247 IF( ln_lc ) THEN ! Langmuir circulation source term added to tke 247 IF( ln_lc ) THEN ! Langmuir circulation source term added to tke (Axell JGR 2002) 248 248 ! !>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> 249 249 ! 250 250 ! !* total energy produce by LC : cumulative sum over jk 251 zpelc(:,:,1) = MAX( rn2b(:,:,1), 0. ) * fsdepw(:,:,1) * fse3w(:,:,1)251 zpelc(:,:,1) = MAX( rn2b(:,:,1), 0._wp ) * fsdepw(:,:,1) * fse3w(:,:,1) 252 252 DO jk = 2, jpk 253 zpelc(:,:,jk) = zpelc(:,:,jk-1) + MAX( rn2b(:,:,jk), 0. ) * fsdepw(:,:,jk) * fse3w(:,:,jk)253 zpelc(:,:,jk) = zpelc(:,:,jk-1) + MAX( rn2b(:,:,jk), 0._wp ) * fsdepw(:,:,jk) * fse3w(:,:,jk) 254 254 END DO 255 255 ! !* finite Langmuir Circulation depth … … 328 328 & / ( fse3t(ji,jj,jk-1) * fse3w(ji,jj,jk ) ) 329 329 ! ! shear prod. at w-point weightened by mask 330 zesh2 = ( avmu(ji-1,jj,jk) + avmu(ji,jj,jk) ) / MAX( 1. e0, umask(ji-1,jj,jk) + umask(ji,jj,jk) ) &331 & + ( avmv(ji,jj-1,jk) + avmv(ji,jj,jk) ) / MAX( 1. e0, vmask(ji,jj-1,jk) + vmask(ji,jj,jk) )330 zesh2 = ( avmu(ji-1,jj,jk) + avmu(ji,jj,jk) ) / MAX( 1._wp , umask(ji-1,jj,jk) + umask(ji,jj,jk) ) & 331 & + ( avmv(ji,jj-1,jk) + avmv(ji,jj,jk) ) / MAX( 1._wp , vmask(ji,jj-1,jk) + vmask(ji,jj,jk) ) 332 332 ! 333 333 zd_up(ji,jj,jk) = zzd_up ! Matrix (zdiag, zd_up, zd_lw) 334 334 zd_lw(ji,jj,jk) = zzd_lw 335 zdiag(ji,jj,jk) = 1. e0- zzd_lw - zzd_up + zfact2 * dissl(ji,jj,jk) * tmask(ji,jj,jk)335 zdiag(ji,jj,jk) = 1._wp - zzd_lw - zzd_up + zfact2 * dissl(ji,jj,jk) * tmask(ji,jj,jk) 336 336 ! 337 337 ! ! right hand side in en … … 384 384 ! ! TKE due to surface and internal wave breaking 385 385 ! !<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<< 386 IF( nn_etau == 1 ) THEN !* penetration throughout the water column386 IF( nn_etau == 1 ) THEN !* penetration below the mixed layer (rn_efr fraction) 387 387 DO jk = 2, jpkm1 388 388 DO jj = 2, jpjm1 389 389 DO ji = fs_2, fs_jpim1 ! vector opt. 390 390 en(ji,jj,jk) = en(ji,jj,jk) + rn_efr * en(ji,jj,1) * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) ) & 391 & * ( 1. e0- fr_i(ji,jj) ) * tmask(ji,jj,jk)392 END DO 393 END DO 394 END DO 395 ELSEIF( nn_etau == 2 ) THEN !* act only at the base of the mixed layer (jk=nmln) 391 & * ( 1._wp - fr_i(ji,jj) ) * tmask(ji,jj,jk) 392 END DO 393 END DO 394 END DO 395 ELSEIF( nn_etau == 2 ) THEN !* act only at the base of the mixed layer (jk=nmln) (rn_efr fraction) 396 396 DO jj = 2, jpjm1 397 397 DO ji = fs_2, fs_jpim1 ! vector opt. 398 398 jk = nmln(ji,jj) 399 399 en(ji,jj,jk) = en(ji,jj,jk) + rn_efr * en(ji,jj,1) * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) ) & 400 & * ( 1. e0- fr_i(ji,jj) ) * tmask(ji,jj,jk)401 END DO 402 END DO 403 ELSEIF( nn_etau == 3 ) THEN !* penetration throughout the water column400 & * ( 1._wp - fr_i(ji,jj) ) * tmask(ji,jj,jk) 401 END DO 402 END DO 403 ELSEIF( nn_etau == 3 ) THEN !* penetration belox the mixed layer (HF variability) 404 404 !CDIR NOVERRCHK 405 405 DO jk = 2, jpkm1 … … 410 410 ztx2 = utau(ji-1,jj ) + utau(ji,jj) 411 411 zty2 = vtau(ji ,jj-1) + vtau(ji,jj) 412 ztau = 0.5 * SQRT( ztx2 * ztx2 + zty2 * zty2 )! module of the mean stress413 zdif = taum(ji,jj) - ztau ! mean of the module - moduleof the mean414 zdif = r n_sclhft * MAX( 0.e0, zdif + rn_addhft )! apply some modifications...412 ztau = 0.5_wp * SQRT( ztx2 * ztx2 + zty2 * zty2 ) ! module of the mean stress 413 zdif = taum(ji,jj) - ztau ! mean of modulus - modulus of the mean 414 zdif = rhftau_scl * MAX( 0._wp, zdif + rhftau_add ) ! apply some modifications... 415 415 en(ji,jj,jk) = en(ji,jj,jk) + zbbrau * zdif * EXP( -fsdepw(ji,jj,jk) / htau(ji,jj) ) & 416 & * ( 1. e0- fr_i(ji,jj) ) * tmask(ji,jj,jk)416 & * ( 1._wp - fr_i(ji,jj) ) * tmask(ji,jj,jk) 417 417 END DO 418 418 END DO … … 439 439 !! mxl = sqrt(2*en) / N 440 440 !! where N is the brunt-vaisala frequency, with a minimum value set 441 !! to r n_lmin (rn_lmin0) in the interior (surface) ocean.441 !! to rmxl_min (rn_mxl0) in the interior (surface) ocean. 442 442 !! The mixing and dissipative length scale are bound as follow : 443 443 !! nn_mxl=0 : mxl bounded by the distance to surface and bottom. … … 479 479 ! !* Buoyancy length scale: l=sqrt(2*e/n**2) 480 480 ! 481 IF( ln_mxl0 ) THEN ! surface mixing length = F(stress) : l=vkarmn*2.e5*taum/(rau0*g) 482 !!gm this should be useless 483 zmxlm(:,:,1) = 0.e0 484 !!gm end 485 zraug = vkarmn * 2.e5 / ( rau0 * grav ) 486 DO jj = 2, jpjm1 487 DO ji = fs_2, fs_jpim1 ! vector opt. 488 zmxlm(ji,jj,1) = MAX( rn_lmin0, zraug * taum(ji,jj) ) 489 END DO 490 END DO 491 ELSE ! surface set to the minimum value 492 zmxlm(:,:,1) = rn_lmin0 481 IF( ln_mxl0 ) THEN ! surface mixing length = F(stress) : l=vkarmn*2.e5*taum/(rau0*g) 482 zraug = vkarmn * 2.e5_wp / ( rau0 * grav ) 483 zmxlm(:,:,1) = MAX( rn_mxl0, zraug * taum(ji,jj) ) 484 ELSE ! surface set to the minimum value 485 zmxlm(:,:,1) = rn_mxl0 493 486 ENDIF 494 zmxlm(:,:,jpk) = rn_lmin ! bottomset to the interior minium value495 ! 496 !CDIR NOVERRCHK 497 DO jk = 2, jpkm1 ! interior value : l=sqrt(2*e/n**2)487 zmxlm(:,:,jpk) = rmxl_min ! last level set to the interior minium value 488 ! 489 !CDIR NOVERRCHK 490 DO jk = 2, jpkm1 ! interior value : l=sqrt(2*e/n^2) 498 491 !CDIR NOVERRCHK 499 492 DO jj = 2, jpjm1 … … 501 494 DO ji = fs_2, fs_jpim1 ! vector opt. 502 495 zrn2 = MAX( rn2(ji,jj,jk), rsmall ) 503 zmxlm(ji,jj,jk) = MAX( rn_lmin, SQRT( 2. * en(ji,jj,jk) / zrn2 ) ) 504 END DO 505 END DO 506 END DO 507 ! 508 !!gm CAUTION: to be added here the bottom boundary condition on zmxlm 496 zmxlm(ji,jj,jk) = MAX( rmxl_min, SQRT( 2._wp * en(ji,jj,jk) / zrn2 ) ) 497 END DO 498 END DO 499 END DO 500 ! 509 501 ! 510 502 ! !* Physical limits for the mixing length 511 503 ! 512 zmxld(:,:, 1 ) = zmxlm(:,:,1) ! surface set to the minimumvalue513 zmxld(:,:,jpk) = r n_lmin ! bottomset to the minimum value504 zmxld(:,:, 1 ) = zmxlm(:,:,1) ! surface set to the zmxlm value 505 zmxld(:,:,jpk) = rmxl_min ! last level set to the minimum value 514 506 ! 515 507 SELECT CASE ( nn_mxl ) … … 625 617 DO jj = 2, jpjm1 626 618 DO ji = fs_2, fs_jpim1 ! vector opt. 627 zcoef = 0.5 / ( fse3w(ji,jj,jk) * fse3w(ji,jj,jk))619 zcoef = avm(ji,jj,jk) * 2._wp * fse3w(ji,jj,jk) * fse3w(ji,jj,jk) 628 620 ! ! shear 629 621 zdku = avmu(ji-1,jj,jk) * ( un(ji-1,jj,jk-1) - un(ji-1,jj,jk) ) * ( ub(ji-1,jj,jk-1) - ub(ji-1,jj,jk) ) & … … 632 624 & + avmv(ji,jj ,jk) * ( vn(ji,jj ,jk-1) - vn(ji,jj ,jk) ) * ( vb(ji,jj ,jk-1) - vb(ji,jj ,jk) ) 633 625 ! ! local Richardson number 634 zri = MAX( rn2b(ji,jj,jk), 0. ) * avm(ji,jj,jk) / (zcoef * (zdku + zdkv + rn_bshear ))635 zpdlr = MAX( 0.1 , 0.2 / MAX( 0.2 , zri ) )626 zri = MAX( rn2b(ji,jj,jk), 0._wp ) * zcoef / (zdku + zdkv + rn_bshear ) 627 zpdlr = MAX( 0.1_wp, 0.2 / MAX( 0.2 , zri ) ) 636 628 !!gm and even better with the use of the "true" ri_crit=0.22222... (this change the results!) 637 !!gm zpdlr = MAX( 0.1 , ri_crit / MAX( ri_crit , zri ) )629 !!gm zpdlr = MAX( 0.1_wp, ri_crit / MAX( ri_crit , zri ) ) 638 630 avt(ji,jj,jk) = MAX( zpdlr * avt(ji,jj,jk), avtb_2d(ji,jj) * avtb(jk) ) * tmask(ji,jj,jk) 639 631 # if defined key_c1d … … 672 664 INTEGER :: ji, jj, jk ! dummy loop indices 673 665 !! 674 NAMELIST/namzdf_tke/ rn_ediff, rn_ediss , rn_ebb , rn_emin , & 675 & rn_emin0, rn_bshear, nn_mxl , ln_mxl0 , & 676 & rn_lmin , rn_lmin0 , nn_pdl , nn_etau , & 677 & nn_htau , rn_efr , rn_addhft, rn_sclhft, & 678 & ln_lc , rn_lc 666 NAMELIST/namzdf_tke/ rn_ediff, rn_ediss , rn_ebb , rn_emin , & 667 & rn_emin0, rn_bshear, nn_mxl , ln_mxl0 , & 668 & rn_mxl0 , nn_pdl , ln_lc , rn_lc , & 669 & nn_etau , nn_htau , rn_efr 679 670 !!---------------------------------------------------------------------- 680 671 … … 682 673 READ ( numnam, namzdf_tke ) 683 674 684 ri_cri = 2. / ( 2. + rn_ediss / rn_ediff ) ! resulting critical Richardson number 675 ri_cri = 2._wp / ( 2._wp + rn_ediss / rn_ediff ) ! resulting critical Richardson number 676 rmxl_min = 1.e-6_wp / ( rn_ediff * SQRT( rn_emin ) ) ! resulting minimum length to recover molecular viscosity 685 677 686 678 IF(lwp) THEN !* Control print … … 698 690 WRITE(numout,*) ' prandl number flag nn_pdl = ', nn_pdl 699 691 WRITE(numout,*) ' surface mixing length = F(stress) or not ln_mxl0 = ', ln_mxl0 700 WRITE(numout,*) ' surface mixing length minimum value rn_lmin0 = ', rn_lmin0 701 WRITE(numout,*) ' interior mixing length minimum value rn_lmin0 = ', rn_lmin0 692 WRITE(numout,*) ' surface mixing length minimum value rn_mxl0 = ', rn_mxl0 693 WRITE(numout,*) ' flag to take into acc. Langmuir circ. ln_lc = ', ln_lc 694 WRITE(numout,*) ' coef to compute verticla velocity of LC rn_lc = ', rn_lc 702 695 WRITE(numout,*) ' test param. to add tke induced by wind nn_etau = ', nn_etau 703 696 WRITE(numout,*) ' flag for computation of exp. tke profile nn_htau = ', nn_htau 704 697 WRITE(numout,*) ' fraction of en which pene. the thermocline rn_efr = ', rn_efr 705 WRITE(numout,*) ' add offset applied to HF tau rn_addhft = ', rn_addhft706 WRITE(numout,*) ' scale factor applied to HF tau rn_sclhft = ', rn_sclhft707 WRITE(numout,*) ' flag to take into acc. Langmuir circ. ln_lc = ', ln_lc708 WRITE(numout,*) ' coef to compute verticla velocity of LC rn_lc = ', rn_lc709 698 WRITE(numout,*) 710 699 WRITE(numout,*) ' critical Richardson nb with your parameters ri_cri = ', ri_cri … … 712 701 713 702 ! !* Check of some namelist values 714 IF( nn_mxl < 0 .OR. nn_mxl > 3 ) CALL ctl_stop( 'bad flag: nn_mxl is 0, 1 or 2 ' ) 715 IF( nn_pdl < 0 .OR. nn_pdl > 1 ) CALL ctl_stop( 'bad flag: nn_pdl is 0 or 1 ' ) 716 IF( nn_htau < 0 .OR. nn_htau > 2 ) CALL ctl_stop( 'bad flag: nn_htau is 0, 1 or 2 ' ) 717 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 ' ) 718 IF( rn_sclhft == 0. .AND. nn_etau == 3 ) CALL ctl_stop( 'force null HF tau to penetrate the thermocline...' ) 719 IF( .NOT. lhftau .AND. nn_etau == 3 ) CALL ctl_stop( 'bad flag: nn_etau == 3 must be used with HF tau' ) 720 703 IF( nn_mxl < 0 .OR. nn_mxl > 3 ) CALL ctl_stop( 'bad flag: nn_mxl is 0, 1 or 2 ' ) 704 IF( nn_pdl < 0 .OR. nn_pdl > 1 ) CALL ctl_stop( 'bad flag: nn_pdl is 0 or 1 ' ) 705 IF( nn_htau < 0 .OR. nn_htau > 1 ) CALL ctl_stop( 'bad flag: nn_htau is 0, 1 or 2 ' ) 706 #if ! key_coupled 707 IF( nn_etau == 3 ) CALL ctl_stop( 'nn_etau == 3 : HF taum only known in coupled mode' ) 708 #endif 709 710 IF( ln_mxl0 ) THEN 711 IF(lwp) WRITE(numout,*) ' use a surface mixing length = F(stress) : set rn_mxl0 = rmxl_min' 712 rn_mxl0 = rmxl_min 713 ENDIF 714 721 715 IF( nn_etau == 2 ) CALL zdf_mxl( nit000 ) ! Initialization of nmln 722 716 … … 724 718 IF( nn_etau /= 0 ) THEN 725 719 SELECT CASE( nn_htau ) ! Choice of the depth of penetration 726 CASE( 0 ) ! constant depth penetration (here 10 meters) 727 htau(:,:) = 10.e0 728 CASE( 1 ) ! F(latitude) : 0.5m to 30m at high lat. 729 DO jj = 1, jpj 730 DO ji = 1, jpi 731 htau(ji,jj) = MAX( 0.5, 3./4. * MIN( 40., 60.*ABS( SIN( rpi/180. * gphit(ji,jj) ) ) ) ) 732 END DO 733 END DO 734 CASE( 2 ) ! constant depth penetration (here 30 meters) 735 htau(:,:) = 30.e0 720 CASE( 0 ) ! constant depth penetration (here 10 meters) 721 htau(:,:) = 10._wp 722 CASE( 1 ) ! F(latitude) : 0.5m to 30m poleward of 40 degrees 723 htau(:,:) = MAX( 0.5_wp, MIN( 30._wp, 45._wp* ABS( SIN( rpi/180._wp * gphit(:,:) ) ) ) ) 736 724 END SELECT 737 725 ENDIF … … 744 732 avmv(:,:,jk) = avmb(jk) * vmask(:,:,jk) 745 733 END DO 746 dissl(:,:,:) = 1.e-12 734 dissl(:,:,:) = 1.e-12_wp 747 735 ! !* read or initialize all required files 748 736 CALL tke_rst( nit000, 'READ' )
Note: See TracChangeset
for help on using the changeset viewer.