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! ---
8subroutine ts_show_regions(ucell,na_u,xa,N_Elec,Elecs)
9  use precision, only : dp
10  use units, only : Ang
11  use parallel, only : IONode
12  use fdf, only : fdf_get
13
14  use m_ts_electype
15  use m_ts_method, only : r_aBuf, r_aC, na_Buf, atom_type
16  use m_ts_method, only : TYP_BUFFER, TYP_DEVICE
17  use m_region
18
19  implicit none
20! ********************
21! * INPUT variables  *
22! ********************
23  real(dp), intent(in) :: ucell(3,3)
24  integer, intent(in)  :: na_u
25  real(dp), intent(in) :: xa(3,na_u)
26  integer, intent(in) :: N_Elec
27  type(Elec), intent(in) :: Elecs(N_Elec)
28
29! ********************
30! * LOCAL variables  *
31! ********************
32  integer :: ia, i, ia_mid
33  type(tRgn) :: rgn, rtmp
34
35  if ( .not. IONode ) return
36
37  if ( .not. fdf_get('TS.Atoms.Print',.false.) ) then
38
39    write(*,'(/,a)') 'transiesta: Regions of atoms:'
40
41    if ( r_aBuf%n > 0 ) then
42      call rgn_copy(r_aBuf, rgn)
43      call rgn_sort(rgn)
44      rgn%name = 'Buffer'
45      call rgn_print(rgn, name='//', seq_max=12, indent=3)
46    end if
47
48    ! Prepare the device region
49    call rgn_copy(r_aC, rtmp)
50    do i = 1 , N_Elec
51      ia = Elecs(i)%idx_a
52      ia_mid = ia + TotUsedAtoms(Elecs(i)) - 1
53      call rgn_range(rgn, ia, ia_mid)
54      rgn%name = 'Elec.'//trim(Elecs(i)%name)
55      call rgn_print(rgn, name='##', seq_max=12, indent=3)
56      ! Remove the electrode from the device
57      call rgn_remove(rtmp, rgn, rtmp)
58    end do
59    ! Sort the region
60    call rgn_sort(rtmp)
61    rtmp%name = 'Device'
62    call rgn_print(rtmp, name='--', seq_max=12, indent=3)
63
64    call rgn_delete(rgn, rtmp)
65
66    write(*,*) ! new-line
67
68    return
69
70  end if
71
72  ! mid-point of device
73  ia_mid = (na_u - na_Buf-sum(TotUsedAtoms(Elecs))+1) / 2
74
75  write(*,'(/,a)') 'transiesta: Atomic coordinates and regions (Ang):'
76  ia = 1
77  do while ( ia <= na_u )
78     select case ( atom_type(ia) )
79     case ( TYP_BUFFER )
80        i = ia
81        do while ( i <= na_u )
82           ! this assures that we do not have a memory leak
83           if ( atom_type(i) /= TYP_BUFFER ) exit
84           i = i + 1
85        end do
86        i = i - ia
87        ! steps ia counter
88        call out_REGION(ia,i,'Buffer','/')
89     case ( TYP_DEVICE )
90        ia_mid = ia_mid - 1
91        if ( ia_mid == 0 ) then
92           write(*,'(tr1,i5,tr2,3(tr2,f12.7),tr8,a)') ia,xa(:,ia)/Ang,'Device'
93        else
94           write(*,'(tr1,i5,tr2,3(tr2,f12.7))') ia,xa(:,ia)/Ang
95        end if
96        ia = ia + 1
97     case default
98        ! electrode position
99        do i = 1 , N_Elec
100           if ( ia == Elecs(i)%idx_a ) then
101              call out_REGION(ia,TotUsedAtoms(Elecs(i)), &
102                   trim(name(Elecs(i)))//' electrode','#')
103              exit
104           end if
105        end do
106     end select
107  end do
108
109  write(*,*)
110
111contains
112
113  subroutine out_REGION(ia,NA,name,marker,first,last)
114    integer, intent(inout) :: ia
115    integer, intent(in)    :: NA
116    character(len=*), intent(in) :: name
117    character(len=1), intent(in) :: marker
118    character(len=1), intent(in), optional :: first, last
119    integer :: i, mid
120    logical :: lfirst, llast
121    lfirst = present(first)
122    llast = present(last)
123
124    if ( NA < 1 ) return
125    mid = (NA+1) / 2
126    write(*,'(tr7,a)') repeat(marker,46)
127    do i = 1, NA
128       if ( i == 1 .and. lfirst ) then
129          write(*,'(tr1,i5,tr1,a1,3(tr2,f12.7),tr2,a1,tr1,a1)') &
130               ia,marker,xa(:,ia)/Ang,marker,first
131       else if ( i == mid ) then
132          write(*,'(tr1,i5,tr1,a1,3(tr2,f12.7),tr2,a1,tr5,a)') &
133               ia,marker,xa(:,ia)/Ang,marker,trim(name)
134       else if ( i == NA .and. llast ) then
135          write(*,'(tr1,i5,tr1,a1,3(tr2,f12.7),tr2,a1,tr1,a1)') &
136               ia,marker,xa(:,ia)/Ang,marker,last
137       else
138          write(*,'(tr1,i5,tr1,a1,3(tr2,f12.7),tr2,a1)') &
139               ia,marker,xa(:,ia)/Ang,marker
140       end if
141       ia = ia + 1
142    end do
143    write(*,'(tr7,a)') repeat(marker,46)
144
145  end subroutine out_REGION
146
147end subroutine ts_show_regions
148
149
150