Changeset 879 for LMDZ4/trunk/libf/phylmd/physiq.F
- Timestamp:
- Jan 17, 2008, 10:20:32 PM (16 years ago)
- File:
-
- 1 edited
-
LMDZ4/trunk/libf/phylmd/physiq.F (modified) (18 diffs)
Legend:
- Unmodified
- Added
- Removed
-
LMDZ4/trunk/libf/phylmd/physiq.F
r878 r879 108 108 LOGICAL ok_gust ! pour activer l'effet des gust sur flux surface 109 109 PARAMETER (ok_gust=.FALSE.) 110 integer iflag_radia ! active ou non le rayonnement (MPL) 111 save iflag_radia 110 112 c====================================================================== 111 113 LOGICAL check ! Verifier la conservation du modele en eau … … 131 133 REAL fluxg(klon) !flux turbulents ocean-atmosphere 132 134 REAL amn, amx 135 INTEGER igout 133 136 c====================================================================== 134 137 c Clef controlant l'activation du cycle diurne: … … 222 225 REAL u(klon,klev) 223 226 REAL v(klon,klev) 224 REAL t(klon,klev) 227 REAL t(klon,klev),theta(klon,klev) 225 228 REAL qx(klon,klev,nqmax) 226 229 … … 784 787 cym SAVE wd ! sb 785 788 789 REAL,allocatable,save :: sigd(:) 790 c 791 c================================================================================================= 792 cCR04.12.07: on ajoute les nouvelles variables du nouveau schema de convection avec poches froides 793 c Variables liées à la poche froide (jyg) 794 795 REAL mip(klon,klev) ! mass flux shed by the adiab ascent at each level 796 REAL Vprecip(klon,klev) ! precipitation vertical profile 797 REAL,allocatable,save :: cin(:) ! CIN 798 REAL,allocatable,save :: ftd(:,:) ! differential heating between wake and environment 799 REAL,allocatable,save :: fqd(:,:) ! differential moistening between wake and environment 800 c34EK 801 c -- Variables de controle de ALE et ALP 802 REAL,allocatable,save :: ALE(:) ! Energie disponible pour soulevement : utilisee par la 803 c convection d'Emanuel pour le declenchement et la regulation 804 REAL,allocatable,save :: ALP(:) ! Puissance disponible pour soulevement ! 805 c 806 REAL wape_prescr, fip_prescr 807 INTEGER it_wape_prescr 808 SAVE wape_prescr, fip_prescr, it_wape_prescr 809 c 810 c variables supplementaires de concvl 811 REAL Tconv(klon,klev) 812 REAL ment(klon,klev,klev),sij(klon,klev,klev) 813 REAL dd_t(klon,klev),dd_q(klon,klev) 814 cnouvelles variables pour le couplage convection-couche limite 815 real,allocatable,save :: Ale_bl(:) 816 real,allocatable,save :: Alp_bl(:) 817 integer,allocatable,save :: lalim_conv(:) 818 real,allocatable,save :: wght_th(:,:) 819 real alp_bl_prescr 820 save alp_bl_prescr 821 real ale_bl_prescr 822 save ale_bl_prescr 823 real ale_wake(klon) 824 real alp_wake(klon) 825 cRC 826 c Variables liées à la poche froide (jyg et rr) 827 c Version diagnostique pour l'instant : pas de rétroaction sur la convection 828 829 REAL t_wake(klon,klev),q_wake(klon,klev) ! wake pour la convection 830 831 REAL,allocatable,save :: wake_deltat(:,:) ! Wake : ecart de temperature avec la 832 c zone non perturbee 833 REAL,allocatable,save :: wake_deltaq(:,:) ! Wake : ecart d'humidite avec la 834 c zone non perturbee 835 REAL wake_dth(klon,klev) ! wake : temp pot difference 836 837 REAL wake_d_deltat_gw(klon,klev)! wake : delta T tendency due to Gravity Wave (/s) 838 REAL wake_omgbdth(klon,klev) ! Wake : flux of Delta_Theta transported by LS omega 839 REAL wake_dp_omgb(klon,klev) ! Wake : vertical gradient of large scale omega 840 REAL wake_dtKE(klon,klev) ! Wake : differential heating (wake - unpertubed) CONV 841 REAL wake_dqKE(klon,klev) ! Wake : differential moistening (wake - unpertubed) CONV 842 REAL wake_dtPBL(klon,klev) ! Wake : differential heating (wake - unpertubed) PBL 843 REAL wake_dqPBL(klon,klev) ! Wake : differential moistening (wake - unpertubed) PBL 844 REAL wake_omg(klon,klev) ! Wake : velocity difference (wake - unpertubed) 845 REAL wake_ddeltat(klon,klev),wake_ddeltaq(klon,klev) 846 REAL wake_dp_deltomg(klon,klev) ! Wake : gradient vertical de wake_omg 847 REAL wake_spread(klon,klev) ! spreading term in wake_delt 848 REAL,allocatable,save :: wake_Cstar(:) ! vitesse d'etalement de la poche 849 cpourquoi y'a pas de save?? 850 REAL wake_h(klon) ! Wake : hauteur de la poche froide 851 REAL,allocatable,save :: wake_s(:) ! Wake : fraction surfacique occupee par 852 c la poche froide 853 INTEGER wake_k(klon) ! Wake sommet 854 c 855 REAL t_undi(klon,klev) ! temperature moyenne dans la zone non perturbee 856 REAL q_undi(klon,klev) ! humidite moyenne dans la zone non perturbee 857 c 858 REAL wake_pe(klon) ! Wake potential energy - WAPE 859 REAL,allocatable,save :: wake_fip(:) ! Gust Front Impinging power - ALP 860 861 REAL wake_gfl(klon) ! Gust Front Length 862 REAL wake_dens(klon) 863 c 864 REAL,allocatable,save :: dt_wake(:,:) 865 REAL,allocatable,save :: dq_wake(:,:) ! LS tendencies due to wake 866 c 867 REAL dt_dwn(klon,klev) 868 REAL dq_dwn(klon,klev) 869 REAL wdt_PBL(klon,klev) 870 REAL udt_PBL(klon,klev) 871 REAL wdq_PBL(klon,klev) 872 REAL udq_PBL(klon,klev) 873 REAL M_dwn(klon,klev) 874 REAL M_up(klon,klev) 875 REAL dt_a(klon,klev) 876 REAL dq_a(klon,klev) 877 c 878 cRR:fin declarations poches froides 879 c======================================================================================================= 880 786 881 c Variables locales pour la couche limite (al1): 787 882 c … … 909 1004 EXTERNAL phyredem ! ecrire l'etat de redemarrage de la physique 910 1005 EXTERNAL radlwsw ! rayonnements solaire et infrarouge 911 EXTERNAL suphe c! initialiser certaines constantes1006 EXTERNAL suphel ! initialiser certaines constantes 912 1007 EXTERNAL transp ! transport total de l'eau et de l'energie 913 1008 EXTERNAL ecribina ! ecrire le fichier binaire global … … 1378 1473 c 1379 1474 c====================================================================== 1475 ! Ecriture eventuelle d'un profil verticale en entree de la physique. 1476 ! Utilise notamment en 1D mais peut etre active egalement en 3D 1477 ! en imposant la valeur de igout. 1478 c====================================================================== 1479 1480 if (klon.eq.1) then 1481 print*,'WARNING !!!! omega=0' 1482 omega=0. 1483 igout=1 1484 write(lunout,*) 'DEBUT DE PHYSIQ !!!!!!!!!!!!!!!!!!!!' 1485 write(lunout,*) 1486 s 'nlon,nlev,nqmax,debut,lafin,rjourvrai,gmtime,pdtphys' 1487 write(lunout,*) 1488 s nlon,nlev,nqmax,debut,lafin,rjourvrai,gmtime,pdtphys 1489 1490 write(lunout,*) 'papers, play, phi, u, v, t, omega' 1491 do k=1,nlev 1492 write(lunout,*) paprs(igout,k),pplay(igout,k),pphi(igout,k), 1493 s u(igout,k),v(igout,k),t(igout,k),omega(igout,k) 1494 enddo 1495 write(lunout,*) 'ovap (g/kg), oliq (g/kg)' 1496 do k=1,nlev 1497 write(lunout,*) qx(igout,k,1)*1000,qx(igout,k,2)*1000. 1498 enddo 1499 endif 1500 1501 c====================================================================== 1380 1502 1381 1503 cym => necessaire pour iflag_con != 2 … … 1417 1539 allocate( zuthe(klon),zvthe(klon)) 1418 1540 allocate( alb_neig(klon)) 1419 allocate( ema_workcbmf(klon)) 1541 allocate( ema_workcbmf(klon)) 1542 cCR:nvelles variables convection/poches froides 1543 allocate( sigd(klon)) 1544 allocate( cin(klon)) 1545 allocate( ftd(klon,klev)) 1546 allocate( fqd(klon,klev)) 1547 allocate( ALE(klon)) 1548 allocate( ALP(klon)) 1549 allocate( Ale_bl(klon)) 1550 allocate( Alp_bl(klon)) 1551 allocate( lalim_conv(klon)) 1552 allocate( wght_th(klon,klev)) 1553 allocate( wake_deltat(klon,klev)) 1554 allocate( wake_deltaq(klon,klev)) 1555 allocate( wake_Cstar(klon)) 1556 allocate( wake_s(klon)) 1557 allocate( wake_fip(klon)) 1558 allocate( dt_wake(klon,klev)) 1559 allocate( dq_wake(klon,klev)) 1560 cfinCR 1420 1561 allocate( ema_cbmf(klon)) 1421 1562 allocate( ema_pcb(klon)) … … 1528 1669 ENDIF 1529 1670 IF (debut) THEN 1530 CALL suphe c! initialiser constantes et parametres phys.1671 CALL suphel ! initialiser constantes et parametres phys. 1531 1672 ENDIF 1532 1673 … … 1543 1684 C 1544 1685 !rv 1686 cCRinitialisation de wght_th et lalim_conv pour la definition de la couche alimentation 1687 cde la convection a partir des caracteristiques du thermique 1688 wght_th(:,:)=1. 1689 lalim_conv(:)=1 1690 cRC 1545 1691 u10m(:,:)=0. 1546 1692 v10m(:,:)=0. … … 1593 1739 call conf_phys(ocean, ok_veget, ok_journe, ok_mensuel, 1594 1740 . ok_instan, ok_hf, seuil_inversion, 1595 . fact_cldcon, facttemps,ok_newmicro, 1741 . fact_cldcon, facttemps,ok_newmicro,iflag_radia, 1596 1742 . iflag_cldcon,iflag_ratqs,ratqsbas,ratqshaut, 1597 1743 . ok_ade, ok_aie, 1598 1744 . bl95_b0, bl95_b1, 1599 . iflag_thermals,nsplit_thermals) 1745 . iflag_thermals,nsplit_thermals, 1746 cnv flags pour la convection et les poches froides 1747 . iflag_coupl,iflag_clos,iflag_wake) 1748 1749 print*,'iflag_coupl,iflag_clos,iflag_wake', 1750 . iflag_coupl,iflag_clos,iflag_wake 1600 1751 1601 1752 c … … 1711 1862 ENDDO 1712 1863 cIM15/11/02 rajout initialisation ibas_con,itop_con cf. SB =>END 1864 c=============================================================================== 1865 cCR:04.12.07: initialisations poches froides 1866 c Controle de ALE et ALP pour la fermeture convective (jyg) 1867 if (iflag_wake.eq.1) then 1868 CALL ini_wake(0.,0.,it_wape_prescr,wape_prescr,fip_prescr 1869 s ,alp_bl_prescr, ale_bl_prescr) 1870 c 11/09/06 rajout initialisation ALE et ALP du wake et PBL(YU) 1871 c print*,'apres ini_wake iflag_cldcon=', iflag_cldcon 1872 endif 1873 1874 do i = 1,klon 1875 wake_s(i) = 0. 1876 wake_fip(i) = 0. 1877 wake_cstar(i) = 0. 1878 DO k=1,klev 1879 wake_deltat(i,k)=0. 1880 wake_deltaq(i,k)=0. 1881 ENDDO 1882 enddo 1883 c================================================================================ 1713 1884 1714 1885 ENDIF … … 2230 2401 c supprimer les calculs / ftra. 2231 2402 ntra = 1 2403 2404 c===================================================================================== 2405 cajout pour la parametrisation des poches froides: 2406 ccalcul de t_wake et t_undi: si pas de poches froides, t_wake=t_undi=t_seri 2407 do k=1,klev 2408 do i=1,klon 2409 if (iflag_wake.eq.1) then 2410 t_wake(i,k) = t_seri(i,k) 2411 . +(1-wake_s(i))*wake_deltat(i,k) 2412 q_wake(i,k) = q_seri(i,k) 2413 . +(1-wake_s(i))*wake_deltaq(i,k) 2414 t_undi(i,k) = t_seri(i,k) 2415 . -wake_s(i)*wake_deltat(i,k) 2416 q_undi(i,k) = q_seri(i,k) 2417 . -wake_s(i)*wake_deltaq(i,k) 2418 else 2419 t_wake(i,k) = t_seri(i,k) 2420 q_wake(i,k) = q_seri(i,k) 2421 t_undi(i,k) = t_seri(i,k) 2422 q_undi(i,k) = q_seri(i,k) 2423 endif 2424 enddo 2425 enddo 2426 2427 cc-- Calcul de l'energie disponible ALE (J/kg) et de la puissance disponible ALP (W/m2) 2428 cc-- pour le soulevement des particules dans le modele convectif 2429 c 2430 do i = 1,klon 2431 ALE(i) = 0. 2432 ALP(i) = 0. 2433 enddo 2434 c 2435 ccalcul de ale_wake et alp_wake 2436 do i = 1,klon 2437 if (iflag_wake.eq.1) then 2438 ale_wake(i) = 0.5*wake_cstar(i)**2 2439 alp_wake(i) = wake_fip(i) 2440 else 2441 ale_wake(i) = 0. 2442 alp_wake(i) = 0. 2443 endif 2444 enddo 2445 ccombinaison avec ale et alp de couche limite: constantes si pas de couplage, valeurs calculees 2446 cdans le thermique sinon 2447 if (iflag_coupl.eq.0) then 2448 print*,'ALE et ALP imposes' 2449 do i = 1,klon 2450 con ne couple que ale 2451 c ALE(i) = max(ale_wake(i),Ale_bl(i)) 2452 ALE(i) = max(ale_wake(i),ale_bl_prescr) 2453 con ne couple que alp 2454 c ALP(i) = alp_wake(i) + Alp_bl(i) 2455 ALP(i) = alp_wake(i) + alp_bl_prescr 2456 enddo 2457 else 2458 print*,'ALE et ALP couples au thermique' 2459 do i = 1,klon 2460 ALE(i) = max(ale_wake(i),Ale_bl(i)) 2461 ALP(i) = alp_wake(i) + Alp_bl(i) 2462 c write(20,*)'ALE',ALE(i),Ale_bl(i),ale_wake(i) 2463 c write(21,*)'ALP',ALP(i),Alp_bl(i),alp_wake(i) 2464 enddo 2465 endif 2466 2467 cfin calcul ale et alp 2468 c================================================================================================= 2469 2470 2232 2471 c sb, oct02: 2233 2472 c Schema de convection modularise et vectorise: … … 2236 2475 IF (ok_cvl) THEN ! new driver for convectL 2237 2476 2238 CALL concvl (iflag_con, 2239 . dtime,paprs,pplay,t_seri,q_seri, 2240 . u_seri,v_seri,tr_seri,ntra, 2477 CALL concvl (iflag_con,iflag_clos, 2478 . dtime,paprs,pplay,t_undi,q_undi, 2479 . t_wake,q_wake, 2480 . u_seri,v_seri,tr_seri,nbtr, 2481 . ALE,ALP, 2241 2482 . ema_work1,ema_work2, 2242 2483 . d_t_con,d_q_con,d_u_con,d_v_con,d_tr, 2243 . rain_con, snow_con, ibas_con, itop_con, 2484 . rain_con, snow_con, ibas_con, itop_con, sigd, 2244 2485 . upwd,dnwd,dnwd0, 2245 . Ma, cape,tvp,iflagctrl,2486 . Ma,mip,Vprecip,cape,cin,tvp,Tconv,iflagctrl, 2246 2487 . pbase,bbase,dtvpdt1,dtvpdq1,dplcldt,dplcldr,qcondc,wd, 2247 . pmflxr,pmflxs, 2248 . da,phi,mp)2488 . pmflxr,pmflxs,da,phi,mp, 2489 . ftd,fqd,lalim_conv,wght_th) 2249 2490 2250 2491 cIM cf. FH … … 2314 2555 ENDDO 2315 2556 DO i = 1, klon 2316 ema_pct(i) = paprs(i,itop_con(i)) 2557 2558 ! L'idicage de itop_con peut cacher un pb potentiel 2559 ! FH sous la dictee de JYG, CR 2560 ema_pct(i) = paprs(i,itop_con(i)+1) 2561 2317 2562 if (itop_con(i).gt.klev-3) then 2318 2563 print*,'La convection monte trop haut ' … … 2397 2642 ENDIF 2398 2643 zx_ajustq=.FALSE. 2644 2645 c 2646 c============================================================================= 2647 cRR:Evolution de la poche froide: on ne fait pas de separation wake/env 2648 cpour la couche limite diffuse pour l instant 2649 c 2650 if (iflag_wake.eq.1) then 2651 DO k=1,klev 2652 DO i=1,klon 2653 dt_dwn(i,k) = ftd(i,k) 2654 wdt_PBL(i,k) = 0. 2655 dq_dwn(i,k) = fqd(i,k) 2656 wdq_PBL(i,k) = 0. 2657 M_dwn(i,k) = dnwd0(i,k) 2658 M_up(i,k) = upwd(i,k) 2659 dt_a(i,k) = d_t_con(i,k)/dtime - ftd(i,k) 2660 udt_PBL(i,k) = 0. 2661 dq_a(i,k) = d_q_con(i,k)/dtime - fqd(i,k) 2662 udq_PBL(i,k) = 0. 2663 ENDDO 2664 ENDDO 2665 c 2666 ccalcul caracteristiques de la poche froide 2667 call calWAKE (paprs,pplay,dtime 2668 : ,t_seri,q_seri,omega,ibas_con 2669 : ,dt_dwn,dq_dwn,M_dwn,M_up 2670 : ,dt_a,dq_a,sigd 2671 : ,wdt_PBL,wdq_PBL 2672 : ,udt_PBL,udq_PBL 2673 o ,wake_deltat,wake_deltaq,wake_dth 2674 o ,wake_h,wake_s,wake_dens 2675 o ,wake_pe,wake_fip,wake_gfl 2676 o ,dt_wake,dq_wake 2677 o ,wake_k, t_undi,q_undi 2678 o ,wake_omgbdth,wake_dp_omgb 2679 o ,wake_dtKE,wake_dqKE 2680 o ,wake_dtPBL,wake_dqPBL 2681 o ,wake_omg,wake_dp_deltomg 2682 o ,wake_spread,wake_Cstar,wake_d_deltat_gw 2683 o ,wake_ddeltat,wake_ddeltaq) 2684 c 2685 cIncrementation des tendances de la poche froide 2686 DO k = 1, klev 2687 DO i = 1, klon 2688 t_seri(i,k) = t_seri(i,k) + dt_wake(i,k)*dtime 2689 q_seri(i,k) = q_seri(i,k) + dq_wake(i,k)*dtime 2690 ENDDO 2691 ENDDO 2692 2693 endif 2694 c print*,'apres callwake iflag_cldcon=', iflag_cldcon 2399 2695 c 2400 2696 c=================================================================== … … 2441 2737 s ,u_seri,v_seri,t_seri,q_seri,zqsat,debut 2442 2738 s ,d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs 2443 s ,fm_therm,entr_therm,zqasc,clwcon0th,lmax_th,ratqscth, 2444 s ratqsdiff,zqsatth) 2739 s ,fm_therm,entr_therm,zqasc,clwcon0th,lmax_th,ratqscth 2740 s ,ratqsdiff,zqsatth 2741 con rajoute ale et alp, et les caracteristiques de la couche alim 2742 s ,Ale_bl,Alp_bl,lalim_conv,wght_th) 2445 2743 endif 2446 2744 … … 3016 3314 ENDIF 3017 3315 itaprad = itaprad + 1 3316 3317 if (iflag_radia.eq.0) then 3318 print *,'--------------------------------------------------' 3319 print *,'>>>> ATTENTION rayonnement desactive pour ce cas' 3320 print *,'>>>> heat et cool mis a zero ' 3321 print *,'--------------------------------------------------' 3322 heat=0. 3323 cool=0. 3324 endif 3325 3326 3018 3327 c 3019 3328 c Ajouter la tendance des rayonnements (tous les pas) … … 3420 3729 ENDDO 3421 3730 c 3731 !========================================================================== 3732 ! Sorties des tendances pour un point particulier 3733 ! a utiliser en 1D, avec igout=1 ou en 3D sur un point particulier 3734 ! pour le debug 3735 ! La valeur de igout est attribuee plus haut dans le programme 3736 !========================================================================== 3737 3738 if (klon.eq.1) then 3739 write(lunout,*) 'FIN DE PHYSIQ !!!!!!!!!!!!!!!!!!!!' 3740 write(lunout,*) 3741 s 'nlon,nlev,nqmax,debut,lafin,rjourvrai,gmtime,pdtphys' 3742 write(lunout,*) 3743 s nlon,nlev,nqmax,debut,lafin,rjourvrai,gmtime,pdtphys 3744 write(lunout,*) 'd_t_dyn,d_t_con,d_t_lsc,d_t_ajsb,d_t_ajs,d_t_eva' 3745 do k=1,nlev 3746 write(lunout,*) d_t_dyn(igout,k),d_t_con(igout,k), 3747 s d_t_lsc(igout,k),d_t_ajsb(igout,k),d_t_ajs(igout,k), 3748 s d_t_eva(igout,k) 3749 enddo 3750 write(lunout,*) 'cool,heat' 3751 do k=1,nlev 3752 write(lunout,*) cool(igout,k),heat(igout,k) 3753 enddo 3754 3755 write(lunout,*) 'd_t_oli,d_t_vdf,d_t_oro,d_t_lif,d_t_ec' 3756 do k=1,nlev 3757 write(lunout,*) d_t_oli(igout,k),d_t_vdf(igout,k), 3758 s d_t_oro(igout,k),d_t_lif(igout,k),d_t_ec(igout,k) 3759 enddo 3760 3761 write(lunout,*) 'd_ps ',d_ps(igout) 3762 write(lunout,*) 'd_u, d_v, d_t, d_qx1, d_qx2 ' 3763 do k=1,nlev 3764 write(lunout,*) d_u(igout,k),d_v(igout,k),d_t(igout,k), 3765 s d_qx(igout,k,1),d_qx(igout,k,2) 3766 enddo 3767 endif 3768 3769 !========================================================================== 3770 3771 c============================================================ 3772 c Calcul de la temperature potentielle 3773 c============================================================ 3774 DO k = 1, klev 3775 DO i = 1, klon 3776 theta(i,k)=t(i,k)*(100000./pplay(i,k))**(RD/RCPD) 3777 ENDDO 3778 ENDDO 3779 c 3780 3422 3781 c 22.03.04 BEG 3423 3782 c=============================================================
Note: See TracChangeset
for help on using the changeset viewer.
