TABLE OF CONTENTS


ABINIT/m_io_redirect [ Modules ]

[ Top ] [ Modules ]

NAME

  m_io_redirect

FUNCTION

  management of output and log files when parallelisation on cells is activated

COPYRIGHT

  Copyright (C) 2001-2018 ABINIT group (FJ,MT)
  This file is distributed under the terms of the
  GNU General Public License, see ~ABINIT/Infos/copyright
  or http://www.gnu.org/copyleft/gpl.txt .

SOURCE

17 #if defined HAVE_CONFIG_H
18 #include "config.h"
19 #endif
20 
21 #include "abi_common.h"
22 
23 module m_io_redirect
24 
25  use defs_basis
26  use m_abicore
27  use m_errors
28 
29  use m_xmpi,        only : xmpi_comm_rank, xmpi_barrier
30  use m_libpaw_tools,only : libpaw_write_comm_set
31  use m_io_tools,    only : delete_file, flush_unit, open_file, get_unit
32  use m_fstrings,    only : int2char4
33 
34  implicit none
35 
36  private
37 
38  character(len=fnlen),allocatable,save :: fillog(:),filout(:)
39 
40 !Procedures ------------------------------------
41  public :: localfilnam    ! Define new appended name
42  public :: localwrfile    ! write local files
43  public :: localrdfile    ! read local files
44  public :: localredirect  ! redirect communicators for local files
45 
46 contains

m_io_redirect/localfilnam [ Functions ]

[ Top ] [ m_io_redirect ] [ Functions ]

NAME

 localfilnam

FUNCTION

 Define new appended name

INPUTS

   commspace= communicator between cells
   commworld= communicator of the world
   filnam(5)=character strings giving file names
   nam= name of the extension of the files
   nfil= number of files to be named

PARENTS

      dfpt_looppert,gstateimg

CHILDREN

      abi_io_redirect,libpaw_write_comm_set

SOURCE

 71  subroutine localfilnam(commspace,commspace1,commworld,filnam,nam,nfil)
 72 
 73 
 74 !This section has been created automatically by the script Abilint (TD).
 75 !Do not modify the following lines by hand.
 76 #undef ABI_FUNC
 77 #define ABI_FUNC 'localfilnam'
 78 !End of the abilint section
 79 
 80  implicit none
 81  
 82 !Arguments ------------------------------------
 83 !scalars
 84  integer,intent(in) :: commspace,commspace1,commworld,nfil
 85  character(len=4) :: nam
 86 !arrays
 87  character(len=fnlen),intent(in) :: filnam(:)
 88 
 89 !Local variables-------------------------------
 90 !scalars
 91  integer :: ii,ios,me,me_loc
 92  character(len=10) :: appen,tag
 93 
 94 ! *********************************************************************
 95  
 96  if (nfil<=1) return
 97 
 98  ABI_ALLOCATE(filout,(nfil))
 99  ABI_ALLOCATE(fillog,(nfil))
100  me=xmpi_comm_rank(commworld)
101  me_loc=xmpi_comm_rank(commspace1)
102  call int2char4(me,tag)
103  ABI_CHECK((tag(1:1)/='#'),'Bug: string length too short!')
104  do ii=1,nfil
105    call appdig(ii,'',appen)
106    filout(ii)=trim(filnam(2))//nam//trim(appen)
107    if (me_loc==0) then
108      fillog(ii)=trim(filnam(5))//nam//trim(appen)//"_LOG"
109    else
110      fillog(ii)=trim(filnam(5))//nam//trim(appen)//"_LOG_"//trim(tag)
111    end if
112  end do
113 
114  if (me==0) then
115    do ii=1,nfil
116      call delete_file(filout(ii),ios)
117    end do
118  end if
119 
120  call xmpi_barrier(commspace)
121 
122  end subroutine localfilnam

m_io_redirect/localrdfile [ Functions ]

[ Top ] [ m_io_redirect ] [ Functions ]

NAME

   localrdfile

FUNCTION

   read local (for each cell) output and log files
   to gather in one output or log file

INPUTS

   commspace= communicator between cells
   commworld= communicator of the world
   compute_all= flag to activate the reading on all cells
   [dyn(nfil)]= --optional-- indicate if this cell is to be taken into account
   nfil= number of files to be named
   paral= flag to activate parallelisation
   prtvol= flag to activate printing

PARENTS

      dfpt_looppert,gstateimg

CHILDREN

      abi_io_redirect,libpaw_write_comm_set

SOURCE

228  subroutine localrdfile(commspace,commworld,compute_all,nfil,paral,prtvol,dyn)
229 
230 
231 !This section has been created automatically by the script Abilint (TD).
232 !Do not modify the following lines by hand.
233 #undef ABI_FUNC
234 #define ABI_FUNC 'localrdfile'
235 !End of the abilint section
236 
237  implicit none
238  
239 !Arguments ------------------------------------
240 !scalars
241  integer,intent(in) :: commspace,commworld,nfil,paral,prtvol
242  logical,intent(in) :: compute_all
243 !arrays
244  integer,intent(in),target,optional :: dyn(nfil)
245 
246 !Local variables ------------------------------
247 !scalars
248  integer :: ii,ios,me,temp_unit
249  logical :: eof
250  character(len=500) :: msg
251  character(len=fnlen) :: line
252 !arrays
253  integer,pointer :: my_dyn(:)
254 
255 ! *********************************************************************
256 
257  if (nfil<=1) return
258  
259  me=xmpi_comm_rank(commworld)
260  call xmpi_barrier(commspace) !waiting for all files to be written and close
261 
262  if (prtvol<=0.and.(paral==1)) then
263    if (me==0) then
264      if (present(dyn)) then
265        my_dyn => dyn
266      else
267        ABI_ALLOCATE(my_dyn,(nfil))
268        my_dyn(:)=1
269      end if
270      do ii=1,nfil
271        if (open_file(filout(ii),msg,newunit=temp_unit,status='old') == 0) then
272          if (my_dyn(ii)==1.or.compute_all) then
273            eof=.false.
274            do while (.not.eof)
275              read(temp_unit,fmt='(a)',err=111,end=111,iostat=ios) line
276              if (ios==0) write(ab_out,'(a)') trim(line)
277              goto 112
278              111              eof=.true.
279              112              continue
280            end do
281          end if
282          close(unit=temp_unit,status='delete')
283        end if
284      end do
285      call flush_unit(ab_out)
286      if (.not.present(dyn)) then
287        ABI_DEALLOCATE(my_dyn)
288      end if
289    end if
290    call xmpi_barrier(commspace)
291  end if
292 
293  if (paral==1.and.do_write_log) then
294    if (me==0) then
295      do ii=1,nfil
296        if (open_file(fillog(ii),msg,newunit=temp_unit,status='old') == 0) then
297          eof=.false.
298          do while (.not.eof)
299            read(temp_unit,fmt='(a)',err=113,end=113,iostat=ios) line
300            if (ios==0) write(std_out,'(a)') trim(line)
301            goto 114
302            113            eof=.true.
303            114            continue
304          end do
305          close(unit=temp_unit,status='delete')
306        end if
307      end do
308    end if
309  end if
310 
311  call xmpi_barrier(commspace)
312 
313  ABI_DEALLOCATE(filout)
314  ABI_DEALLOCATE(fillog)
315 
316  end subroutine localrdfile

m_io_redirect/localredirect [ Functions ]

[ Top ] [ m_io_redirect ] [ Functions ]

NAME

   localredirect

FUNCTION

   Close output units ; restore defaults

INPUTS

   commspace= communicator between cells
   commworld= communicator of the world
   nfil= number of files to be named
   paral= flag to activate parallelisation
   prtvol= flag to activate printing  

PARENTS

      dfpt_looppert,gstateimg

CHILDREN

      abi_io_redirect,libpaw_write_comm_set

SOURCE

341  subroutine localredirect(commspace,commworld,nfil,paral,prtvol)
342 
343 
344 !This section has been created automatically by the script Abilint (TD).
345 !Do not modify the following lines by hand.
346 #undef ABI_FUNC
347 #define ABI_FUNC 'localredirect'
348 !End of the abilint section
349 
350  implicit none
351  
352 !Arguments ------------------------------------
353 !scalars
354  integer, intent(in) :: commspace,commworld,nfil,paral,prtvol
355 
356 !Local variables ------------------------------
357  integer :: me
358 
359 ! *********************************************************************
360 
361  if (nfil<=1) return
362 
363  me=xmpi_comm_rank(commspace)
364  if (prtvol>0) then
365    if (do_write_log.and.me==0) close(unit=ab_out)
366  else if (paral==1) then
367    if (me==0) close(unit=ab_out)
368  end if
369  if (paral==1.and.me==0.and.do_write_log) close(unit=std_out)
370  call abi_io_redirect(new_ab_out=ab_out_default,new_std_out=std_out_default,new_io_comm=commworld)
371  call libpaw_write_comm_set(commworld)
372 
373  end subroutine localredirect

m_io_redirect/localwrfile [ Functions ]

[ Top ] [ m_io_redirect ] [ Functions ]

NAME

   localwrfile

FUNCTION

   Write local (for each cell) output and log files

INPUTS

   commspace= communicator inside cells
   ii=number of the iteration on the cells
   nfil= number of files to be named
   paral= flag to activate parallelisation
   prtvol= flag to activate printing

PARENTS

      dfpt_looppert,gstateimg

CHILDREN

      abi_io_redirect,libpaw_write_comm_set

SOURCE

147  subroutine localwrfile(commspace,ii,nfil,paral,prtvol)
148 
149 
150 !This section has been created automatically by the script Abilint (TD).
151 !Do not modify the following lines by hand.
152 #undef ABI_FUNC
153 #define ABI_FUNC 'localwrfile'
154 !End of the abilint section
155 
156  implicit none
157 
158 !Arguments ------------------------------------
159  integer, intent(in) :: commspace,ii,nfil,paral,prtvol
160 
161 !Local variables ------------------------------
162  integer :: me
163  character(len=500) :: msg
164 
165 ! *********************************************************************
166 
167  if (nfil<=1) return
168 
169  me=xmpi_comm_rank(commspace)
170  if (paral==1) then
171    call abi_io_redirect(new_io_comm=commspace)
172    call libpaw_write_comm_set(commspace)
173  end if
174  if (prtvol>0) then
175    if (do_write_log) then
176      call abi_io_redirect(new_ab_out=get_unit())
177      if (me==0) then 
178        if (open_file(NULL_FILE,msg,unit=ab_out,status='unknown') /= 0) then
179          MSG_ERROR(msg)
180        end if
181      end if
182    else
183      call abi_io_redirect(new_ab_out=std_out)
184    end if
185  else if (paral==1) then
186    call abi_io_redirect(new_ab_out=get_unit())
187    if (me==0) then 
188      if (open_file(filout(ii),msg,unit=ab_out,status='unknown') /= 0) then
189        MSG_ERROR(msg)
190      end if
191    end if
192  end if
193  if (paral==1.and.me==0.and.do_write_log) then
194    call abi_io_redirect(new_std_out=get_unit())
195    if (open_file(fillog(ii),msg,unit=std_out,status='unknown') /= 0) then
196      MSG_ERROR(msg)
197    end if
198  end if
199 
200  end subroutine localwrfile