source: codes/icosagcm/trunk/src/dynetat0_hz_mod.f90 @ 150

Last change on this file since 150 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: 10.1 KB
Line 
1MODULE dynetat0_hz_mod 
2  USE genmod
3  USE icosa
4  USE caldyn_gcm_mod 
5  IMPLICIT NONE
6  PRIVATE
7
8       PUBLIC  etat0
9       INTEGER,SAVE::ncell
10       TYPE(t_field),POINTER:: f_iu(:)
11       TYPE(t_field),POINTER:: f_iv(:) 
12       TYPE(t_field),POINTER:: f_iue(:)
13       TYPE(t_field),POINTER:: f_ive(:) 
14       REAL(rstd),POINTER :: iu(:,:),iv(:,:)
15       REAL(rstd),POINTER :: iue(:,:),ive(:,:) 
16
17 
18CONTAINS
19
20    SUBROUTINE etat0(f_ps,f_phis,f_theta_rhodz,f_u,f_q) 
21  USE icosa
22  USE caldyn_mod
23  USE write_field 
24  USE maxicosa
25        IMPLICIT NONE
26     TYPE(t_domain),POINTER :: d 
27     TYPE(t_field),POINTER:: f_ps(:)
28     TYPE(t_field),POINTER:: f_phis(:)
29     TYPE(t_field),POINTER:: f_u(:)
30     TYPE(t_field),POINTER:: f_q(:)
31     TYPE(t_field),POINTER:: f_theta_rhodz(:) 
32     TYPE(t_field),POINTER::  f_buf_i3(:), f_buf1_i(:), f_buf2_i(:)
33     REAL(rstd),POINTER :: ps(:)
34     REAL(rstd),POINTER :: phis(:)
35     REAL(rstd),POINTER :: theta_rhodz(:,:)
36     REAL(rstd),POINTER :: u(:,:) 
37     REAL(rstd),POINTER :: q(:,:,:)
38     REAL(rstd):: maxff,minff,maxuu,minuu
39     INTEGER :: ind
40 
41     CALL allocate_field(f_iu,field_t,type_real,llm) 
42     CALL allocate_field(f_iv,field_t,type_real,llm)
43     CALL allocate_field(f_iue,field_u,type_real,llm) 
44     CALL allocate_field(f_ive,field_u,type_real,llm)
45     CALL allocate_field(f_u,field_u,type_real,llm) 
46     CALL allocate_field(f_buf1_i,field_t,type_real,llm) 
47     CALL allocate_field(f_buf2_i,field_t,type_real,llm) 
48     CALL allocate_field(f_buf_i3,field_u,type_real,3,llm) 
49 
50     PRINT*,"IN NETCDF READ"
51!------------------------------------zero
52        DO ind=1,ndomain
53      CALL swap_dimensions(ind)
54      CALL swap_geometry(ind)
55          iu = f_iu(ind) 
56          iv = f_iv(ind) 
57      iue = f_iue(ind) 
58         ive = f_ive(ind)       
59           u = f_u(ind) 
60          iu = 0.0 
61          iv = 0.0 
62           u = 0.0 
63         iue = 0.0     
64         ive = 0.0 
65     END DO 
66!--------------------------------------------
67     ncell = 0
68     DO ind=1,ndomain
69      CALL swap_dimensions(ind)
70      CALL swap_geometry(ind)
71      d => domain_glo(ind)
72      ps=f_ps(ind)
73      phis=f_phis(ind)
74      theta_rhodz=f_theta_rhodz(ind)
75      q=f_q(ind)
76      iu=f_iu(ind) 
77      iv=f_iv(ind) 
78      CALL compute_dynetat0(ind,d,ps,phis,theta_rhodz,iu,iv,q)
79      ENDDO
80
81     CALL transfert_request(f_ps,req_i1)
82        CALL transfert_request(f_phis,req_i1)
83        CALL transfert_request(f_theta_rhodz,req_i1)
84        CALL transfert_request(f_q,req_i1)
85        CALL transfert_request(f_iu,req_i1)
86        CALL transfert_request(f_iv,req_i1)
87!------------------------------------------
88        DO ind=1,ndomain
89      CALL swap_dimensions(ind)
90      CALL swap_geometry(ind)
91         u=f_u(ind)
92         iu=f_iu(ind) 
93         iv=f_iv(ind) 
94         iue=f_iue(ind) 
95         ive=f_ive(ind) 
96         CALL compute_dynetatu(iu,iv,iue,ive,u)
97        ENDDO
98!------------- OUTPUT OF Variables
99        CALL un2ulonlat(f_u,f_buf_i3,f_buf1_i,f_buf2_i) 
100        CALL writefield("buf1",f_buf1_i)
101  END SUBROUTINE etat0
102
103!==================================================================
104  SUBROUTINE compute_dynetat0(ind,d,ps,phis,theta_rhodz,iu,iv,q) 
105   use icosa
106   use netcdf
107   use wind_mod 
108   USE disvert_mod
109        IMPLICIT NONE
110   TYPE(t_domain),POINTER :: d 
111   CHARACTER*20::dimname 
112   REAL(rstd), INTENT(OUT) :: ps(iim*jjm)
113   REAL(rstd), INTENT(OUT) :: phis(iim*jjm)
114   REAL(rstd), INTENT(OUT) :: theta_rhodz(iim*jjm,llm)
115   REAL(rstd), INTENT(OUT) :: q(iim*jjm,llm,nqtot)
116   REAL(rstd),ALLOCATABLE :: mass(:,:)   ! mass   
117   REAL(rstd),ALLOCATABLE :: rhodz(:,:)   ! mass density 
118   REAL(rstd),ALLOCATABLE :: theta(:,:) 
119   REAL(rstd),ALLOCATABLE :: p(:,:)  ! pression
120   REAL(rstd),POINTER :: iu(:,:),iv(:,:)
121   REAL(rstd),POINTER :: icops(:)
122   REAL(rstd),POINTER :: icophis(:)
123   REAL(rstd),POINTER :: icou(:,:),icov(:,:)
124   REAL(rstd),POINTER :: icotheta(:,:)
125   REAL(rstd),POINTER :: icoq(:,:,:)
126   
127   INTEGER length,iq,ind,l
128   PARAMETER (length = 1)
129   CHARACTER(LEN=200):: iqq ! tableau des parametres du run
130   INTEGER::ierr,nid,ncid,nvarid,dimid,nind
131   INTEGER::ncells 
132   INTEGER::halo_size,i,j,k,ij
133   LOGICAL::single 
134   INTEGER::nDims,nVars,nGlobalAtts,unlimdimid
135   INTEGER:: len
136   
137!       OPEN NETCDF FILE
138         ierr = NF90_OPEN ("icosa_hz30.nc",NF90_NOWRITE,nid)
139      IF (ierr .NE. NF90_NOERR) THEN
140        write(*,*)'dynetat0: with file icosa_hz30.nc'
141        write(*,*)' ierr = ', ierr
142        STOP
143      ENDIF
144
145     ierr= nf90_inquire(nid,nDims,nVars,nGlobalAtts,unlimdimid)
146        IF (ierr .NE. NF90_NOERR) THEN
147        write(*,*)'Problem in inquire'
148        write(*,*)' ierr = ', ierr
149        STOP
150      ENDIF
151
152        PRINT*,"nDims,nVars,nGlobalAtts,unlimdimid"
153        PRINT*,nDims,nVars,nGlobalAtts,unlimdimid
154
155
156         ierr = NF90_INQ_DIMID(nid,"ncells",dimid)
157         IF (ierr .NE. NF90_NOERR ) THEN
158           write(*,*)'ncells is not present in hzdy_30.nc'
159        write(*,*)' ierr = ', ierr
160           STOP
161         ENDIF
162           
163         ierr = nf90_inquire_dimension(nid,dimid,dimname,ncells)
164         IF (ierr .NE. NF90_NOERR ) THEN
165           write(*,*)'ncells  in hzdy_30.nc'
166        write(*,*)' ierr = ', ierr
167           STOP
168         ENDIF
169
170          ALLOCATE(icops(ncells))
171          ALLOCATE(icophis(ncells))
172          ALLOCATE(icou(ncells,llm))
173          ALLOCATE(icov(ncells,llm))
174          ALLOCATE(icotheta(ncells,llm))
175          ALLOCATE(icoq(ncells,llm,nqtot))
176       ALLOCATE(p(iim*jjm,llm+1))
177          ALLOCATE(theta(iim*jjm,llm)) 
178       ALLOCATE(mass(iim*jjm,llm))   ! mass   
179       ALLOCATE(rhodz(iim*jjm,llm))   ! mass density   
180!============================================================
181      ierr = NF90_INQ_VARID(nid, "phis", nvarid)
182      IF (ierr .NE. NF90_NOERR) THEN
183        write(*,*)"dynetat0: phis is absent"
184           write(*,*)' ierr = ', ierr
185         STOP       
186      ENDIF
187
188      ierr = NF90_GET_VAR(nid, nvarid, icophis)
189      IF (ierr .NE. NF90_NOERR) THEN
190         write(*,*)"dynetat0: PROBLEM IN PHIS"
191         STOP
192      ENDIF
193!==============================================================
194          ierr = NF90_INQ_VARID(nid, "ps", nvarid)
195      IF (ierr .NE. NF90_NOERR) THEN
196        write(*,*)"dynetat0: ps is absent"
197           write(*,*)' ierr = ', ierr
198         STOP       
199      ENDIF
200
201      ierr = NF90_GET_VAR(nid, nvarid, icops)
202      IF (ierr .NE. NF90_NOERR) THEN
203         write(*,*)"dynetat0: PROBLEM IN PS"
204         STOP
205      ENDIF
206!================================================================
207          ierr = NF90_INQ_VARID(nid, "theta", nvarid)
208      IF (ierr .NE. NF90_NOERR) THEN
209        write(*,*)"dynetat0: teta is not available in start.nc"
210           write(*,*)' ierr = ', ierr
211         STOP       
212      ENDIF
213
214      ierr = NF90_GET_VAR(nid, nvarid,icotheta)
215      IF (ierr .NE. NF90_NOERR) THEN
216         write(*,*)"dynetat0: PROBLEM IN Teta"
217         STOP
218      ENDIF
219!================================================================
220        DO iq = 1,nqtot   
221                write(iqq,*)INT(iq)
222                iqq=ADJUSTL(iqq) 
223        ierr = NF90_INQ_VARID(nid,"q"//iqq, nvarid)
224      IF (ierr .NE. NF90_NOERR) THEN
225        write(*,*)"dynetat0: ","q"//iqq,"not here"
226           write(*,*)' ierr = ', ierr
227         STOP       
228      ENDIF
229
230      ierr = NF90_GET_VAR(nid, nvarid,icoq(:,:,iq))
231      IF (ierr .NE. NF90_NOERR) THEN
232         write(*,*)"dynetat0: PROBLEM IN Q"
233         STOP
234      ENDIF
235        END DO
236!===============================================================
237         ierr = NF90_INQ_VARID(nid, "ulon", nvarid)
238      IF (ierr .NE. NF90_NOERR) THEN
239        write(*,*)"dynetat0: ulon is not in file"
240           write(*,*)' ierr = ', ierr
241         STOP       
242      ENDIF
243       
244      ierr = NF90_GET_VAR(nid, nvarid,icou)
245      IF (ierr .NE. NF90_NOERR) THEN
246         write(*,*)"dynetat0: PROBLEM IN ucov"
247         STOP
248      ENDIF
249!================================================================
250          ierr = NF90_INQ_VARID(nid, "ulat", nvarid)
251      IF (ierr .NE. NF90_NOERR) THEN
252        write(*,*)"dynetat0: ulat is not available in start.nc"
253           write(*,*)' ierr = ', ierr
254         STOP       
255      ENDIF
256
257      ierr = NF90_GET_VAR(nid, nvarid,icov)
258      IF (ierr .NE. NF90_NOERR) THEN
259         write(*,*)"dynetat0: PROBLEM IN vlat"
260         STOP
261      ENDIF
262!================================================================
263                        iu = 0.0 ; iv = 0.0 
264             DO j=d%jj_begin,d%jj_end
265              DO i=d%ii_begin,d%ii_end
266                  k=d%iim*(j-1)+i
267                 IF (d%assign_domain(i,j)==ind ) THEN
268                     ncell=ncell+1
269                                 phis(k)= icophis(ncell) 
270                        ps(k)= icops(ncell) 
271                        theta(k,:) = icotheta(ncell,:) 
272                        q(k,:,:)= icoq(ncell,:,:) 
273                        iu(k,:) = icou(ncell,:)
274                                 iv(k,:) = icov(ncell,:)
275                           ENDIF
276                ENDDO
277              ENDDO
278
279    DO    l    = 1, llm+1
280      DO j=jj_begin,jj_end
281        DO i=ii_begin,ii_end
282          ij=(j-1)*iim+i
283          p(ij,l) = ap(l) + bp(l) * ps(ij)
284        ENDDO
285      ENDDO
286    ENDDO
287
288   DO l = 1, llm
289     DO j=jj_begin,jj_end
290       DO i=ii_begin,ii_end
291         ij=(j-1)*iim+i
292         mass(ij,l) = ( p(ij,l) - p(ij,l+1) )*Ai(ij)/g
293         rhodz(ij,l) = mass(ij,l) / Ai(ij)
294       ENDDO
295     ENDDO
296   ENDDO
297
298    DO    l    = 1, llm
299      DO j=jj_begin,jj_end
300        DO i=ii_begin,ii_end
301          ij=(j-1)*iim+i
302          theta_rhodz(ij,l) = theta(ij,l)*rhodz(ij,l)
303        ENDDO
304      ENDDO
305    ENDDO
306
307          DEALLOCATE(icops)
308          DEALLOCATE(icophis)
309          DEALLOCATE(icou)
310          DEALLOCATE(icov)
311          DEALLOCATE(icotheta)
312          DEALLOCATE(p)
313          DEALLOCATE(theta) 
314       DEALLOCATE(mass)   ! mass   
315       DEALLOCATE(rhodz)   !
316        END SUBROUTINE compute_dynetat0
317!==================================================================
318          SUBROUTINE compute_dynetatu(iu,iv,iue,ive,u) 
319   use icosa
320   use wind_mod 
321        IMPLICIT NONE
322   CHARACTER*20::dimname 
323   REAL(rstd),INTENT(OUT):: u(3*iim*jjm,llm)
324   REAL(rstd) :: iu(iim*jjm,llm),iv(iim*jjm,llm)
325   REAL(rstd) :: iue(3*iim*jjm,llm),ive(3*iim*jjm,llm)
326   INTEGER::halo_size,i,j,k,ij,l
327
328
329  Do l = 1, llm
330   DO j=jj_begin-1,jj_end+1
331      DO i=ii_begin-1,ii_end+1
332         k=iim*(j-1)+i
333           iue(k+u_right,l)=0.5*(iu(k,l)+iu(k+t_right,l))
334        iue(k+u_lup,l)  =0.5*(iu(k,l)+iu(k+t_lup,l))
335        iue(k+u_ldown,l)=0.5*(iu(k,l)+iu(k+t_ldown,l)) 
336!------------------------------------------------------
337           ive(k+u_right,l)=0.5*(iv(k,l)+iv(k+t_right,l))
338        ive(k+u_lup,l)  =0.5*(iv(k,l)+iv(k+t_lup,l))
339        ive(k+u_ldown,l)=0.5*(iv(k,l)+iv(k+t_ldown,l)) 
340         END DO
341    END DO
342  END DO
343        CALL compute_wind_perp_from_lonlat_compound(iue,ive,u) 
344
345        END SUBROUTINE compute_dynetatu
346
347
348 END MODULE dynetat0_hz_mod 
Note: See TracBrowser for help on using the repository browser.