Context Navigation

iophy.F90 @ 1661

Last change on this file since 1661 was 1545, checked in by emillour, 9 years ago

Venus and Titan GCMs:

Adaptation wrt previous changes for Titan and Venus where
longitude and latitude arrays (in phycommon/geometry_mod) were overwritten
with values from startphy.nc files, where values are given in degrees.
For the sake of homegeneity with other physics package, revert to "default"
behaviour: longitude/latitude are in radians and longitude_deg/latitude_deg
are in degrees.
Also added checking in phyetat0 that the longitude/latitude read in the
restartphy.nc files match the ones provided by the dynamics.

File size: 13.1 KB

Line
1	!
2	! $Header$
3	!
4	module iophy
5
6	! abd REAL,private,allocatable,dimension(:),save :: io_lat
7	! abd REAL,private,allocatable,dimension(:),save :: io_lon
8	REAL,allocatable,dimension(:),save :: io_lat
9	REAL,allocatable,dimension(:),save :: io_lon
10	INTEGER, save :: phys_domain_id
11	INTEGER, save :: npstn
12	INTEGER, allocatable, dimension(:), save :: nptabij
13
14	INTERFACE histwrite_phy
15	MODULE PROCEDURE histwrite2d_phy,histwrite3d_phy
16	END INTERFACE
17
18	INTERFACE histbeg_phy_all
19	MODULE PROCEDURE histbeg_phy,histbeg_phy_points
20	END INTERFACE
21
22
23	contains
24
25	subroutine init_iophy_new(rlat,rlon)
26	USE dimphy, only: klon
27	USE mod_phys_lmdz_para
28	USE mod_grid_phy_lmdz, only: nbp_lon, nbp_lat, klon_glo
29	USE ioipsl
30	implicit none
31	real,dimension(klon),intent(in) :: rlon ! longitudes, in degrees
32	real,dimension(klon),intent(in) :: rlat ! latitudes, in degrees
33
34	REAL,dimension(klon_glo) :: rlat_glo
35	REAL,dimension(klon_glo) :: rlon_glo
36
37	INTEGER,DIMENSION(2) :: ddid
38	INTEGER,DIMENSION(2) :: dsg
39	INTEGER,DIMENSION(2) :: dsl
40	INTEGER,DIMENSION(2) :: dpf
41	INTEGER,DIMENSION(2) :: dpl
42	INTEGER,DIMENSION(2) :: dhs
43	INTEGER,DIMENSION(2) :: dhe
44	INTEGER :: i
45
46	CALL gather(rlat,rlat_glo)
47	CALL bcast(rlat_glo)
48	CALL gather(rlon,rlon_glo)
49	CALL bcast(rlon_glo)
50
51	!$OMP MASTER
52	ALLOCATE(io_lat(nbp_lat))
53	IF (klon_glo == 1) THEN
54	io_lat(1)=rlat_glo(1)
55	ELSE
56	io_lat(1)=rlat_glo(1)
57	io_lat(nbp_lat)=rlat_glo(klon_glo)
58	DO i=2,nbp_lat-1
59	io_lat(i)=rlat_glo(2+(i-2)*nbp_lon)
60	ENDDO
61	ENDIF
62
63	ALLOCATE(io_lon(nbp_lon))
64	IF (klon_glo == 1) THEN
65	io_lon(1)=rlon_glo(1)
66	ELSE
67	io_lon(1:nbp_lon)=rlon_glo(2:nbp_lon+1)
68	ENDIF
69
70	! POUR VENUS, LA CARTE EST A L ENVERS !!
71	io_lat = -1.*io_lat
72	io_lon = -1.*io_lon
73
74	ddid=(/ 1,2 /)
75	dsg=(/ nbp_lon, nbp_lat /)
76	dsl=(/ nbp_lon, jj_nb /)
77	dpf=(/ 1,jj_begin /)
78	dpl=(/ nbp_lon, jj_end /)
79	dhs=(/ ii_begin-1,0 /)
80	IF (mpi_rank==mpi_size-1) THEN
81	dhe=(/0,0/)
82	ELSE
83	dhe=(/ nbp_lon-ii_end,0 /)
84	ENDIF
85
86	call flio_dom_set(mpi_size,mpi_rank,ddid,dsg,dsl,dpf,dpl,dhs,dhe, &
87	'APPLE',phys_domain_id)
88
89	!$OMP END MASTER
90
91	end subroutine init_iophy_new
92
93	subroutine init_iophy(lat,lon)
94	USE dimphy
95	USE mod_phys_lmdz_para
96	use ioipsl, only: flio_dom_set
97	USE mod_grid_phy_lmdz, ONLY : nbp_lon, nbp_lat
98	implicit none
99	real,dimension(nbp_lon),intent(in) :: lon
100	real,dimension(nbp_lat),intent(in) :: lat
101
102	INTEGER,DIMENSION(2) :: ddid
103	INTEGER,DIMENSION(2) :: dsg
104	INTEGER,DIMENSION(2) :: dsl
105	INTEGER,DIMENSION(2) :: dpf
106	INTEGER,DIMENSION(2) :: dpl
107	INTEGER,DIMENSION(2) :: dhs
108	INTEGER,DIMENSION(2) :: dhe
109
110	!$OMP MASTER
111	allocate(io_lat(nbp_lat))
112	io_lat(:)=lat(:)
113	allocate(io_lon(nbp_lon))
114	io_lon(:)=lon(:)
115
116	ddid=(/ 1,2 /)
117	dsg=(/ nbp_lon, nbp_lat /)
118	dsl=(/ nbp_lon, jj_nb /)
119	dpf=(/ 1,jj_begin /)
120	dpl=(/ nbp_lon, jj_end /)
121	dhs=(/ ii_begin-1,0 /)
122	if (mpi_rank==mpi_size-1) then
123	dhe=(/0,0/)
124	else
125	dhe=(/ nbp_lon-ii_end,0 /)
126	endif
127
128	call flio_dom_set(mpi_size,mpi_rank,ddid,dsg,dsl,dpf,dpl,dhs,dhe, &
129	'APPLE',phys_domain_id)
130
131	!$OMP END MASTER
132
133	end subroutine init_iophy
134
135	subroutine histbeg_phy(name,itau0,zjulian,dtime,nhori,nid_day)
136	USE mod_phys_lmdz_para, only: jj_begin, jj_end, jj_nb, is_sequential
137	USE mod_grid_phy_lmdz, ONLY : nbp_lon, nbp_lat
138	use ioipsl, only: histbeg
139	use write_field
140	implicit none
141
142	character(), intent(IN) :: name
143	integer, intent(in) :: itau0
144	real,intent(in) :: zjulian
145	real,intent(in) :: dtime
146	integer,intent(out) :: nhori
147	integer,intent(out) :: nid_day
148
149	!$OMP MASTER
150	if (is_sequential) then
151	call histbeg(name,nbp_lon,io_lon, jj_nb,io_lat(jj_begin:jj_end), &
152	1,nbp_lon,1,jj_nb,itau0, zjulian, dtime, nhori, nid_day)
153	else
154	call histbeg(name,nbp_lon,io_lon, jj_nb,io_lat(jj_begin:jj_end), &
155	1,nbp_lon,1,jj_nb,itau0, zjulian, dtime, nhori, nid_day,phys_domain_id)
156	endif
157	!$OMP END MASTER
158
159	end subroutine histbeg_phy
160
161	subroutine histbeg_phy_points(rlon,rlat,pim,tabij,ipt,jpt, &
162	plon,plat,plon_bounds,plat_bounds, &
163	nname,itau0,zjulian,dtime,nnhori,nnid_day)
164	USE dimphy, only: klon
165	USE mod_phys_lmdz_para
166	USE mod_grid_phy_lmdz, only: klon_glo, nbp_lon, nbp_lat
167	use ioipsl, only: histbeg
168
169	implicit none
170
171	real,dimension(klon),intent(in) :: rlon
172	real,dimension(klon),intent(in) :: rlat
173	integer, intent(in) :: itau0
174	real,intent(in) :: zjulian
175	real,intent(in) :: dtime
176	integer, intent(in) :: pim
177	integer, intent(out) :: nnhori
178	character(len=20), intent(in) :: nname
179	INTEGER, intent(out) :: nnid_day
180	integer :: i
181	REAL,dimension(klon_glo) :: rlat_glo
182	REAL,dimension(klon_glo) :: rlon_glo
183	INTEGER, DIMENSION(pim), intent(in) :: tabij
184	REAL,dimension(pim), intent(in) :: plat, plon
185	INTEGER,dimension(pim), intent(in) :: ipt, jpt
186	REAL,dimension(pim,2), intent(out) :: plat_bounds, plon_bounds
187
188	INTEGER, SAVE :: tabprocbeg, tabprocend
189	!$OMP THREADPRIVATE(tabprocbeg, tabprocend)
190	INTEGER :: ip
191	INTEGER, PARAMETER :: nip=1
192	INTEGER :: npproc
193	REAL, allocatable, dimension(:) :: npplat, npplon
194	REAL, allocatable, dimension(:,:) :: npplat_bounds, npplon_bounds
195	REAL, dimension(nbp_lon,nbp_lat) :: zx_lon, zx_lat
196
197	CALL gather(rlat,rlat_glo)
198	CALL bcast(rlat_glo)
199	CALL gather(rlon,rlon_glo)
200	CALL bcast(rlon_glo)
201
202	!$OMP MASTER
203	DO i=1,pim
204
205	! print*,'CFMIP_iophy i tabij lon lat',i,tabij(i),plon(i),plat(i)
206
207	plon_bounds(i,1)=rlon_glo(tabij(i)-1)
208	plon_bounds(i,2)=rlon_glo(tabij(i)+1)
209	if(plon_bounds(i,2).LE.0..AND.plon_bounds(i,1).GE.0.) THEN
210	if(rlon_glo(tabij(i)).GE.0.) THEN
211	plon_bounds(i,2)=-1*plon_bounds(i,2)
212	endif
213	endif
214	if(plon_bounds(i,2).GE.0..AND.plon_bounds(i,1).LE.0.) THEN
215	if(rlon_glo(tabij(i)).LE.0.) THEN
216	plon_bounds(i,2)=-1*plon_bounds(i,2)
217	endif
218	endif
219	!
220	IF ( tabij(i).LE.nbp_lon) THEN
221	plat_bounds(i,1)=rlat_glo(tabij(i))
222	ELSE
223	plat_bounds(i,1)=rlat_glo(tabij(i)-nbp_lon)
224	ENDIF
225	plat_bounds(i,2)=rlat_glo(tabij(i)+nbp_lon)
226	!
227	! print*,'CFMIP_iophy point i lon lon_bds',i,plon_bounds(i,1),rlon_glo(tabij(i)),plon_bounds(i,2)
228	! print*,'CFMIP_iophy point i lat lat_bds',i,plat_bounds(i,1),rlat_glo(tabij(i)),plat_bounds(i,2)
229	!
230	ENDDO
231	if (is_sequential) then
232
233	npstn=pim
234	IF(.NOT. ALLOCATED(nptabij)) THEN
235	ALLOCATE(nptabij(pim))
236	ENDIF
237	DO i=1,pim
238	nptabij(i)=tabij(i)
239	ENDDO
240
241	CALL gr_fi_ecrit(1,klon,nbp_lon,nbp_lat,rlon_glo,zx_lon)
242	if ((nbp_lon*nbp_lat).gt.1) then
243	DO i = 1, nbp_lon
244	zx_lon(i,1) = rlon_glo(i+1)
245	zx_lon(i,nbp_lat) = rlon_glo(i+1)
246	ENDDO
247	endif
248	CALL gr_fi_ecrit(1,klon,nbp_lon,nbp_lat,rlat_glo,zx_lat)
249
250	DO i=1,pim
251	! print*,'CFMIP_iophy i tabij lon lat',i,tabij(i),plon(i),plat(i)
252
253	plon_bounds(i,1)=zx_lon(ipt(i)-1,jpt(i))
254	plon_bounds(i,2)=zx_lon(ipt(i)+1,jpt(i))
255
256	if (ipt(i).EQ.1) then
257	plon_bounds(i,1)=zx_lon(nbp_lon,jpt(i))
258	plon_bounds(i,2)=360.+zx_lon(ipt(i)+1,jpt(i))
259	endif
260
261	if (ipt(i).EQ.nbp_lon) then
262	plon_bounds(i,2)=360.+zx_lon(1,jpt(i))
263	endif
264
265	plat_bounds(i,1)=zx_lat(ipt(i),jpt(i)-1)
266	plat_bounds(i,2)=zx_lat(ipt(i),jpt(i)+1)
267
268	if (jpt(i).EQ.1) then
269	plat_bounds(i,1)=zx_lat(ipt(i),1)+0.001
270	plat_bounds(i,2)=zx_lat(ipt(i),1)-0.001
271	endif
272
273	if (jpt(i).EQ.nbp_lat) then
274	plat_bounds(i,1)=zx_lat(ipt(i),nbp_lat)+0.001
275	plat_bounds(i,2)=zx_lat(ipt(i),nbp_lat)-0.001
276	endif
277	!
278	! print*,'CFMIP_iophy point i lon lon_bds',i,plon_bounds(i,1),rlon(tabij(i)),plon_bounds(i,2)
279	! print*,'CFMIP_iophy point i lat lat_bds',i,plat_bounds(i,1),rlat(tabij(i)),plat_bounds(i,2)
280	!
281	ENDDO
282	! print*,'iophy is_sequential nname, nnhori, nnid_day=',trim(nname),nnhori,nnid_day
283	call histbeg(nname,pim,plon,plon_bounds, &
284	plat,plat_bounds, &
285	itau0, zjulian, dtime, nnhori, nnid_day)
286	else
287	npproc=0
288	DO ip=1, pim
289	tabprocbeg=klon_mpi_begin
290	tabprocend=klon_mpi_end
291	IF(tabij(ip).GE.tabprocbeg.AND.tabij(ip).LE.tabprocend) THEN
292	npproc=npproc+1
293	npstn=npproc
294	ENDIF
295	ENDDO
296	! print*,'CFMIP_iophy mpi_rank npstn',mpi_rank,npstn
297	IF(.NOT. ALLOCATED(nptabij)) THEN
298	ALLOCATE(nptabij(npstn))
299	ALLOCATE(npplon(npstn), npplat(npstn))
300	ALLOCATE(npplon_bounds(npstn,2), npplat_bounds(npstn,2))
301	ENDIF
302	npproc=0
303	DO ip=1, pim
304	IF(tabij(ip).GE.tabprocbeg.AND.tabij(ip).LE.tabprocend) THEN
305	npproc=npproc+1
306	nptabij(npproc)=tabij(ip)
307	! print*,'mpi_rank npproc ip plon plat tabij=',mpi_rank,npproc,ip, &
308	! plon(ip),plat(ip),tabij(ip)
309	npplon(npproc)=plon(ip)
310	npplat(npproc)=plat(ip)
311	npplon_bounds(npproc,1)=plon_bounds(ip,1)
312	npplon_bounds(npproc,2)=plon_bounds(ip,2)
313	npplat_bounds(npproc,1)=plat_bounds(ip,1)
314	npplat_bounds(npproc,2)=plat_bounds(ip,2)
315	!!!
316	!!! print qui sert a reordonner les points stations selon l'ordre CFMIP
317	!!! ne pas enlever
318	print*,'iophy_mpi rank ip lon lat',mpi_rank,ip,plon(ip),plat(ip)
319	!!!
320	ENDIF
321	ENDDO
322	call histbeg(nname,npstn,npplon,npplon_bounds, &
323	npplat,npplat_bounds, &
324	itau0,zjulian,dtime,nnhori,nnid_day,phys_domain_id)
325	endif
326	!$OMP END MASTER
327
328	end subroutine histbeg_phy_points
329
330	subroutine histwrite2d_phy(nid,lpoint,name,itau,field)
331	USE dimphy, only: klon
332	USE mod_phys_lmdz_para
333	USE ioipsl, only: histwrite
334	USE mod_grid_phy_lmdz, ONLY : nbp_lon, nbp_lat
335	implicit none
336
337	integer,intent(in) :: nid
338	logical,intent(in) :: lpoint
339	character(), intent(IN) :: name
340	integer, intent(in) :: itau
341	real,dimension(:),intent(in) :: field
342	REAL,dimension(klon_mpi) :: buffer_omp
343	INTEGER, allocatable, dimension(:) :: index2d
344	REAL :: Field2d(nbp_lon,jj_nb)
345
346	integer :: ip
347	real,allocatable,dimension(:) :: fieldok
348
349	IF (size(field)/=klon) CALL abort_gcm('iophy::histwrite2d','Field first dimension not equal to klon',1)
350
351	CALL Gather_omp(field,buffer_omp)
352	!$OMP MASTER
353	CALL grid1Dto2D_mpi(buffer_omp,Field2d)
354	if(.NOT.lpoint) THEN
355	ALLOCATE(index2d(nbp_lon*jj_nb))
356	ALLOCATE(fieldok(nbp_lon*jj_nb))
357	CALL histwrite(nid,name,itau,Field2d,nbp_lon*jj_nb,index2d)
358	else
359	ALLOCATE(fieldok(npstn))
360	ALLOCATE(index2d(npstn))
361
362	if(is_sequential) then
363	! klon_mpi_begin=1
364	! klon_mpi_end=klon
365	DO ip=1, npstn
366	fieldok(ip)=buffer_omp(nptabij(ip))
367	ENDDO
368	else
369	DO ip=1, npstn
370	! print*,'histwrite2d is_sequential npstn ip name nptabij',npstn,ip,name,nptabij(ip)
371	IF(nptabij(ip).GE.klon_mpi_begin.AND. &
372	nptabij(ip).LE.klon_mpi_end) THEN
373	fieldok(ip)=buffer_omp(nptabij(ip)-klon_mpi_begin+1)
374	ENDIF
375	ENDDO
376	endif
377	CALL histwrite(nid,name,itau,fieldok,npstn,index2d)
378	!
379	endif
380	deallocate(index2d)
381	deallocate(fieldok)
382	!$OMP END MASTER
383	end subroutine histwrite2d_phy
384
385	subroutine histwrite3d_phy(nid,lpoint,name,itau,field)
386	USE dimphy, only: klon
387	USE mod_phys_lmdz_para
388	USE mod_grid_phy_lmdz, ONLY : nbp_lon, nbp_lat
389	use ioipsl, only: histwrite
390	implicit none
391
392	integer,intent(in) :: nid
393	logical,intent(in) :: lpoint
394	character(), intent(IN) :: name
395	integer, intent(in) :: itau
396	real,dimension(:,:),intent(in) :: field ! --> field(klon,:)
397	REAL,dimension(klon_mpi,size(field,2)) :: buffer_omp
398	REAL :: Field3d(nbp_lon,jj_nb,size(field,2))
399	INTEGER :: ip, n, nlev
400	INTEGER, ALLOCATABLE, dimension(:) :: index3d
401	real,allocatable, dimension(:,:) :: fieldok
402
403	IF (size(field,1)/=klon) CALL abort_gcm('iophy::histwrite3d','Field first dimension not equal to klon',1)
404	nlev=size(field,2)
405
406	! print*,'hist3d_phy mpi_rank npstn=',mpi_rank,npstn
407
408	! DO ip=1, npstn
409	! print*,'hist3d_phy mpi_rank nptabij',mpi_rank,nptabij(ip)
410	! ENDDO
411
412	CALL Gather_omp(field,buffer_omp)
413	!$OMP MASTER
414	CALL grid1Dto2D_mpi(buffer_omp,field3d)
415	if(.NOT.lpoint) THEN
416	ALLOCATE(index3d(nbp_lonjj_nbnlev))
417	ALLOCATE(fieldok(nbp_lon*jj_nb,nlev))
418	CALL histwrite(nid,name,itau,Field3d,nbp_lonjj_nbnlev,index3d)
419	else
420	nlev=size(field,2)
421	ALLOCATE(index3d(npstn*nlev))
422	ALLOCATE(fieldok(npstn,nlev))
423
424	if(is_sequential) then
425	! klon_mpi_begin=1
426	! klon_mpi_end=klon
427	DO n=1, nlev
428	DO ip=1, npstn
429	fieldok(ip,n)=buffer_omp(nptabij(ip),n)
430	ENDDO
431	ENDDO
432	else
433	DO n=1, nlev
434	DO ip=1, npstn
435	IF(nptabij(ip).GE.klon_mpi_begin.AND. &
436	nptabij(ip).LE.klon_mpi_end) THEN
437	fieldok(ip,n)=buffer_omp(nptabij(ip)-klon_mpi_begin+1,n)
438	ENDIF
439	ENDDO
440	ENDDO
441	endif
442	CALL histwrite(nid,name,itau,fieldok,npstn*nlev,index3d)
443	endif
444	deallocate(index3d)
445	deallocate(fieldok)
446	!$OMP END MASTER
447	end subroutine histwrite3d_phy
448
449	end module iophy

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Original Format