source: LMDZ6/trunk/libf/phylmd/dyn1d/1D_interp_cases.h @ 3543

Last change on this file since 3543 was 3541, checked in by fhourdin, 5 years ago

Gros nettoyage en cours sur le 1D.
Le nouveau lmdz1d.F90 est une coquille vide qui choisit entre
old_lmdz1d.F90 (l'ancien lmdz1d.F90) et scm.F90 (le nouveau au format standard).
Plusieur fichiers sont renommés old_truc, le truc étant au format standard,
nettoyé des anciens formats.
Le 1DUTILS.h est lui même séparé entre des routines génériques venant remplacer
notamment des routines de dyn3d (la vocation d'origine de 1DUTILS.h) et
les routiles de lecture spécifiques mises dans old_1DUTILS.h
On perdra un peu de l'utilister de svn au moment de cette grosse bascule.
Mais les old_ sont faits pour ne plus bouger, et les versions standard
sont en pleine évolution.
Fredho

  • Property svn:keywords set to Id
File size: 7.5 KB
Line 
1
2         print*,'FORCING CASE forcing_case2'
3        print*,                                                             &
4     & '#### ITAP,day,day1,(day-day1)*86400,(day-day1)*86400/pdt_cas=',     &
5     &    daytime,day1,(daytime-day1)*86400.,                               &
6     &    (daytime-day1)*86400/pdt_cas
7
8! time interpolation:
9        CALL interp_case_time_std(daytime,day1,annee_ref                                       &
10!    &       ,year_ini_cas,day_ju_ini_cas,nt_cas,pdt_cas,nlev_cas                           &
11     &       ,nt_cas,nlev_cas                                                               &
12     &       ,ts_cas,ps_cas,plev_cas,t_cas,th_cas,thv_cas,thl_cas,qv_cas,ql_cas,qi_cas      &
13     &       ,u_cas,v_cas,ug_cas,vg_cas,vitw_cas,omega_cas,du_cas,hu_cas,vu_cas             &
14     &       ,dv_cas,hv_cas,vv_cas,dt_cas,ht_cas,vt_cas,dtrad_cas                           &
15     &       ,dq_cas,hq_cas,vq_cas,dth_cas,hth_cas,vth_cas,lat_cas,sens_cas,ustar_cas       &
16     &       ,uw_cas,vw_cas,q1_cas,q2_cas,tke_cas                                           &
17!
18     &       ,ts_prof_cas,plev_prof_cas,t_prof_cas,theta_prof_cas,thv_prof_cas  &
19     &       ,thl_prof_cas,qv_prof_cas,ql_prof_cas,qi_prof_cas                              &
20     &       ,u_prof_cas,v_prof_cas,ug_prof_cas,vg_prof_cas,vitw_prof_cas,omega_prof_cas    &
21     &       ,du_prof_cas,hu_prof_cas,vu_prof_cas                                           &
22     &       ,dv_prof_cas,hv_prof_cas,vv_prof_cas,dt_prof_cas,ht_prof_cas,vt_prof_cas       &
23     &       ,dtrad_prof_cas,dq_prof_cas,hq_prof_cas,vq_prof_cas                            &
24     &       ,dth_prof_cas,hth_prof_cas,vth_prof_cas,lat_prof_cas                           &
25     &       ,sens_prof_cas,ustar_prof_cas,uw_prof_cas,vw_prof_cas,q1_prof_cas,q2_prof_cas,tke_prof_cas)
26
27             ts_cur = ts_prof_cas
28!            psurf=plev_prof_cas(1)
29             psurf=ps_prof_cas
30
31! vertical interpolation:
32      CALL interp2_case_vertical(play,nlev_cas,plev_prof_cas                                              &
33     &         ,t_prof_cas,theta_prof_cas,thv_prof_cas,thl_prof_cas                                          &
34     &         ,qv_prof_cas,ql_prof_cas,qi_prof_cas,u_prof_cas,v_prof_cas                                 &
35     &         ,ug_prof_cas,vg_prof_cas,vitw_prof_cas,omega_prof_cas                                      &
36     &         ,du_prof_cas,hu_prof_cas,vu_prof_cas,dv_prof_cas,hv_prof_cas,vv_prof_cas                   &
37     &         ,dt_prof_cas,ht_prof_cas,vt_prof_cas,dtrad_prof_cas,dq_prof_cas,hq_prof_cas,vq_prof_cas    &
38     &         ,dth_prof_cas,hth_prof_cas,vth_prof_cas                                                    &
39!
40     &         ,t_mod_cas,theta_mod_cas,thv_mod_cas,thl_mod_cas,qv_mod_cas,ql_mod_cas,qi_mod_cas          &
41     &         ,u_mod_cas,v_mod_cas,ug_mod_cas,vg_mod_cas,w_mod_cas,omega_mod_cas                         &
42     &         ,du_mod_cas,hu_mod_cas,vu_mod_cas,dv_mod_cas,hv_mod_cas,vv_mod_cas                         &
43     &         ,dt_mod_cas,ht_mod_cas,vt_mod_cas,dtrad_mod_cas,dq_mod_cas,hq_mod_cas,vq_mod_cas           &
44     &         ,dth_mod_cas,hth_mod_cas,vth_mod_cas,mxcalc)
45
46
47      DO l=1,llm
48      teta(l)=temp(l)*(100000./play(l))**(rd/rcpd)
49      ENDDO
50!calcul de l'advection verticale a partir du omega
51!Calcul des gradients verticaux
52!initialisation
53      d_t_z(:)=0.
54      d_th_z(:)=0.
55      d_q_z(:)=0.
56      d_u_z(:)=0.
57      d_v_z(:)=0.
58      d_t_dyn_z(:)=0.
59      d_th_dyn_z(:)=0.
60      d_q_dyn_z(:)=0.
61      d_u_dyn_z(:)=0.
62      d_v_dyn_z(:)=0.
63      if (1==0) then
64         DO l=2,llm-1
65          d_t_z(l)=(temp(l+1)-temp(l-1))/(play(l+1)-play(l-1))
66          d_th_z(l)=(teta(l+1)-teta(l-1))/(play(l+1)-play(l-1))
67          d_q_z(l)=(q(l+1,1)-q(l-1,1))/(play(l+1)-play(l-1))
68          d_u_z(l)=(u(l+1)-u(l-1))/(play(l+1)-play(l-1))
69          d_v_z(l)=(v(l+1)-v(l-1))/(play(l+1)-play(l-1))
70         ENDDO
71      else
72         DO l=2,llm-1
73            IF (omega(l)>0.) THEN
74             d_t_z(l)=(temp(l+1)-temp(l))/(play(l+1)-play(l))
75             d_th_z(l)=(teta(l+1)-teta(l))/(play(l+1)-play(l))
76             d_q_z(l)=(q(l+1,1)-q(l,1))/(play(l+1)-play(l))
77             d_u_z(l)=(u(l+1)-u(l))/(play(l+1)-play(l))
78             d_v_z(l)=(v(l+1)-v(l))/(play(l+1)-play(l))
79            ELSE
80             d_t_z(l)=(temp(l-1)-temp(l))/(play(l-1)-play(l))
81             d_th_z(l)=(teta(l-1)-teta(l))/(play(l-1)-play(l))
82             d_q_z(l)=(q(l-1,1)-q(l,1))/(play(l-1)-play(l))
83             d_u_z(l)=(u(l-1)-u(l))/(play(l-1)-play(l))
84             d_v_z(l)=(v(l-1)-v(l))/(play(l-1)-play(l))
85            ENDIF
86         ENDDO
87      endif
88      d_t_z(1)=d_t_z(2)
89      d_t_z(1)=d_t_z(2)
90      d_th_z(1)=d_th_z(2)
91      d_q_z(1)=d_q_z(2)
92      d_u_z(1)=d_u_z(2)
93      d_v_z(1)=d_v_z(2)
94      d_t_z(llm)=d_t_z(llm-1)
95      d_th_z(llm)=d_th_z(llm-1)
96      d_q_z(llm)=d_q_z(llm-1)
97      d_u_z(llm)=d_u_z(llm-1)
98      d_v_z(llm)=d_v_z(llm-1)
99
100! TRAVAIL : PRENDRE DES NOTATIONS COHERENTES POUR W
101      do l = 1, llm
102! Modif w_mod_cas -> omega_mod_cas (MM+MPL 20170309)
103       omega(l) = -w_mod_cas(l)*play(l)*rg/(rd*temp(l))
104      enddo
105
106!Calcul de l advection verticale
107! Modif w_mod_cas -> omega_mod_cas (MM+MPL 20170310)
108      d_t_dyn_z(:)=omega(:)*d_t_z(:)
109      d_th_dyn_z(:)=omega(:)*d_th_z(:)
110      d_q_dyn_z(:)=omega(:)*d_q_z(:)
111      d_u_dyn_z(:)=omega(:)*d_u_z(:)
112      d_v_dyn_z(:)=omega(:)*d_v_z(:)
113
114!geostrophic wind
115      if (forc_geo.eq.1) then
116        do l=1,llm
117        ug(l) = ug_mod_cas(l) 
118        vg(l) = vg_mod_cas(l) 
119        enddo
120      endif
121!wind nudging
122      if (nudging_u.gt.0.) then
123        do l=1,llm
124           u(l)=u(l)+timestep*(u_mod_cas(l)-u(l))/(nudge_u)
125        enddo
126!     else
127!       do l=1,llm
128!          u(l) = u_mod_cas(l) 
129!       enddo
130      endif
131
132      if (nudging_v.gt.0.) then
133        do l=1,llm
134           v(l)=v(l)+timestep*(v_mod_cas(l)-v(l))/(nudge_v)
135        enddo
136!     else
137!       do l=1,llm
138!          v(l) = v_mod_cas(l) 
139!       enddo
140      endif
141
142      if (nudging_w.gt.0.) then
143        do l=1,llm
144           w(l)=w(l)+timestep*(w_mod_cas(l)-w(l))/(nudge_w)
145        enddo
146 !    else
147 !      do l=1,llm
148 !         w(l) = w_mod_cas(l) 
149 !      enddo
150      endif
151
152!nudging of q and temp
153      if (nudging_t.gt.0.) then
154        do l=1,llm
155           temp(l)=temp(l)+timestep*(t_mod_cas(l)-temp(l))/(nudge_t)
156        enddo
157      endif
158      if (nudging_q.gt.0.) then
159        do l=1,llm
160           q(l,1)=q(l,1)+timestep*(q_mod_cas(l)-q(l,1))/(nudge_q)
161        enddo
162      endif
163
164      do l = 1, llm
165! Modif w_mod_cas -> omega_mod_cas (MM+MPL 20170309)
166       omega2(l)= omega(l)/rg*airefi ! flxmass_w calcule comme ds physiq
167       alpha = rd*temp(l)*(1.+(rv/rd-1.)*q(l,1))/play(l)
168
169!calcul advections
170        d_u_adv(l)=du_mod_cas(l)
171        d_v_adv(l)=dv_mod_cas(l)
172        d_t_adv(l)=alpha*omega(l)/rcpd+dt_mod_cas(l)
173        d_q_adv(l,1)=dq_mod_cas(l)
174
175        if (forc_w==1) then
176           d_q_adv(l,1)=d_q_adv(l,1)-d_q_dyn_z(l)
177           d_t_adv(l)=d_t_adv(l)-d_t_dyn_z(l)
178           d_v_adv(l)=d_v_adv(l)-d_v_dyn_z(l)
179           d_u_adv(l)=d_u_adv(l)-d_u_dyn_z(l)
180        endif
181         
182        if (trad.eq.1) then
183           tend_rayo=1
184           dt_cooling(l) = dtrad_mod_cas(l)
185!          print *,'dt_cooling=',dt_cooling(l)
186        else
187           dt_cooling(l) = 0.0
188        endif
189      enddo
190
191! Faut-il multiplier par -1 ? (MPL 20160713)
192      IF(ok_flux_surf) THEN
193       fsens=-1.*sens_prof_cas
194       flat=-1.*lat_prof_cas
195       print *,'1D_interp: sens,flat',fsens,flat
196      ENDIF
197!
198      IF (ok_prescr_ust) THEN
199       ust=ustar_prof_cas
200       print *,'ust=',ust
201      ENDIF
Note: See TracBrowser for help on using the repository browser.