1! ---
2! Copyright (C) 1996-2016	The SIESTA group
3!  This file is distributed under the terms of the
4!  GNU General Public License: see COPYING in the top directory
5!  or http://www.gnu.org/copyleft/gpl.txt .
6! See Docs/Contributors.txt for a list of contributors.
7! ---
8
9! This code has been fully implemented by:
10!  Nick Papior, 2015
11module m_ts_io_version
12
13  use precision, only : dp
14  implicit none
15
16contains
17
18  subroutine write_TSHS_0(filename, &
19       onlyS, Gamma, TSGamma, &
20       ucell, na_u, no_l, no_u, no_s, maxnh, nspin,  &
21       kscell, kdispl, &
22       xa, iza, lasto, &
23       numh, listhptr, listh, xij, &
24       H, S, Ef, &
25       Qtot, Temp, &
26       istep, ia1)
27
28    use geom_helper, only : ucorb
29
30! **********************
31! * INPUT variables    *
32! **********************
33    character(len=*), intent(in) :: filename
34    logical, intent(in) :: onlyS
35    logical, intent(in) :: Gamma, TSGamma
36    real(dp), intent(in) :: ucell(3,3)
37    integer, intent(in) :: na_u, no_l, no_u, no_s, maxnh, nspin
38    real(dp), intent(in) :: xa(3,na_u)
39    integer, intent(in) :: iza(na_u)
40    integer, intent(in) :: numh(no_l), listhptr(no_l)
41    integer, intent(in) :: listh(maxnh)
42    real(dp), intent(in) :: xij(3,maxnh)
43    integer, intent(in) :: lasto(0:na_u)
44    real(dp), intent(in) :: H(maxnh,nspin), S(maxnh)
45    real(dp), intent(in) :: Ef
46    integer, intent(in) :: kscell(3,3)
47    real(dp), intent(in) :: kdispl(3)
48    real(dp), intent(in) :: Qtot,Temp
49    integer, intent(in) :: istep, ia1
50
51! ************************
52! * LOCAL variables      *
53! ************************
54    integer :: indxuo(no_s)
55    integer :: iu
56    integer :: ispin, i,j
57
58    external :: io_assign, io_close
59
60    ! Open file
61    call io_assign( iu )
62    open( iu, file=filename, form='unformatted', status='unknown' )
63
64    ! Write Dimensions Information
65    write(iu) na_u, no_u, no_s, nspin, maxnh
66
67    ! Write Geometry information
68    write(iu) xa
69    write(iu) iza
70    write(iu) ucell
71
72    ! Write k-point samplung information
73    write(iu) Gamma
74    ! Is this only an S containing object?
75    write(iu) onlyS
76    write(iu) TSGamma
77    write(iu) kscell
78    write(iu) kdispl
79    write(iu) istep, ia1
80
81    write(iu) lasto
82
83    if ( .not. Gamma ) then
84       do i = 1 , no_s
85          indxuo(i) = ucorb(i,no_u)
86       end do
87       write(iu) indxuo
88    endif
89
90    write(iu) numh
91
92    ! Write Electronic Structure Information
93    write(iu) Qtot,Temp
94    write(iu) Ef
95
96    ! Write listh
97    do i = 1 , no_u
98       write(iu) listh(listhptr(i)+1:listhptr(i)+numh(i))
99    end do
100
101    ! Write Overlap matrix
102    do i = 1 , no_u
103       write(iu) S(listhptr(i)+1:listhptr(i)+numh(i))
104    end do
105
106    if ( .not. onlyS ) then
107       ! Write Hamiltonian
108       do ispin = 1 , nspin
109          do i = 1 , no_u
110             write(iu) H(listhptr(i)+1:listhptr(i)+numh(i),ispin)
111          end do
112       end do
113    end if
114
115    if ( .not. Gamma ) then
116       do i = 1 , no_u
117          write(iu) (xij(j,listhptr(i)+1:listhptr(i)+numh(i)),j=1,3)
118       end do
119    end if
120
121    ! Close file
122    call io_close( iu )
123
124  end subroutine write_TSHS_0
125
126end module m_ts_io_version
127
128
129