TABLE OF CONTENTS


ABINIT/work_var_lotf [ Modules ]

[ Top ] [ Modules ]

NAME

 work_var_lotf

FUNCTION

COPYRIGHT

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

SOURCE

15 #if defined HAVE_CONFIG_H
16 #include "config.h"
17 #endif
18 
19 #include "abi_common.h"
20 
21 module work_var_lotf
22 
23  use defs_basis
24  use defs_param_lotf
25  use m_errors
26  use m_abicore
27 
28  implicit none
29 
30  public ::             &
31    work_var_set,       &
32    work_var_dealloc,   &
33    cutoff_init,        &
34    smallfit
35 
36  public
37 
38  !--Control variables
39  ! !--Atomflags variables
40  integer,allocatable ::  ifixed(:) !--MMANCINI what is its utility
41 
42  !--Quantflags variables
43  integer,allocatable :: iq(:)
44 
45  !--Cutoff variables
46  real(dp) :: rcut,rcut_nbl
47  real(dp) :: rcrust
48 
49 contains

work_var_lotf/cutoff_init [ Functions ]

[ Top ] [ work_var_lotf ] [ Functions ]

NAME

  cutoff_init

FUNCTION

INPUTS

PARENTS

      m_lotf

CHILDREN

      dist_pbc

SOURCE

150  subroutine cutoff_init()
151   use pbc_lotf,only : pbc_bb_contract
152   !Local ---------------------------
153 
154 !This section has been created automatically by the script Abilint (TD).
155 !Do not modify the following lines by hand.
156 #undef ABI_FUNC
157 #define ABI_FUNC 'cutoff_init'
158 !End of the abilint section
159 
160   real(dp)  :: bl(3),blmin
161   character(len=500) :: message
162 
163 ! *************************************************************************
164 
165    if (lotfvar%classic==5 .OR. lotfvar%classic==6) then
166      rcut = 4.0d0 / 0.529177 !--check consistency with glue parameters (mind d and its "limits"!)
167      rcut_nbl = rcut + rcrust
168    end if
169 
170    if (lotfvar%me==1.AND.lotfvar%classic==5) then
171      write(message,'(2(a,f12.6,a))')&
172 &     'GLUE radial cutoff used: ',rcut,ch10,&
173 &     'GLUE NEIGBOURS cutoff used: ',rcut_nbl,ch10
174      call wrtout(std_out,message,'COLL')
175    end if
176 
177   !--Cut-off check respect to cell size :
178    bl(:) = pbc_bb_contract()
179 
180    if (lotfvar%me==1) then
181      write(message,'(3a,3f12.6,a)')&
182 &     'LENGTH OF REC. CELL VECTORS : ',ch10,&
183 &     ' bl1, bl2, bl3 : ',   bl(:),ch10
184      call wrtout(std_out,message,'COLL')
185    end if
186    blmin = minval(bl)
187 
188    if (rcut_nbl*blmin  > half) then
189      write(message,'(2a,2(a,f12.6))')&
190 &     'LOTF: cut off too large : ',ch10,&
191 &     ' cut-off (A) is ', rcut_nbl ,  ' min. allowed : ',half/blmin
192      MSG_ERROR(message)
193    end if
194  end subroutine cutoff_init

work_var_lotf/smallfit [ Functions ]

[ Top ] [ work_var_lotf ] [ Functions ]

NAME

 smallfit

FUNCTION

INPUTS

PARENTS

      m_lotf

CHILDREN

      dist_pbc

SOURCE

213  subroutine smallfit(tau0,ndum)
214   use bond_lotf,only : tafit
215   USE pbc_lotf,only : dist_pbc,r2
216 
217 !This section has been created automatically by the script Abilint (TD).
218 !Do not modify the following lines by hand.
219 #undef ABI_FUNC
220 #define ABI_FUNC 'smallfit'
221 !End of the abilint section
222 
223   implicit none
224 
225   !Arguments ------------------------
226   integer,intent(out) :: ndum
227   real(dp),intent(in):: tau0(3,lotfvar%natom)
228 
229   !Local ---------------------------
230   integer  ::  i, j, jat, nquant
231   real(dp) :: rag_fit
232 
233 ! *************************************************************************
234 
235    rag_fit = zero
236 
237    nquant = lotfvar%natom !--to remember old notation with nquant
238 
239    ndum = 0
240    do i = 1, lotfvar%natom
241      if(.not.tafit(i)) then
242        do j=1,nquant
243          jat = iq(j)
244          call dist_pbc(tau0(:,i),tau0(:,jat))
245          if (r2  <  rag_fit) then
246            tafit(i) = .true.
247            ndum = ndum + 1
248            cycle
249          end if
250        end do !nqtot
251      end if !tafit
252    end do
253 
254  end subroutine smallfit
255 
256 end module work_var_lotf

work_var_lotf/work_var_dealloc [ Functions ]

[ Top ] [ work_var_lotf ] [ Functions ]

NAME

 work_var_dealloc

FUNCTION

  deallocate variable

INPUTS

PARENTS

      m_lotf

CHILDREN

      dist_pbc

SOURCE

119  subroutine work_var_dealloc()
120 
121 ! *************************************************************************
122 
123 !This section has been created automatically by the script Abilint (TD).
124 !Do not modify the following lines by hand.
125 #undef ABI_FUNC
126 #define ABI_FUNC 'work_var_dealloc'
127 !End of the abilint section
128 
129    ABI_DEALLOCATE(iq)
130    ABI_DEALLOCATE(ifixed)
131  end subroutine work_var_dealloc

work_var_lotf/work_var_set [ Functions ]

[ Top ] [ work_var_lotf ] [ Functions ]

NAME

 work_var_set

FUNCTION

  set some internal variable of lotf

INPUTS

  natom=number of atoms

PARENTS

      m_lotf

CHILDREN

      dist_pbc

SOURCE

70  subroutine work_var_set()
71 
72 
73 !This section has been created automatically by the script Abilint (TD).
74 !Do not modify the following lines by hand.
75 #undef ABI_FUNC
76 #define ABI_FUNC 'work_var_set'
77 !End of the abilint section
78 
79   implicit none
80   !Local-----------------------------
81   integer :: iat
82 
83 ! *************************************************************************
84 
85   !--ifixed from ATOMFLAGS is initialized :
86    ABI_ALLOCATE(ifixed,(lotfvar%natom))
87    ifixed(:) = 1
88 
89   !--FINDS  FITTED ATOMS
90   ! ABI_ALLOCATE(tquant,(lotfvar%natom))
91   ! tquant(:) = .true.
92   !  nquant = lotfvar%natom
93   !nqxx   = lotfvar%natom
94 
95    ABI_ALLOCATE(iq,(lotfvar%natom))
96    iq(:)=(/(iat,iat=1,lotfvar%natom)/)
97 
98  end subroutine work_var_set