Ignore:
Timestamp:
Apr 16, 2014, 7:16:58 AM (10 years ago)
Author:
fhourdin
Message:

Passage au format libre pour inclure les ficheirs du 1D dans lmdz1d.F90

File:
1 edited

Legend:

Unmodified
Added
Removed

  • LMDZ5/trunk/libf/phylmd/1D_interp_cases.h

    r2017 r2019  
    44      if (forcing_GCSSold) then
    55
    6        call get_uvd(it,timestep,fich_gcssold_ctl,fich_gcssold_dat,
    7      :               ht_gcssold,hq_gcssold,hw_gcssold,
    8      :               hu_gcssold,hv_gcssold,
    9      :               hthturb_gcssold,hqturb_gcssold,Ts_gcssold,
    10      :               imp_fcg_gcssold,ts_fcg_gcssold,
    11      :               Tp_fcg_gcssold,Turb_fcg_gcssold)
     6       call get_uvd(it,timestep,fich_gcssold_ctl,fich_gcssold_dat,              &
     7     &               ht_gcssold,hq_gcssold,hw_gcssold,                          &
     8     &               hu_gcssold,hv_gcssold,                                     &
     9     &               hthturb_gcssold,hqturb_gcssold,Ts_gcssold,                 &
     10     &               imp_fcg_gcssold,ts_fcg_gcssold,                            &
     11     &               Tp_fcg_gcssold,Turb_fcg_gcssold)
    1212       if (prt_level.ge.1) then
    1313         print *,' get_uvd -> hqturb_gcssold ',it,hqturb_gcssold
    … …  
    3838
    3939       if (prt_level.ge.1) then
    40         print*,
    41      : '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/dt_toga=',
    42      :    day,day1,(day-day1)*86400.,(day-day1)*86400/dt_toga
     40        print*,                                                             &
     41     & '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/dt_toga=',     &
     42     &    day,day1,(day-day1)*86400.,(day-day1)*86400/dt_toga
    4343       endif
    4444
    4545! time interpolation:
    46         CALL interp_toga_time(daytime,day1,annee_ref
    47      i             ,year_ini_toga,day_ju_ini_toga,nt_toga,dt_toga
    48      i             ,nlev_toga,ts_toga,plev_toga,t_toga,q_toga,u_toga
    49      i             ,v_toga,w_toga,ht_toga,vt_toga,hq_toga,vq_toga
    50      o             ,ts_prof,plev_prof,t_prof,q_prof,u_prof,v_prof,w_prof
    51      o             ,ht_prof,vt_prof,hq_prof,vq_prof)
     46        CALL interp_toga_time(daytime,day1,annee_ref                        &
     47     &             ,year_ini_toga,day_ju_ini_toga,nt_toga,dt_toga           &
     48     &             ,nlev_toga,ts_toga,plev_toga,t_toga,q_toga,u_toga        &
     49     &             ,v_toga,w_toga,ht_toga,vt_toga,hq_toga,vq_toga           &
     50     &             ,ts_prof,plev_prof,t_prof,q_prof,u_prof,v_prof,w_prof    &
     51     &             ,ht_prof,vt_prof,hq_prof,vq_prof)
    5252
    5353        if (type_ts_forcing.eq.1) ts_cur = ts_prof ! SST used in read_tsurf1d
    5454
    5555! vertical interpolation:
    56       CALL interp_toga_vertical(play,nlev_toga,plev_prof
    57      :         ,t_prof,q_prof,u_prof,v_prof,w_prof
    58      :         ,ht_prof,vt_prof,hq_prof,vq_prof
    59      :         ,t_mod,q_mod,u_mod,v_mod,w_mod
    60      :         ,ht_mod,vt_mod,hq_mod,vq_mod,mxcalc)
     56      CALL interp_toga_vertical(play,nlev_toga,plev_prof                    &
     57     &         ,t_prof,q_prof,u_prof,v_prof,w_prof                          &
     58     &         ,ht_prof,vt_prof,hq_prof,vq_prof                             &
     59     &         ,t_mod,q_mod,u_mod,v_mod,w_mod                               &
     60     &         ,ht_mod,vt_mod,hq_mod,vq_mod,mxcalc)
    6161
    6262! large-scale forcing :
    … …  
    8484      if (forcing_twpice) then
    8585
    86         print*,
    87      : '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/dt_twpi=',
    88      :    daytime,day1,(daytime-day1)*86400.,
    89      :    (daytime-day1)*86400/dt_twpi
     86        print*,                                                             &
     87     & '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/dt_twpi=',     &
     88     &    daytime,day1,(daytime-day1)*86400.,                               &
     89     &    (daytime-day1)*86400/dt_twpi
    9090
    9191! time interpolation:
    92         CALL interp_toga_time(daytime,day1,annee_ref
    93      i       ,year_ini_twpi,day_ju_ini_twpi,nt_twpi,dt_twpi,nlev_twpi
    94      i       ,ts_twpi,plev_twpi,t_twpi,q_twpi,u_twpi,v_twpi,w_twpi
    95      i       ,ht_twpi,vt_twpi,hq_twpi,vq_twpi
    96      o       ,ts_proftwp,plev_proftwp,t_proftwp,q_proftwp,u_proftwp
    97      o       ,v_proftwp,w_proftwp
    98      o       ,ht_proftwp,vt_proftwp,hq_proftwp,vq_proftwp)
     92        CALL interp_toga_time(daytime,day1,annee_ref                        &
     93     &       ,year_ini_twpi,day_ju_ini_twpi,nt_twpi,dt_twpi,nlev_twpi       &
     94     &       ,ts_twpi,plev_twpi,t_twpi,q_twpi,u_twpi,v_twpi,w_twpi          &
     95     &       ,ht_twpi,vt_twpi,hq_twpi,vq_twpi                               &
     96     &       ,ts_proftwp,plev_proftwp,t_proftwp,q_proftwp,u_proftwp         &
     97     &       ,v_proftwp,w_proftwp                                           &
     98     &       ,ht_proftwp,vt_proftwp,hq_proftwp,vq_proftwp)
    9999
    100100! vertical interpolation:
    101       CALL interp_toga_vertical(play,nlev_twpi,plev_proftwp
    102      :         ,t_proftwp,q_proftwp,u_proftwp,v_proftwp,w_proftwp
    103      :         ,ht_proftwp,vt_proftwp,hq_proftwp,vq_proftwp
    104      :         ,t_mod,q_mod,u_mod,v_mod,w_mod
    105      :         ,ht_mod,vt_mod,hq_mod,vq_mod,mxcalc)
     101      CALL interp_toga_vertical(play,nlev_twpi,plev_proftwp                 &
     102     &         ,t_proftwp,q_proftwp,u_proftwp,v_proftwp,w_proftwp           &
     103     &         ,ht_proftwp,vt_proftwp,hq_proftwp,vq_proftwp                 &
     104     &         ,t_mod,q_mod,u_mod,v_mod,w_mod                               &
     105     &         ,ht_mod,vt_mod,hq_mod,vq_mod,mxcalc)
    106106
    107107
    108108!calcul de l'advection verticale a partir du omega
    109 cCalcul des gradients verticaux
    110 cinitialisation
     109!Calcul des gradients verticaux
     110!initialisation
    111111      d_t_z(:)=0.
    112112      d_q_z(:)=0.
    … …  
    114114      d_q_dyn_z(:)=0.
    115115      DO l=2,llm-1
    116        d_t_z(l)=(temp(l+1)-temp(l-1))
    117      &          /(play(l+1)-play(l-1))
    118        d_q_z(l)=(q(l+1,1)-q(l-1,1))
    119      &          /(play(l+1)-play(l-1))
     116       d_t_z(l)=(temp(l+1)-temp(l-1))/(play(l+1)-play(l-1))
     117       d_q_z(l)=(q(l+1,1)-q(l-1,1))/(play(l+1)-play(l-1))
    120118      ENDDO
    121119      d_t_z(1)=d_t_z(2)
    … …  
    124122      d_q_z(llm)=d_q_z(llm-1)
    125123
    126 cCalcul de l advection verticale
     124!Calcul de l advection verticale
    127125      d_t_dyn_z(:)=w_mod(:)*d_t_z(:)
    128126      d_q_dyn_z(:)=w_mod(:)*d_q_z(:)
    … …  
    133131!           if (phi(l).gt.5000.) then
    134132        if (phi(l).gt.0.) then
    135         u(l)=u(l)
    136      .      +timestep*(u_mod(l)-u(l))/(2.*3600.)
    137         v(l)=v(l)
    138      .      +timestep*(v_mod(l)-v(l))/(2.*3600.)
     133        u(l)=u(l)+timestep*(u_mod(l)-u(l))/(2.*3600.)
     134        v(l)=v(l)+timestep*(v_mod(l)-v(l))/(2.*3600.)
    139135           endif   
    140136        else
    … …  
    150146           if ((zz(l).le.16000.).and.(zz(l).gt.15000.)) then
    151147             zfact=(zz(l)-15000.)/1000.
    152         q(l,1)=q(l,1)
    153      .      +timestep*(q_mod(l)-q(l,1))/(6.*3600.)*zfact
    154         temp(l)=temp(l)
    155      .      +timestep*(t_mod(l)-temp(l))/(6.*3600.)*zfact
     148        q(l,1)=q(l,1)+timestep*(q_mod(l)-q(l,1))/(6.*3600.)*zfact
     149        temp(l)=temp(l)+timestep*(t_mod(l)-temp(l))/(6.*3600.)*zfact
    156150           else if (zz(l).gt.16000.) then
    157         q(l,1)=q(l,1)
    158      .      +timestep*(q_mod(l)-q(l,1))/(6.*3600.)
    159         temp(l)=temp(l)
    160      .      +timestep*(t_mod(l)-temp(l))/(6.*3600.)
     151        q(l,1)=q(l,1)+timestep*(q_mod(l)-q(l,1))/(6.*3600.)
     152        temp(l)=temp(l)+timestep*(t_mod(l)-temp(l))/(6.*3600.)
    161153           endif
    162154        enddo   
    … …  
    189181       if (forcing_amma) then
    190182
    191         print*,
    192      : '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/dt_amma=',
    193      :    daytime,day1,(daytime-day1)*86400.,
    194      :    (daytime-day1)*86400/dt_amma
     183        print*,                                                             &
     184     & '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/dt_amma=',     &
     185     &    daytime,day1,(daytime-day1)*86400.,                               &
     186     &    (daytime-day1)*86400/dt_amma
    195187
    196188! time interpolation using TOGA interpolation routine
    197         CALL interp_amma_time(daytime,day1,annee_ref
    198      i       ,year_ini_amma,day_ju_ini_amma,nt_amma,dt_amma,nlev_amma
    199      i       ,vitw_amma,ht_amma,hq_amma,lat_amma,sens_amma
    200      o       ,vitw_profamma,ht_profamma,hq_profamma,lat_profamma
    201      :       ,sens_profamma)
     189        CALL interp_amma_time(daytime,day1,annee_ref                        &
     190     &       ,year_ini_amma,day_ju_ini_amma,nt_amma,dt_amma,nlev_amma       &
     191     &       ,vitw_amma,ht_amma,hq_amma,lat_amma,sens_amma                  &
     192     &       ,vitw_profamma,ht_profamma,hq_profamma,lat_profamma            &
     193     &       ,sens_profamma)
    202194
    203195      print*,'apres interpolation temporelle AMMA'
    … …  
    213205! vertical interpolation using TOGA interpolation routine:
    214206!      write(*,*)'avant interp vert', t_proftwp
    215       CALL interp_toga_vertical(play,nlev_amma,plev_amma
    216      :         ,th_profamma,q_profamma,u_profamma,v_profamma
    217      :         ,vitw_profamma
    218      :         ,ht_profamma,vt_profamma,hq_profamma,vq_profamma
    219      :         ,t_mod,q_mod,u_mod,v_mod,w_mod
    220      :         ,ht_mod,vt_mod,hq_mod,vq_mod,mxcalc)
     207      CALL interp_toga_vertical(play,nlev_amma,plev_amma                      &
     208     &         ,th_profamma,q_profamma,u_profamma,v_profamma                 &
     209     &         ,vitw_profamma                                               &
     210     &         ,ht_profamma,vt_profamma,hq_profamma,vq_profamma             &
     211     &         ,t_mod,q_mod,u_mod,v_mod,w_mod                               &
     212     &         ,ht_mod,vt_mod,hq_mod,vq_mod,mxcalc)
    221213       write(*,*) 'Profil initial forcing AMMA interpole'
    222214
    223215
    224216!calcul de l'advection verticale a partir du omega
    225 cCalcul des gradients verticaux
    226 cinitialisation
     217!Calcul des gradients verticaux
     218!initialisation
    227219      do l=1,llm
    228220      d_t_z(l)=0.
    … …  
    231223
    232224      DO l=2,llm-1
    233        d_t_z(l)=(temp(l+1)-temp(l-1))
    234      &          /(play(l+1)-play(l-1))
    235        d_q_z(l)=(q(l+1,1)-q(l-1,1))
    236      &          /(play(l+1)-play(l-1))
     225       d_t_z(l)=(temp(l+1)-temp(l-1))/(play(l+1)-play(l-1))
     226       d_q_z(l)=(q(l+1,1)-q(l-1,1))/(play(l+1)-play(l-1))
    237227      ENDDO
    238228      d_t_z(1)=d_t_z(2)
    … …  
    250240!        d_th_adv(l) = alpha*omega(l)/rcpd+ht_mod(l)-omega(l)*d_t_z(l)
    251241!attention: on impose dth
    252         d_th_adv(l) = alpha*omega(l)/rcpd+
     242        d_th_adv(l) = alpha*omega(l)/rcpd+                                  &
    253243     &         ht_mod(l)*(play(l)/pzero)**rkappa-omega(l)*d_t_z(l)
    254244!        d_th_adv(l) = 0.
    … …  
    275265!       call lstendH(llm,omega,dt_dyn,dq_dyn,du_dyn, dv_dyn,
    276266!     :  q,temp,u,v,play)
    277        call lstendH(llm,nqtot,omega,dt_dyn,dq_dyn,
    278      :  q,temp,u,v,play)
     267       call lstendH(llm,nqtot,omega,dt_dyn,dq_dyn,q,temp,u,v,play)
    279268
    280269        do l=1,llm
    … …  
    290279      if (forcing_armcu) then
    291280
    292         print*,
    293      : '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/dt_armcu=',
    294      :    day,day1,(day-day1)*86400.,(day-day1)*86400/dt_armcu
     281        print*,                                                             &
     282     & '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/dt_armcu=',    &
     283     &    day,day1,(day-day1)*86400.,(day-day1)*86400/dt_armcu
    295284
    296285! time interpolation:
    297286! ATTENTION, cet appel ne convient pas pour TOGA !!
    298287! revoir 1DUTILS.h et les arguments
    299       CALL interp_armcu_time(daytime,day1,annee_ref
    300      i            ,year_ini_armcu,day_ju_ini_armcu,nt_armcu,dt_armcu
    301      i            ,nlev_armcu,sens_armcu,flat_armcu,adv_theta_armcu
    302      i            ,rad_theta_armcu,adv_qt_armcu,sens_prof,flat_prof
    303      i            ,adv_theta_prof,rad_theta_prof,adv_qt_prof)
     288      CALL interp_armcu_time(daytime,day1,annee_ref                         &
     289     &            ,year_ini_armcu,day_ju_ini_armcu,nt_armcu,dt_armcu        &
     290     &            ,nlev_armcu,sens_armcu,flat_armcu,adv_theta_armcu          &
     291     &            ,rad_theta_armcu,adv_qt_armcu,sens_prof,flat_prof         &
     292     &            ,adv_theta_prof,rad_theta_prof,adv_qt_prof)
    304293
    305294! vertical interpolation:
    … …  
    346335      if (forcing_sandu) then
    347336
    348         print*,
    349      : '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/dt_sandu=',
    350      :    day,day1,(day-day1)*86400.,(day-day1)*86400/dt_sandu
     337        print*,                                                             &
     338     & '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/dt_sandu=',    &
     339     &    day,day1,(day-day1)*86400.,(day-day1)*86400/dt_sandu
    351340
    352341! time interpolation:
    353342! ATTENTION, cet appel ne convient pas pour TOGA !!
    354343! revoir 1DUTILS.h et les arguments
    355       CALL interp_sandu_time(daytime,day1,annee_ref
    356      i             ,year_ini_sandu,day_ju_ini_sandu,nt_sandu,dt_sandu
    357      i             ,nlev_sandu
    358      i             ,ts_sandu,ts_prof)
     344      CALL interp_sandu_time(daytime,day1,annee_ref                         &
     345     &             ,year_ini_sandu,day_ju_ini_sandu,nt_sandu,dt_sandu       &
     346     &             ,nlev_sandu                                              &
     347     &             ,ts_sandu,ts_prof)
    359348
    360349        if (type_ts_forcing.eq.1) ts_cur = ts_prof ! SST used in read_tsurf1d
    361350
    362351! vertical interpolation:
    363       CALL interp_sandu_vertical(play,nlev_sandu,plev_profs
    364      :         ,t_profs,thl_profs,q_profs,u_profs,v_profs,w_profs
    365      :         ,omega_profs,o3mmr_profs
    366      :         ,t_mod,thl_mod,q_mod,u_mod,v_mod,w_mod
    367      :         ,omega_mod,o3mmr_mod,mxcalc)
    368 !calcul de l'advection verticale
    369 cCalcul des gradients verticaux
    370 cinitialisation
    371       d_t_z(:)=0.
    372       d_q_z(:)=0.
    373       d_t_dyn_z(:)=0.
    374       d_q_dyn_z(:)=0.
    375 ! schema centre
    376 !     DO l=2,llm-1
    377 !      d_t_z(l)=(temp(l+1)-temp(l-1))
    378 !    &          /(play(l+1)-play(l-1))
    379 !      d_q_z(l)=(q(l+1,1)-q(l-1,1))
    380 !    &          /(play(l+1)-play(l-1))
    381 ! schema amont
    382       DO l=2,llm-1
    383        d_t_z(l)=(temp(l+1)-temp(l))/(play(l+1)-play(l))
    384        d_q_z(l)=(q(l+1,1)-q(l,1))/(play(l+1)-play(l))
    385 !     print *,'l temp2 temp0 play2 play0 omega_mod',
    386 !    & temp(l+1),temp(l-1),play(l+1),play(l-1),omega_mod(l)
    387       ENDDO
    388       d_t_z(1)=d_t_z(2)
    389       d_q_z(1)=d_q_z(2)
    390       d_t_z(llm)=d_t_z(llm-1)
    391       d_q_z(llm)=d_q_z(llm-1)
    392 
    393 !  calcul de l advection verticale
    394 ! Confusion w (m/s) et omega (Pa/s) !!
    395       d_t_dyn_z(:)=omega_mod(:)*d_t_z(:)
    396       d_q_dyn_z(:)=omega_mod(:)*d_q_z(:)
    397 !     do l=1,llm
    398 !      print *,'d_t_dyn omega_mod d_t_z d_q_dyn d_q_z',
    399 !    :l,d_t_dyn_z(l),omega_mod(l),d_t_z(l),d_q_dyn_z(l),d_q_z(l)
    400 !     enddo
    401 
    402 
    403 ! large-scale forcing : pour le cas Sandu ces forcages sont la SST
    404 ! et une divergence constante -> profil de omega
    405       tsurf = ts_prof
    406       write(*,*) 'SST suivante: ',tsurf
    407       do l = 1, llm
    408        omega(l) = omega_mod(l)
    409        omega2(l)= omega(l)/rg*airefi ! flxmass_w calcule comme ds physiq
    410 
    411        alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
    412 !
    413 !      d_th_adv(l) = 0.0
    414 !      d_q_adv(l,1) = 0.0
    415 !CR:test advection=0
    416 !calcul de l'advection verticale
    417         d_th_adv(l) = alpha*omega(l)/rcpd-d_t_dyn_z(l)
    418 !        print*,'temp adv',l,-d_t_dyn_z(l)
    419         d_q_adv(l,1) = -d_q_dyn_z(l)
    420 !        print*,'q adv',l,-d_q_dyn_z(l)
    421        dt_cooling(l) = 0.0
    422       enddo
    423       endif ! forcing_sandu
    424 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    425 !---------------------------------------------------------------------
    426 ! Interpolation forcing in time and onto model levels
    427 !---------------------------------------------------------------------
    428       if (forcing_astex) then
    429 
    430         print*,
    431      : '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/dt_astex=',
    432      :    day,day1,(day-day1)*86400.,(day-day1)*86400/dt_astex
    433 
    434 ! time interpolation:
    435 ! ATTENTION, cet appel ne convient pas pour TOGA !!
    436 ! revoir 1DUTILS.h et les arguments
    437       CALL interp_astex_time(daytime,day1,annee_ref
    438      i             ,year_ini_astex,day_ju_ini_astex,nt_astex,dt_astex
    439      i             ,nlev_astex,div_astex,ts_astex,ug_astex,vg_astex
    440      i             ,ufa_astex,vfa_astex,div_prof,ts_prof,ug_prof,vg_prof
    441      i             ,ufa_prof,vfa_prof)
    442 
    443         if (type_ts_forcing.eq.1) ts_cur = ts_prof ! SST used in read_tsurf1d
    444 
    445 ! vertical interpolation:
    446       CALL interp_astex_vertical(play,nlev_astex,plev_profa
    447      :         ,t_profa,thl_profa,qv_profa,ql_profa,qt_profa
    448      :         ,u_profa,v_profa,w_profa,tke_profa,o3mmr_profa
    449      :         ,t_mod,thl_mod,qv_mod,ql_mod,qt_mod,u_mod,v_mod,w_mod
    450      :         ,tke_mod,o3mmr_mod,mxcalc)
     352      CALL interp_sandu_vertical(play,nlev_sandu,plev_profs                 &
     353     &         ,t_profs,thl_profs,q_profs,u_profs,v_profs,w_profs           &
     354     &         ,omega_profs,o3mmr_profs                                     &
     355     &         ,t_mod,thl_mod,q_mod,u_mod,v_mod,w_mod                       &
     356     &         ,omega_mod,o3mmr_mod,mxcalc)
    451357!calcul de l'advection verticale
    452358!Calcul des gradients verticaux
    … …  
    476382!  calcul de l advection verticale
    477383! Confusion w (m/s) et omega (Pa/s) !!
    478       d_t_dyn_z(:)=w_mod(:)*d_t_z(:)
    479       d_q_dyn_z(:)=w_mod(:)*d_q_z(:)
     384      d_t_dyn_z(:)=omega_mod(:)*d_t_z(:)
     385      d_q_dyn_z(:)=omega_mod(:)*d_q_z(:)
    480386!     do l=1,llm
    481387!      print *,'d_t_dyn omega_mod d_t_z d_q_dyn d_q_z',
    … …  
    484390
    485391
    486 ! large-scale forcing : pour le cas Astex ces forcages sont la SST
    487 ! la divergence,ug,vg,ufa,vfa
     392! large-scale forcing : pour le cas Sandu ces forcages sont la SST
     393! et une divergence constante -> profil de omega
    488394      tsurf = ts_prof
    489395      write(*,*) 'SST suivante: ',tsurf
    490396      do l = 1, llm
    491        omega(l) = w_mod(l)
     397       omega(l) = omega_mod(l)
    492398       omega2(l)= omega(l)/rg*airefi ! flxmass_w calcule comme ds physiq
    493399
    … …  
    504410       dt_cooling(l) = 0.0
    505411      enddo
     412      endif ! forcing_sandu
     413!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     414!---------------------------------------------------------------------
     415! Interpolation forcing in time and onto model levels
     416!---------------------------------------------------------------------
     417      if (forcing_astex) then
     418
     419        print*,                                                             &
     420     & '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/dt_astex=',    &
     421     &    day,day1,(day-day1)*86400.,(day-day1)*86400/dt_astex
     422
     423! time interpolation:
     424! ATTENTION, cet appel ne convient pas pour TOGA !!
     425! revoir 1DUTILS.h et les arguments
     426      CALL interp_astex_time(daytime,day1,annee_ref                         &
     427     &             ,year_ini_astex,day_ju_ini_astex,nt_astex,dt_astex       &
     428     &             ,nlev_astex,div_astex,ts_astex,ug_astex,vg_astex         &
     429     &             ,ufa_astex,vfa_astex,div_prof,ts_prof,ug_prof,vg_prof    &
     430     &             ,ufa_prof,vfa_prof)
     431
     432        if (type_ts_forcing.eq.1) ts_cur = ts_prof ! SST used in read_tsurf1d
     433
     434! vertical interpolation:
     435      CALL interp_astex_vertical(play,nlev_astex,plev_profa                 &
     436     &         ,t_profa,thl_profa,qv_profa,ql_profa,qt_profa                &
     437     &         ,u_profa,v_profa,w_profa,tke_profa,o3mmr_profa               &
     438     &         ,t_mod,thl_mod,qv_mod,ql_mod,qt_mod,u_mod,v_mod,w_mod        &
     439     &         ,tke_mod,o3mmr_mod,mxcalc)
     440!calcul de l'advection verticale
     441!Calcul des gradients verticaux
     442!initialisation
     443      d_t_z(:)=0.
     444      d_q_z(:)=0.
     445      d_t_dyn_z(:)=0.
     446      d_q_dyn_z(:)=0.
     447! schema centre
     448!     DO l=2,llm-1
     449!      d_t_z(l)=(temp(l+1)-temp(l-1))
     450!    &          /(play(l+1)-play(l-1))
     451!      d_q_z(l)=(q(l+1,1)-q(l-1,1))
     452!    &          /(play(l+1)-play(l-1))
     453! schema amont
     454      DO l=2,llm-1
     455       d_t_z(l)=(temp(l+1)-temp(l))/(play(l+1)-play(l))
     456       d_q_z(l)=(q(l+1,1)-q(l,1))/(play(l+1)-play(l))
     457!     print *,'l temp2 temp0 play2 play0 omega_mod',
     458!    & temp(l+1),temp(l-1),play(l+1),play(l-1),omega_mod(l)
     459      ENDDO
     460      d_t_z(1)=d_t_z(2)
     461      d_q_z(1)=d_q_z(2)
     462      d_t_z(llm)=d_t_z(llm-1)
     463      d_q_z(llm)=d_q_z(llm-1)
     464
     465!  calcul de l advection verticale
     466! Confusion w (m/s) et omega (Pa/s) !!
     467      d_t_dyn_z(:)=w_mod(:)*d_t_z(:)
     468      d_q_dyn_z(:)=w_mod(:)*d_q_z(:)
     469!     do l=1,llm
     470!      print *,'d_t_dyn omega_mod d_t_z d_q_dyn d_q_z',
     471!    :l,d_t_dyn_z(l),omega_mod(l),d_t_z(l),d_q_dyn_z(l),d_q_z(l)
     472!     enddo
     473
     474
     475! large-scale forcing : pour le cas Astex ces forcages sont la SST
     476! la divergence,ug,vg,ufa,vfa
     477      tsurf = ts_prof
     478      write(*,*) 'SST suivante: ',tsurf
     479      do l = 1, llm
     480       omega(l) = w_mod(l)
     481       omega2(l)= omega(l)/rg*airefi ! flxmass_w calcule comme ds physiq
     482
     483       alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
     484!
     485!      d_th_adv(l) = 0.0
     486!      d_q_adv(l,1) = 0.0
     487!CR:test advection=0
     488!calcul de l'advection verticale
     489        d_th_adv(l) = alpha*omega(l)/rcpd-d_t_dyn_z(l)
     490!        print*,'temp adv',l,-d_t_dyn_z(l)
     491        d_q_adv(l,1) = -d_q_dyn_z(l)
     492!        print*,'q adv',l,-d_q_dyn_z(l)
     493       dt_cooling(l) = 0.0
     494      enddo
    506495      endif ! forcing_astex
    507496!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
Note: See TracChangeset for help on using the changeset viewer.