source: trunk/LMDZ.UNIVERSAL/libf/phygeneric/iophy.F90 @ 958

Last change on this file since 958 was 870, checked in by aslmd, 12 years ago

LMDZ.UNIVERSAL. Added compiling for large domains on Gnome. Wrote a better README.
File size: 12.5 KB
Line 
1!
2! $Header$
3!
4module 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
23contains
24
25  subroutine init_iophy_new(rlat,rlon)
26  USE dimphy
27  USE mod_phys_lmdz_para
28  USE mod_grid_phy_lmdz
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
86  USE mod_phys_lmdz_para
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
128  USE mod_phys_lmdz_para
129  use ioipsl
130  implicit none
131  include 'dimensions.h'
132   
133    character*(*), intent(IN) :: name
134    integer, intent(in) :: itau0
135    real,intent(in) :: zjulian
136    real,intent(in) :: dtime
137    integer,intent(out) :: nhori
138    integer,intent(out) :: nid_day
139
140!$OMP MASTER   
141    if (is_sequential) then
142      call histbeg(name,iim,io_lon, jj_nb,io_lat(jj_begin:jj_end), &
143                   1,iim,1,jj_nb,itau0, zjulian, dtime, nhori, nid_day)
144    else
145      call histbeg(name,iim,io_lon, jj_nb,io_lat(jj_begin:jj_end), &
146                   1,iim,1,jj_nb,itau0, zjulian, dtime, nhori, nid_day,phys_domain_id)
147    endif
148!$OMP END MASTER
149 
150  end subroutine histbeg_phy
151
152  subroutine histbeg_phy_points(rlon,rlat,pim,tabij,ipt,jpt, &
153             plon,plat,plon_bounds,plat_bounds, &
154             nname,itau0,zjulian,dtime,nnhori,nnid_day)
155  USE dimphy
156  USE mod_phys_lmdz_para
157  USE mod_grid_phy_lmdz
158  use ioipsl
159  implicit none
160  include 'dimensions.h'
161
162    real,dimension(klon),intent(in) :: rlon
163    real,dimension(klon),intent(in) :: rlat
164    integer, intent(in) :: itau0
165    real,intent(in) :: zjulian
166    real,intent(in) :: dtime
167    integer, intent(in) :: pim
168    integer, intent(out) :: nnhori
169    character(len=20), intent(in) :: nname
170    INTEGER, intent(out) :: nnid_day
171    integer :: i
172    REAL,dimension(klon_glo)        :: rlat_glo
173    REAL,dimension(klon_glo)        :: rlon_glo
174    INTEGER, DIMENSION(pim), intent(in)  :: tabij
175    REAL,dimension(pim), intent(in) :: plat, plon
176    INTEGER,dimension(pim), intent(in) :: ipt, jpt
177    REAL,dimension(pim,2), intent(out) :: plat_bounds, plon_bounds
178
179    INTEGER, SAVE :: tabprocbeg, tabprocend
180!$OMP THREADPRIVATE(tabprocbeg, tabprocend)
181    INTEGER :: ip
182    INTEGER, PARAMETER :: nip=1
183    INTEGER :: npproc
184    REAL, allocatable, dimension(:) :: npplat, npplon
185    REAL, allocatable, dimension(:,:) :: npplat_bounds, npplon_bounds
186    INTEGER, PARAMETER :: jjmp1=jjm+1-1/jjm
187    REAL, dimension(iim,jjmp1) :: zx_lon, zx_lat
188
189    CALL gather(rlat,rlat_glo)
190    CALL bcast(rlat_glo)
191    CALL gather(rlon,rlon_glo)
192    CALL bcast(rlon_glo)
193
194!$OMP MASTER
195    DO i=1,pim
196
197!    print*,'CFMIP_iophy i tabij lon lat',i,tabij(i),plon(i),plat(i)
198
199     plon_bounds(i,1)=rlon_glo(tabij(i)-1)
200     plon_bounds(i,2)=rlon_glo(tabij(i)+1)
201     if(plon_bounds(i,2).LE.0..AND.plon_bounds(i,1).GE.0.) THEN
202      if(rlon_glo(tabij(i)).GE.0.) THEN
203       plon_bounds(i,2)=-1*plon_bounds(i,2)
204      endif
205     endif
206     if(plon_bounds(i,2).GE.0..AND.plon_bounds(i,1).LE.0.) THEN
207      if(rlon_glo(tabij(i)).LE.0.) THEN
208       plon_bounds(i,2)=-1*plon_bounds(i,2)
209      endif
210     endif
211!
212     IF ( tabij(i).LE.iim) THEN
213      plat_bounds(i,1)=rlat_glo(tabij(i))
214     ELSE
215      plat_bounds(i,1)=rlat_glo(tabij(i)-iim)
216     ENDIF
217     plat_bounds(i,2)=rlat_glo(tabij(i)+iim)
218!
219!    print*,'CFMIP_iophy point i lon lon_bds',i,plon_bounds(i,1),rlon_glo(tabij(i)),plon_bounds(i,2)
220!    print*,'CFMIP_iophy point i lat lat_bds',i,plat_bounds(i,1),rlat_glo(tabij(i)),plat_bounds(i,2)
221!
222    ENDDO
223    if (is_sequential) then
224
225     npstn=pim
226     IF(.NOT. ALLOCATED(nptabij)) THEN
227      ALLOCATE(nptabij(pim))
228     ENDIF
229     DO i=1,pim
230      nptabij(i)=tabij(i)
231     ENDDO
232
233       CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlon_glo,zx_lon)
234       if ((iim*jjm).gt.1) then
235       DO i = 1, iim
236         zx_lon(i,1) = rlon_glo(i+1)
237         zx_lon(i,jjmp1) = rlon_glo(i+1)
238       ENDDO
239       endif
240       CALL gr_fi_ecrit(1,klon,iim,jjmp1,rlat_glo,zx_lat)
241
242    DO i=1,pim
243!    print*,'CFMIP_iophy i tabij lon lat',i,tabij(i),plon(i),plat(i)
244
245     plon_bounds(i,1)=zx_lon(ipt(i)-1,jpt(i))
246     plon_bounds(i,2)=zx_lon(ipt(i)+1,jpt(i))
247
248     if (ipt(i).EQ.1) then
249      plon_bounds(i,1)=zx_lon(iim,jpt(i))
250      plon_bounds(i,2)=360.+zx_lon(ipt(i)+1,jpt(i))
251     endif
252 
253     if (ipt(i).EQ.iim) then
254      plon_bounds(i,2)=360.+zx_lon(1,jpt(i))
255     endif
256
257     plat_bounds(i,1)=zx_lat(ipt(i),jpt(i)-1)
258     plat_bounds(i,2)=zx_lat(ipt(i),jpt(i)+1)
259
260     if (jpt(i).EQ.1) then
261      plat_bounds(i,1)=zx_lat(ipt(i),1)+0.001
262      plat_bounds(i,2)=zx_lat(ipt(i),1)-0.001
263     endif
264 
265     if (jpt(i).EQ.jjmp1) then
266      plat_bounds(i,1)=zx_lat(ipt(i),jjmp1)+0.001
267      plat_bounds(i,2)=zx_lat(ipt(i),jjmp1)-0.001
268     endif
269!
270!    print*,'CFMIP_iophy point i lon lon_bds',i,plon_bounds(i,1),rlon(tabij(i)),plon_bounds(i,2)
271!    print*,'CFMIP_iophy point i lat lat_bds',i,plat_bounds(i,1),rlat(tabij(i)),plat_bounds(i,2)
272!
273    ENDDO
274!    print*,'iophy is_sequential nname, nnhori, nnid_day=',trim(nname),nnhori,nnid_day
275     call histbeg(nname,pim,plon,plon_bounds, &
276                           plat,plat_bounds, &
277                           itau0, zjulian, dtime, nnhori, nnid_day)
278    else
279     npproc=0
280     DO ip=1, pim
281      tabprocbeg=klon_mpi_begin
282      tabprocend=klon_mpi_end
283      IF(tabij(ip).GE.tabprocbeg.AND.tabij(ip).LE.tabprocend) THEN
284       npproc=npproc+1
285       npstn=npproc
286      ENDIF
287     ENDDO
288!    print*,'CFMIP_iophy mpi_rank npstn',mpi_rank,npstn
289     IF(.NOT. ALLOCATED(nptabij)) THEN
290      ALLOCATE(nptabij(npstn))
291      ALLOCATE(npplon(npstn), npplat(npstn))
292      ALLOCATE(npplon_bounds(npstn,2), npplat_bounds(npstn,2))
293     ENDIF
294     npproc=0
295     DO ip=1, pim
296      IF(tabij(ip).GE.tabprocbeg.AND.tabij(ip).LE.tabprocend) THEN
297       npproc=npproc+1
298       nptabij(npproc)=tabij(ip)
299!      print*,'mpi_rank npproc ip plon plat tabij=',mpi_rank,npproc,ip, &
300!      plon(ip),plat(ip),tabij(ip)
301       npplon(npproc)=plon(ip)
302       npplat(npproc)=plat(ip)
303       npplon_bounds(npproc,1)=plon_bounds(ip,1)
304       npplon_bounds(npproc,2)=plon_bounds(ip,2)
305       npplat_bounds(npproc,1)=plat_bounds(ip,1)
306       npplat_bounds(npproc,2)=plat_bounds(ip,2)
307!!!
308!!! print qui sert a reordonner les points stations selon l'ordre CFMIP
309!!! ne pas enlever
310        print*,'iophy_mpi rank ip lon lat',mpi_rank,ip,plon(ip),plat(ip)
311!!!
312      ENDIF
313     ENDDO
314     call histbeg(nname,npstn,npplon,npplon_bounds, &
315                            npplat,npplat_bounds, &
316                            itau0,zjulian,dtime,nnhori,nnid_day,phys_domain_id)
317    endif
318!$OMP END MASTER
319
320  end subroutine histbeg_phy_points
321 
322  subroutine histwrite2d_phy(nid,lpoint,name,itau,field)
323  USE dimphy
324  USE mod_phys_lmdz_para
325  USE ioipsl
326  implicit none
327  include 'dimensions.h'
328   
329    integer,intent(in) :: nid
330    logical,intent(in) :: lpoint
331    character*(*), intent(IN) :: name
332    integer, intent(in) :: itau
333    real,dimension(:),intent(in) :: field
334    REAL,dimension(klon_mpi) :: buffer_omp
335    INTEGER, allocatable, dimension(:) :: index2d
336    REAL :: Field2d(iim,jj_nb)
337
338    integer :: ip
339    real,allocatable,dimension(:) :: fieldok
340
341    IF (size(field)/=klon) CALL abort_gcm('iophy::histwrite2d','Field first dimension not equal to klon',1)
342   
343    CALL Gather_omp(field,buffer_omp)   
344!$OMP MASTER
345    CALL grid1Dto2D_mpi(buffer_omp,Field2d)
346    if(.NOT.lpoint) THEN
347     ALLOCATE(index2d(iim*jj_nb))
348     ALLOCATE(fieldok(iim*jj_nb))
349     CALL histwrite(nid,name,itau,Field2d,iim*jj_nb,index2d)
350    else
351     ALLOCATE(fieldok(npstn))
352     ALLOCATE(index2d(npstn))
353
354     if(is_sequential) then
355!     klon_mpi_begin=1
356!     klon_mpi_end=klon
357      DO ip=1, npstn
358       fieldok(ip)=buffer_omp(nptabij(ip))
359      ENDDO
360     else
361      DO ip=1, npstn
362!     print*,'histwrite2d is_sequential npstn ip name nptabij',npstn,ip,name,nptabij(ip)
363       IF(nptabij(ip).GE.klon_mpi_begin.AND. &
364          nptabij(ip).LE.klon_mpi_end) THEN
365         fieldok(ip)=buffer_omp(nptabij(ip)-klon_mpi_begin+1)
366       ENDIF
367      ENDDO
368     endif
369     CALL histwrite(nid,name,itau,fieldok,npstn,index2d)
370!
371    endif
372    deallocate(index2d)
373    deallocate(fieldok)
374!$OMP END MASTER   
375  end subroutine histwrite2d_phy
376
377  subroutine histwrite3d_phy(nid,lpoint,name,itau,field)
378  USE dimphy
379  USE mod_phys_lmdz_para
380
381  use ioipsl
382  implicit none
383  include 'dimensions.h'
384   
385    integer,intent(in) :: nid
386    logical,intent(in) :: lpoint
387    character*(*), intent(IN) :: name
388    integer, intent(in) :: itau
389    real,dimension(:,:),intent(in) :: field  ! --> field(klon,:)
390    REAL,dimension(klon_mpi,size(field,2)) :: buffer_omp
391    REAL :: Field3d(iim,jj_nb,size(field,2))
392    INTEGER :: ip, n, nlev
393    INTEGER, ALLOCATABLE, dimension(:) :: index3d
394    real,allocatable, dimension(:,:) :: fieldok
395
396    IF (size(field,1)/=klon) CALL abort_gcm('iophy::histwrite3d','Field first dimension not equal to klon',1)
397    nlev=size(field,2)
398
399!   print*,'hist3d_phy mpi_rank npstn=',mpi_rank,npstn
400
401!   DO ip=1, npstn
402!    print*,'hist3d_phy mpi_rank nptabij',mpi_rank,nptabij(ip)
403!   ENDDO
404
405    CALL Gather_omp(field,buffer_omp)
406!$OMP MASTER
407    CALL grid1Dto2D_mpi(buffer_omp,field3d)
408    if(.NOT.lpoint) THEN
409     ALLOCATE(index3d(iim*jj_nb*nlev))
410     ALLOCATE(fieldok(iim*jj_nb,nlev))
411     CALL histwrite(nid,name,itau,Field3d,iim*jj_nb*nlev,index3d)
412    else
413      nlev=size(field,2)
414      ALLOCATE(index3d(npstn*nlev))
415      ALLOCATE(fieldok(npstn,nlev))
416
417      if(is_sequential) then
418!      klon_mpi_begin=1
419!      klon_mpi_end=klon
420       DO n=1, nlev
421       DO ip=1, npstn
422        fieldok(ip,n)=buffer_omp(nptabij(ip),n)
423       ENDDO
424       ENDDO
425      else
426       DO n=1, nlev
427       DO ip=1, npstn
428        IF(nptabij(ip).GE.klon_mpi_begin.AND. &
429         nptabij(ip).LE.klon_mpi_end) THEN
430         fieldok(ip,n)=buffer_omp(nptabij(ip)-klon_mpi_begin+1,n)
431        ENDIF
432       ENDDO
433       ENDDO
434      endif
435      CALL histwrite(nid,name,itau,fieldok,npstn*nlev,index3d)
436    endif
437  deallocate(index3d)
438  deallocate(fieldok)
439!$OMP END MASTER   
440  end subroutine histwrite3d_phy
441 
442end module iophy
Note: See TracBrowser for help on using the repository browser.