Ignore:
Timestamp:
Jun 4, 2007, 4:34:47 PM (17 years ago)
Author:
Laurent Fairhead
Message:

Merge entre la version V3_conv et le HEAD
YM, JG, LF

File:
1 edited

Legend:

Unmodified
Added
Removed
  • LMDZ4/trunk/libf/phylmd/iophy.F90

    r629 r766  
    77  REAL,private,allocatable,dimension(:),save :: io_lat
    88  REAL,private,allocatable,dimension(:),save :: io_lon
     9  INTEGER, save :: phys_domain_id
    910 
    1011  INTERFACE histwrite_phy
    1112    MODULE PROCEDURE histwrite2d_phy,histwrite3d_phy
    1213  END INTERFACE
     14 
     15  REAL,private,allocatable,save,dimension(:,:) :: buffer_omp
    1316
    1417contains
     
    1619  subroutine init_iophy(lat,lon)
    1720  use dimphy
     21  use ioipsl
    1822  implicit none
    1923  include 'dimensions90.h'   
    2024    real,dimension(iim),intent(in) :: lon
    2125    real,dimension(jjm+1),intent(in) :: lat
    22  
     26
     27    INTEGER,DIMENSION(2) :: ddid
     28    INTEGER,DIMENSION(2) :: dsg
     29    INTEGER,DIMENSION(2) :: dsl
     30    INTEGER,DIMENSION(2) :: dpf
     31    INTEGER,DIMENSION(2) :: dpl
     32    INTEGER,DIMENSION(2) :: dhs
     33    INTEGER,DIMENSION(2) :: dhe
     34
     35!$OMP MASTER 
    2336    allocate(io_lat(jjm+1))
    2437    io_lat(:)=lat(:)
     
    3144    ndex2d(:)=0
    3245    ndex3d(:)=0
     46    allocate(buffer_omp(klon_mpi,klev))
     47   
     48    ddid=(/ 1,2 /)
     49    dsg=(/ iim, jjm+1 /)
     50    dsl=(/ iim, jjphy_nb /)
     51    dpf=(/ 1,jjphy_begin /)
     52    dpl=(/ iim, jjphy_end /)
     53    dhs=(/ iiphy_begin-1,0 /)
     54    if (phy_rank==phy_size-1) then
     55      dhe=(/0,0/)
     56    else
     57      dhe=(/ iim-iiphy_end,0 /) 
     58    endif
     59   
     60    call flio_dom_set(phy_size,phy_rank,ddid,dsg,dsl,dpf,dpl,dhs,dhe, &
     61                      'APPLE',phys_domain_id)
    3362
     63!$OMP END MASTER
     64!$OMP FLUSH(buffer_omp)
     65     
    3466  end subroutine init_iophy
    3567 
     
    4779    integer,intent(out) :: nhori
    4880    integer,intent(out) :: nid_day
    49    
     81
     82!$OMP MASTER   
    5083    if (monocpu) then
    5184      call histbeg(name,iim,io_lon, jjphy_nb,io_lat(jjphy_begin:jjphy_end), &
    5285                   1,iim,1,jjphy_nb,itau0, zjulian, dtime, nhori, nid_day)
    5386    else
    54       call histbeg(name//'_'//trim(int2str(phy_rank)),iim,io_lon, jjphy_nb,io_lat(jjphy_begin:jjphy_end), &
    55                    1,iim,1,jjphy_nb,itau0, zjulian, dtime, nhori, nid_day)
     87      call histbeg(name,iim,io_lon, jjphy_nb,io_lat(jjphy_begin:jjphy_end), &
     88                   1,iim,1,jjphy_nb,itau0, zjulian, dtime, nhori, nid_day,phys_domain_id)
    5689    endif
     90!$OMP END MASTER
    5791 
    5892  end subroutine histbeg_phy
     
    68102    integer, intent(in) :: itau
    69103    real,dimension(klon),intent(in) :: field
    70    
    71     CALL phy2dyn(field,tmp_tab2d,1)
     104
     105    CALL GatherField_omp(field,buffer_omp,1)   
     106!$OMP MASTER
     107    CALL phy2dyn(buffer_omp,tmp_tab2d,1)
    72108    CALL histwrite(nid,name,itau,tmp_tab2d,iim*jjphy_nb,ndex2d)
    73    
     109!$OMP END MASTER   
    74110  end subroutine histwrite2d_phy
    75111 
     
    85121    real,dimension(klon,klev),intent(in) :: field
    86122   
    87     CALL phy2dyn(field,tmp_tab3d,klev)
     123    CALL GatherField_omp(field,buffer_omp,klev)
     124!$OMP MASTER
     125    CALL phy2dyn(buffer_omp,tmp_tab3d,klev)
    88126    CALL histwrite(nid,name,itau,tmp_tab3d,iim*jjphy_nb*klev,ndex3d)
    89    
     127!$OMP END MASTER   
    90128  end subroutine histwrite3d_phy
    91129 
     
    96134  include 'dimensions90.h'
    97135 
    98     real,dimension(klon,nlev),intent(in) :: field_phy
     136    real,dimension(klon_mpi,nlev),intent(in) :: field_phy
    99137    real,dimension(iim,jjphy_nb,nlev),intent(out) :: field_dyn
    100138    integer,intent(in) :: nlev
     
    107145        if (jjphy_begin==jjphy_end) then
    108146          field_dyn(:,1,l)=0.
    109           field_dyn(iiphy_begin:iiphy_end,1,l)=field_phy(1:klon,l)
     147          field_dyn(iiphy_begin:iiphy_end,1,l)=field_phy(1:klon_mpi,l)
    110148        else
    111149       
     
    115153         else
    116154           field_dyn(:,1,l)=0.
    117            next=iim-iiphy_begin+1
     155           next=iim-iiphy_begin+2
    118156           field_dyn(iiphy_begin:iim,1,l)=field_phy(1:next-1,l)   
    119157         endif
     
    125163         
    126164          if (jjphy_end==jjm+1) then
    127              field_dyn(:,jjphy_nb,l)=field_phy(klon,l)
     165             field_dyn(:,jjphy_nb,l)=field_phy(klon_mpi,l)
    128166          else
    129167           field_dyn(:,jjphy_nb,l)=0.
Note: See TracChangeset for help on using the changeset viewer.