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