TABLE OF CONTENTS


ABINIT/m_haydock_io [ Modules ]

[ Top ] [ Modules ]

NAME

 m_haydock_io

FUNCTION

  This module provides routines to read the Haydock file

COPYRIGHT

  Copyright (C) 2013-2018 ABINIT group (YG)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

SOURCE

19 #if defined HAVE_CONFIG_H
20 #include "config.h"
21 #endif
22 
23 #include "abi_common.h"
24 
25 MODULE m_haydock_io
26 
27  use defs_basis
28  use m_abicore
29  use m_errors
30 
31  use m_io_tools,       only : open_file
32 
33  implicit none
34 
35  private

m_haydock_io/close_haydock [ Functions ]

[ Top ] [ m_haydock_io ] [ Functions ]

NAME

 close_haydock

FUNCTION

  Closes the haydock file and free the memory related to the file descriptor

INPUTS

  haydock_file = haydock file descriptor

PARENTS

      bsepostproc,m_haydock

CHILDREN

SOURCE

497 subroutine close_haydock(haydock_file)
498 
499 
500 !This section has been created automatically by the script Abilint (TD).
501 !Do not modify the following lines by hand.
502 #undef ABI_FUNC
503 #define ABI_FUNC 'close_haydock'
504 !End of the abilint section
505 
506  implicit none
507 
508 !Arguments ------------------------------------
509  type(haydock_type),intent(inout) :: haydock_file
510 
511 ! *************************************************************************
512 
513  close(haydock_file%unt)
514 
515  if(allocated(haydock_file%qpoints)) then
516    ABI_FREE(haydock_file%qpoints)
517  end if
518 
519  if(allocated(haydock_file%niter)) then
520    ABI_FREE(haydock_file%niter)
521  end if
522 
523 end subroutine close_haydock

m_haydock_io/haydock_type [ Types ]

[ Top ] [ m_haydock_io ] [ Types ]

NAME

  haydock_type

FUNCTION

  The structure defining the content of the haydock file

SOURCE

49  type,public :: haydock_type
50 
51    integer :: version
52    ! Version of the file
53 
54    integer :: hsize
55    ! Size of the hamiltonian
56 
57    integer :: use_coupling
58    ! 1 if coupling is used, 0 if not
59 
60    integer :: op
61    ! Type of the file
62 
63    integer :: nq
64    ! Number of q-points in the file
65 
66    integer :: unt
67    ! Fortran unit number.
68 
69    real(dp) :: broad
70    ! Broadening factor
71 
72    complex(dpc) :: factor
73 
74    real(dp),allocatable :: qpoints(:,:)
75    ! qpoints(3,nq)
76 
77    integer,allocatable :: niter(:)
78    ! niter(nq)
79 
80  end type haydock_type
81 
82  public :: open_haydock             ! Init the Haydock file
83  public :: read_haydock             ! Reads the data associated with 1 q-point
84  public :: read_dim_haydock         ! Read dimensions
85  public :: write_dim_haydock        ! Write dimensions.
86  public :: skip_dim_haydock         ! Skip the section with dimensions.
87  public :: write_haydock            ! Writes data related to a q-point inside the file
88  public :: close_haydock            ! Close the Haydock file and release memory

m_haydock_io/open_haydock [ Functions ]

[ Top ] [ m_haydock_io ] [ Functions ]

NAME

 open_haydock

FUNCTION

 Open the file and initialize haydock file descriptor that will
 be used for later access

INPUTS

  filename = Name of the file to be opened

OUTPUT

  haydock_file = file descriptor for the haydock file

PARENTS

      bsepostproc,m_haydock

CHILDREN

SOURCE

115 subroutine open_haydock(filename, haydock_file)
116 
117 
118 !This section has been created automatically by the script Abilint (TD).
119 !Do not modify the following lines by hand.
120 #undef ABI_FUNC
121 #define ABI_FUNC 'open_haydock'
122 !End of the abilint section
123 
124  implicit none
125 
126 !Arguments ------------------------------------
127 !scalars
128  character(len=*),intent(in) :: filename
129  type(haydock_type),intent(out) :: haydock_file
130 !arrays
131 
132 !Local variables ------------------------------
133 !scalars
134  character(len=500) :: msg
135 
136 !************************************************************************
137 
138  if (open_file(filename,msg,newunit=haydock_file%unt,form="unformatted") /= 0) then
139    MSG_ERROR(msg)
140  end if
141 
142  haydock_file%version = CUR_VERSION
143 
144 end subroutine open_haydock

m_haydock_io/read_dim_haydock [ Functions ]

[ Top ] [ m_haydock_io ] [ Functions ]

NAME

 read_dim_haydock

FUNCTION

  Reads the header dimensions of the haydock file and store them in
  the haydock file descriptor

 INPUT/OUTPUT
  haydock_file = haydock file descriptor

PARENTS

      bsepostproc,m_haydock

CHILDREN

SOURCE

167 subroutine read_dim_haydock(haydock_file)
168 
169 
170 !This section has been created automatically by the script Abilint (TD).
171 !Do not modify the following lines by hand.
172 #undef ABI_FUNC
173 #define ABI_FUNC 'read_dim_haydock'
174 !End of the abilint section
175 
176  implicit none
177 
178 !Arguments ------------------------------------
179 !scalars
180  type(haydock_type),intent(inout) :: haydock_file
181 
182 !Local variables ------------------------------
183 !scalars
184  integer :: niter_file
185  integer :: iq, it
186  character(len=500) :: msg
187 !arrays
188  real(dp) :: q_file(3)
189 
190 !************************************************************************
191 
192  read(haydock_file%unt) haydock_file%version
193 
194  msg = "The haydock file has been produced by an uncompatible version of abinit"
195  ABI_CHECK(haydock_file%version == CUR_VERSION, msg)
196 
197  read(haydock_file%unt) haydock_file%hsize,haydock_file%use_coupling,&
198 &   haydock_file%op,haydock_file%nq,haydock_file%broad
199 
200  ABI_MALLOC(haydock_file%qpoints,(3,haydock_file%nq))
201  ABI_MALLOC(haydock_file%niter,(haydock_file%nq))
202 
203  do iq = 1,haydock_file%nq
204    read(haydock_file%unt) q_file(:)
205    read(haydock_file%unt) niter_file
206 
207    haydock_file%qpoints(:,iq) = q_file(:)
208    haydock_file%niter(iq) = niter_file
209    ! Skip data for this q.
210    do it=1,niter_file
211      read(haydock_file%unt) ! it,aa(it),bb(it)
212    end do
213    read(haydock_file%unt) ! phi_nm1
214    read(haydock_file%unt) ! phi_n
215    read(haydock_file%unt) ! factor
216  end do
217 
218 end subroutine read_dim_haydock

m_haydock_io/read_haydock [ Functions ]

[ Top ] [ m_haydock_io ] [ Functions ]

NAME

 read_haydock

FUNCTION

  Reads the data related to one q-point inside the file accessed by file descriptor

INPUTS

  haydock_file = haydock file descriptor
  q = q-point to be searched inside the file

OUTPUT

  aa = coefficients "a" of the lanczos chain
  bb = coefficients "b" of the lanczos chain
  phi_n = last vector of the lanczos chain
  phi_nm1 = penultimate vector of the lanczos chain
  niter = number of iterations done
  factor = pre-factor used to obtain the green function

NOTES

  niter = 0 if the q-point has not been found

PARENTS

      bsepostproc,m_haydock

CHILDREN

SOURCE

338 subroutine read_haydock(haydock_file, q, aa, bb, phi_n, phi_nm1, niter, factor)
339 
340 
341 !This section has been created automatically by the script Abilint (TD).
342 !Do not modify the following lines by hand.
343 #undef ABI_FUNC
344 #define ABI_FUNC 'read_haydock'
345 !End of the abilint section
346 
347  implicit none
348 
349 !Arguments ------------------------------------
350 !scalars
351  integer,intent(out) :: niter
352  complex(dpc),intent(out) :: factor
353  type(haydock_type),intent(in) :: haydock_file
354 !arrays
355  real(dp),intent(in) :: q(3)
356  real(dp),allocatable,intent(out) :: bb(:)
357  complex(dpc),allocatable,intent(out) :: aa(:),phi_n(:),phi_nm1(:)
358 
359 !Local variables ------------------------------
360 !scalars
361  integer :: iq, it, inn
362  integer :: niter_file
363  logical :: found_q
364  character(len=500) :: msg
365 !arrays
366  real(dp) :: q_file(3)
367 
368 ! *************************************************************************
369 
370  rewind(haydock_file%unt)
371  call skip_dim_haydock(haydock_file)
372 
373  do iq = 1,haydock_file%nq
374    read(haydock_file%unt) q_file(:)
375    read(haydock_file%unt) niter_file
376 
377    if ( ALL(ABS(q_file - q) < tol6) ) then
378       found_q = .TRUE.; EXIT
379    else
380       ! Skip data for this q.
381       do it=1,niter_file
382         read(haydock_file%unt) ! it,aa(it),bb(it)
383       end do
384       read(haydock_file%unt) ! phi_nm1
385       read(haydock_file%unt) ! phi_n
386       read(haydock_file%unt) ! factor
387    end if
388  end do
389 
390  if(found_q) then
391    niter = niter_file
392    ABI_ALLOCATE(aa,(niter))
393    ABI_ALLOCATE(bb,(niter))
394    do inn=1,niter
395      read(haydock_file%unt)it,aa(inn),bb(inn)
396      if (inn/=it) then
397        write(msg,'(2(a,i0))')" Found it_file: ",it," while it should be: ",inn
398        MSG_ERROR(msg)
399      end if
400    end do
401    ABI_MALLOC(phi_nm1,(haydock_file%hsize))
402    ABI_MALLOC(phi_n,(haydock_file%hsize))
403    read(haydock_file%unt)phi_nm1
404    read(haydock_file%unt)phi_n
405    read(haydock_file%unt)factor
406  else
407    niter = 0
408  end if
409 
410 end subroutine read_haydock

m_haydock_io/skip_dim_haydock [ Functions ]

[ Top ] [ m_haydock_io ] [ Functions ]

NAME

 skip_dim_haydock

FUNCTION

  Skip the part of the file reading basic dimensions contained in the header

INPUTS

  haydock_file = haydock file descriptor

OUTPUT

PARENTS

      m_haydock_io

CHILDREN

SOURCE

284 subroutine skip_dim_haydock(haydock_file)
285 
286 
287 !This section has been created automatically by the script Abilint (TD).
288 !Do not modify the following lines by hand.
289 #undef ABI_FUNC
290 #define ABI_FUNC 'skip_dim_haydock'
291 !End of the abilint section
292 
293  implicit none
294 
295 !Arguments ------------------------------------
296 !scalars
297  type(haydock_type),intent(in) :: haydock_file
298 
299 ! *************************************************************************
300 
301  read(haydock_file%unt)
302  read(haydock_file%unt)
303 
304 end subroutine skip_dim_haydock

m_haydock_io/write_dim_haydock [ Functions ]

[ Top ] [ m_haydock_io ] [ Functions ]

NAME

 write_dim_haydock

FUNCTION

  Writes the basic dimensions stored inside the haydock descriptor inside the file

INPUTS

  haydock_file = haydock file descriptor

PARENTS

      m_haydock

CHILDREN

SOURCE

240 subroutine write_dim_haydock(haydock_file)
241 
242 
243 !This section has been created automatically by the script Abilint (TD).
244 !Do not modify the following lines by hand.
245 #undef ABI_FUNC
246 #define ABI_FUNC 'write_dim_haydock'
247 !End of the abilint section
248 
249  implicit none
250 
251 !Arguments ------------------------------------
252  type(haydock_type),intent(in) :: haydock_file
253 
254 ! *************************************************************************
255 
256  write(haydock_file%unt) haydock_file%version
257  write(haydock_file%unt) haydock_file%hsize,haydock_file%use_coupling, &
258 &  haydock_file%op,haydock_file%nq,haydock_file%broad
259 
260 end subroutine write_dim_haydock

m_haydock_io/write_haydock [ Functions ]

[ Top ] [ m_haydock_io ] [ Functions ]

NAME

 write_haydock

FUNCTION

  Writes data related to a q-point inside the file

INPUTS

  haydock_file = haydock file descriptor
  q = q-point to be searched inside the file
  aa = coefficients "a" of the lanczos chain
  bb = coefficients "b" of the lanczos chain
  phi_n = last vector of the lanczos chain
  phi_nm1 = penultimate vector of the lanczos chain
  niter = number of iterations done
  factor = pre-factor used to obtain the green function

PARENTS

      m_haydock

CHILDREN

SOURCE

439 subroutine write_haydock(haydock_file, hsize, q, aa, bb, phi_n, phi_nm1, niter, factor)
440 
441 
442 !This section has been created automatically by the script Abilint (TD).
443 !Do not modify the following lines by hand.
444 #undef ABI_FUNC
445 #define ABI_FUNC 'write_haydock'
446 !End of the abilint section
447 
448  implicit none
449 
450 !Arguments ------------------------------------
451 !scalars
452  integer,intent(in) :: niter,hsize
453  complex(dpc),intent(in) :: factor
454  type(haydock_type),intent(in) :: haydock_file
455 !arrays
456  real(dp),intent(in) :: q(3)
457  real(dp),intent(in) :: bb(niter)
458  complex(dpc),intent(in) :: aa(niter),phi_n(hsize),phi_nm1(hsize)
459 
460 !Local variables -----------------------------
461 !scalars
462  integer :: it
463 
464 ! *************************************************************************
465 
466  write(haydock_file%unt) q
467  write(haydock_file%unt) niter  ! NB if the previous loop completed inn=niter_max+1
468  do it=1,niter        ! if we exited then inn is not incremented by one.
469    write(haydock_file%unt)it,aa(it),bb(it)
470  end do
471  write(haydock_file%unt) phi_nm1
472  write(haydock_file%unt) phi_n
473  write(haydock_file%unt) factor
474 
475 end subroutine write_haydock