source: codes/icosagcm/trunk/src/physics_dry.f90 @ 149

Last change on this file since 149 was 149, checked in by sdubey, 11 years ago
Added few new routines to read NC files and compute diagnostics to r145.
Few routines of dry physics including radiation module, surface process and convective adjustment in new routine phyparam.f90. dynetat to read start files for dynamics. check_conserve routine to compute conservation of quatities like mass, energy etc.etat0_heldsz.f90 for held-suarez test case initial conditions. new Key time_style=lmd or dcmip to use day_step, ndays like in LMDZ
File size: 8.9 KB
Line 
1        MODULE physics_dry_mod
2  USE ICOSA
3  PUBLIC init_physics_dry, physics_dry
4
5CONTAINS
6
7   SUBROUTINE init_physics_dry
8      USE ICOSA !sarvesh
9      USE time_mod !sarvesh
10      USE dimphys_mod
11      USE RADIATION
12
13      IMPLICIT NONE
14      INTEGER::i,j,ij
15!-------------------------------------- ORBITAL PARAMETER----
16    periheli=150.
17    CALL getin('periheli', periheli)
18    aphelie=150.
19    CALL getin('aphelie',aphelie) 
20    coefir=0.08
21    CALL getin('coefir',coefir)
22    coefvis=0.99
23    CALL getin('coefvis',coefvis)
24    obliquit=0.0
25    CALL getin('obliquit',obliquit) 
26    peri_day=0.
27    CALL getin('peri_day',peri_day)
28    year_day=360.
29    CALL getin('year_day',year_day) 
30    callrad=.true.
31    CALL getin('callrad', callrad)
32    calldifv=.true.
33    CALL getin('calldifv', calldifv)
34    calladj=.true.
35    CALL getin('calladj', calladj)
36    callcond=.true.
37    callsoil=.true.
38    CALL getin('callsoil',callsoil)
39    season=.true.
40    CALL getin('season',season)
41    diurnal=.true.
42    CALL getin('diurnal',diurnal)
43    lverbose=.false.
44    CALL getin('lverbose',lverbose)
45    period_sort=1.
46   CALL getin('period_sort',period_sort)
47!    ptimestep=dt
48!   CALL getin('ptimestep',ptimestep)
49
50      print*,'Activation de la physique:'
51      print*,' Rayonnement ',callrad
52      print*,' Diffusion verticale turbulente ', calldifv
53      print*,' Ajustement convectif ',calladj
54      print*,' Sol ',callsoil
55      print*,' Cycle diurne ',diurnal
56!      choice of the frequency of the computation of radiations
57      IF(diurnal) THEN
58         iradia=NINT(daysec/(20.*dt))
59      ELSE
60         iradia=NINT(daysec/(4.*dt))
61      ENDIF
62      iradia=1
63
64      ngridmx=iim*jjm ; nlayermx=llm
65        offset=halo
66
67      ALLOCATE(albedo(ngridmx));ALLOCATE(emissiv(ngridmx))
68      ALLOCATE(inertie(ngridmx));ALLOCATE(z0(ngridmx))
69      ALLOCATE(rnatur(ngridmx));ALLOCATE(tsurf(ngridmx))
70      ALLOCATE(tsoil(ngridmx,nlayermx));ALLOCATE(fluxgrd(ngridmx))
71      ALLOCATE(fluxrad(ngridmx));ALLOCATE(dtrad(ngridmx,llm+1))
72      ALLOCATE(q2(ngridmx,llm+1));ALLOCATE(q2l(ngridmx,llm+1))
73      ALLOCATE(capcal(ngridmx))
74
75      CALL iniorbit(aphelie,periheli,year_day,peri_day,obliquit)
76
77      PRINT*,'unjours',daysec
78      PRINT*,'The radiative transfer is computed each ', iradia,' physical time-step or each ', &
79              iradia*dt,' seconds'
80  END SUBROUTINE init_physics_dry
81
82
83  SUBROUTINE physics_dry( it,jD_cur,jH_cur,f_phis, f_ps, f_theta_rhodz, f_ue, f_q)
84  USE icosa
85  IMPLICIT NONE
86    INTEGER,INTENT(IN)    :: it
87    REAL(rstd),INTENT(IN) :: jD_cur,jH_cur
88    TYPE(t_field),POINTER :: f_phis(:)
89    TYPE(t_field),POINTER :: f_ps(:)
90    TYPE(t_field),POINTER :: f_theta_rhodz(:)
91    TYPE(t_field),POINTER :: f_ue(:)
92    TYPE(t_field),POINTER :: f_q(:)
93 
94    REAL(rstd),POINTER :: phis(:)
95    REAL(rstd),POINTER :: ps(:)
96    REAL(rstd),POINTER :: theta_rhodz(:,:)
97    REAL(rstd),POINTER :: ue(:,:)
98    REAL(rstd),POINTER :: q(:,:,:)
99!    REAL(rstd),POINTER :: precl(:)
100    INTEGER :: ind
101!    LOGICAL:: firstcall,lastcall
102
103    CALL transfert_request(f_ue,req_e1_vect)
104    CALL transfert_request(f_theta_rhodz,req_i1)
105
106    DO ind=1,ndomain
107      CALL swap_dimensions(ind)
108      CALL swap_geometry(ind)
109      phis=f_phis(ind)
110      ps=f_ps(ind)
111      theta_rhodz=f_theta_rhodz(ind)
112      ue=f_ue(ind)
113      q=f_q(ind)
114!      out_i=f_out_i(ind)
115!      precl=f_precl(ind)
116!     print*,"====================================ind",ind,"----------it",it
117      CALL compute_physics_dry(it,jD_cur,jH_cur,phis, ps, theta_rhodz, ue, q(:,:,1))
118     ENDDO     
119
120!    CALL writefield("out_i",f_out_i)
121   
122!    IF (mod(it,itau_out)==0 ) THEN
123!      CALL writefield("precl",f_precl)
124!    ENDIF
125
126  END SUBROUTINE physics_dry
127 
128  SUBROUTINE compute_physics_dry(it,jD_cur,jH_cur,phis, ps, theta_rhodz, ue, q)
129  USE icosa
130  USE pression_mod
131  USE exner_mod
132  USE theta2theta_rhodz_mod
133  USE geopotential_mod
134  USE wind_mod
135  USE PHY
136 
137  IMPLICIT NONE
138    INTEGER::it
139    REAL(rstd) :: jD_cur
140    REAL(rstd) :: jH_cur
141    REAL(rstd) :: phis(iim*jjm)
142    REAL(rstd) :: ps(iim*jjm)
143    REAL(rstd) :: theta_rhodz(iim*jjm,llm)
144    REAL(rstd) :: ue(3*iim*jjm,llm)
145    REAL(rstd) :: q(iim*jjm,llm)
146!    REAL(rstd) :: precl(iim*jjm)
147
148    REAL(rstd) :: p(iim*jjm,llm+1)
149    REAL(rstd) :: pks(iim*jjm)
150    REAL(rstd) :: pk(iim*jjm,llm)
151    REAL(rstd) :: phi(iim*jjm,llm)
152    REAL(rstd) :: T(iim*jjm,llm)
153    REAL(rstd) :: Tfi(iim*jjm,llm)
154    REAL(rstd) :: theta(iim*jjm,llm)
155
156   REAL(rstd) :: uc(iim*jjm,3,llm)
157   REAL(rstd) :: u(iim*jjm,llm)
158   REAL(rstd) :: v(iim*jjm,llm)
159    REAL(rstd) :: ufi(iim*jjm,llm)
160    REAL(rstd) :: vfi(iim*jjm,llm)
161    REAL(rstd) :: qfi(iim*jjm,llm)
162    REAL(rstd) :: utemp(iim*jjm,llm)
163    REAL(rstd) :: vtemp(iim*jjm,llm)
164    REAL(rstd) :: lat(iim*jjm)
165    REAL(rstd) :: lon(iim*jjm)
166    REAL(rstd) :: pmid(iim*jjm,llm)
167    REAL(rstd) :: pint(iim*jjm,llm+1)
168    REAL(rstd) :: pdel(iim*jjm,llm)
169    REAL(rstd) :: plev(iim*jjm,llm+1),play(iim*jjm,llm) 
170    REAL(rstd) :: pkbycp
171    INTEGER :: i,j,l,ij,ig
172       
173!-------------------
174!    LOGICAL:: firstcall,lastcall
175    REAL(rstd) :: dufi(iim*jjm,llm)
176    REAL(rstd) :: dvfi(iim*jjm,llm)
177    REAL(rstd) :: dTfi(iim*jjm,llm)
178    REAL(rstd) :: dpsfi(iim*jjm)
179    REAL(rstd) :: dqfi(iim*jjm,llm)
180!    PRINT *,'Entering in LMD SIMPLE physics'   
181
182
183        offset=halo 
184    CALL compute_pression(ps,p,halo)
185    CALL compute_exner(ps,p,pks,pk,halo)
186    CALL compute_theta_rhodz2theta(ps,theta_rhodz,theta,halo)
187    CALL compute_geopotential(phis,pks,pk,theta,phi,halo)
188    CALL compute_theta_rhodz2temperature(ps,theta_rhodz,T,halo)
189    CALL compute_wind_centered(ue,uc)
190    CALL compute_wind_centered_lonlat_compound(uc, u, v)
191
192     DO j=jj_begin-offset,jj_end+offset
193      DO i=ii_begin-offset,ii_end+offset
194        ij=(j-1)*iim+i
195        CALL xyz2lonlat(xyz_i(ij,:),lon(ij),lat(ij)) 
196      ENDDO
197     ENDDO     
198   
199    DO l=1,llm
200      DO j=jj_begin-offset,jj_end+offset
201        DO i=ii_begin-offset,ii_end+offset
202          ij=(j-1)*iim+i
203!          Tfi(ij,l)=T(ij,l)
204!          ufi(ij,l)=u(ij,l)
205!          vfi(ij,l)=v(ij,l)
206!          qfi(ij,l)=q(ij,l)
207          dTfi(ij,l)=0.0
208          dufi(ij,l)=0.0
209          dvfi(ij,l)=0.0
210          dqfi(ij,l)=0.0
211        ENDDO
212      ENDDO
213    ENDDO   
214        plev(:,:) = p(:,:) 
215        dpsfi=0.0 
216
217   DO l=1,llm
218     DO j=jj_begin-offset,jj_end+offset
219       DO i=ii_begin-offset,ii_end+offset
220          ij=(j-1)*iim+i
221          pkbycp=pk(ij,l)/cpp
222          play(ij,l)=preff*pkbycp**(1./kappa) 
223        ENDDO
224      ENDDO
225     ENDDO 
226
227
228      CALL phyparam_lmd(it,iim*jjm,llm,1,dt,lat,lon,jD_cur,jH_cur, &
229      plev,play,phi,phis,u,v,T,q,dufi,dvfi,dTfi,dqfi,dpsfi) 
230
231      CALL ADDFI(u,v,T,q,ps,dufi,dvfi,dTfi,dqfi,dpsfi)
232
233!       CALL SARCHECKF(llm)
234!       print*,"plev",(maxval(plev(:,l)),l=1,llm+1)
235
236!      CALL phyparam_lmd(it,iim*jjm,llm,1,dt,lat,lon,jD_cur,jH_cur, &
237!      plev,play,phi,phis,ufi,vfi,Tfi,qfi,dufi,dvfi,dTfi,dqfi,dpsfi)
238
239!               Print*,"going ADD FI",it
240!      CALL ADDFI(ufi,vfi,Tfi,qfi,ps,dufi,dvfi,dTfi,dqfi,dpsfi)
241
242!       WRITE(11,*)"ducovfi",maxval(dufi),minval(dufi),it
243!       WRITE(11,*)"ucovfi",maxval(ufi),minval(ufi)
244!       WRITE(11,*)"dtetafi",maxval(dTfi),minval(dTfi)
245
246!=============================================
247!       go to 1234
248    DO l=1,llm
249      DO j=jj_begin-offset,jj_end+offset
250        DO i=ii_begin-offset,ii_end+offset
251          ij=(j-1)*iim+i
252          uc(ij,:,l)=(dufi(ij,l)*elon_i(ij,:)+dvfi(ij,l)*elat_i(ij,:))*dt
253        ENDDO
254      ENDDO
255    ENDDO
256
257   
258    DO l=1,llm
259      DO j=jj_begin-offset,jj_end+offset
260        DO i=ii_begin-offset,ii_end+offset
261          ij=(j-1)*iim+i
262          ue(ij+u_right,l)=ue(ij+u_right,l)+sum( 0.5*(uc(ij,:,l) + uc(ij+t_right,:,l))*ep_e(ij+u_right,:) )
263          ue(ij+u_lup,l)=ue(ij+u_lup,l)+sum( 0.5*(uc(ij,:,l) + uc(ij+t_lup,:,l))*ep_e(ij+u_lup,:) )
264          ue(ij+u_ldown,l)=ue(ij+u_ldown,l)+sum( 0.5*(uc(ij,:,l) + uc(ij+t_ldown,:,l))*ep_e(ij+u_ldown,:) )
265        ENDDO
266      ENDDO
267    ENDDO
268!1234    continue
269   
270!   CALL compute_temperature2theta_rhodz(ps,Tfi,theta_rhodz,halo)
271   CALL compute_temperature2theta_rhodz(ps,T,theta_rhodz,halo)
272
273!       WRITE(13,*)"tetafi",maxval(theta_rhodz),minval(theta_rhodz)
274    RETURN   
275  END SUBROUTINE compute_physics_dry
276
277    SUBROUTINE addfi(ufi,vfi,Tfi,qfi,ps,dufi,dvfi,dTfi,dqfi,dpsfi)
278        USE ICOSA
279        IMPLICIT NONE
280          REAL(rstd) :: dufi(iim*jjm,llm)
281          REAL(rstd) :: dvfi(iim*jjm,llm)
282        REAL(rstd) :: dTfi(iim*jjm,llm)
283        REAL(rstd) :: dpsfi(iim*jjm)
284        REAL(rstd) :: dqfi(iim*jjm,llm)
285        REAL(rstd) :: ufi(iim*jjm,llm)
286        REAL(rstd) :: vfi(iim*jjm,llm)
287        REAL(rstd) :: qfi(iim*jjm,llm)
288          REAL(rstd) :: ps(iim*jjm)
289          REAL(rstd) :: Tfi(iim*jjm,llm)
290          INTEGER::i,j,l,ij,offset
291                offset=halo
292       
293     DO l=1,llm
294      DO j=jj_begin-offset,jj_end+offset
295        DO i=ii_begin-offset,ii_end+offset
296          ij=(j-1)*iim+i
297          Tfi(ij,l)=Tfi(ij,l)+dTfi(ij,l)*dt
298          ufi(ij,l)=ufi(ij,l)+dufi(ij,l)*dt 
299          vfi(ij,l)=vfi(ij,l)+dvfi(ij,l)*dt 
300          qfi(ij,l)=qfi(ij,l)+dqfi(ij,l)*dt 
301        END DO
302      END DO
303     END DO 
304         ps(:)=ps(:) + dpsfi(:)*dt
305        END SUBROUTINE addfi 
306   
307END MODULE physics_dry_mod
Note: See TracBrowser for help on using the repository browser.