TABLE OF CONTENTS


ABINIT/bond_lotf [ Modules ]

[ Top ] [ Modules ]

NAME

 bond_lotf

FUNCTION

  Define BOND variables and the procedure to
  set them.

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

18 #if defined HAVE_CONFIG_H
19 #include "config.h"
20 #endif
21 
22 #include "abi_common.h"
23 
24 module bond_lotf
25 
26  use defs_basis
27  use m_errors
28  use m_abicore
29 
30  implicit none
31 
32  public
33 
34  !--Fittedatoms variables
35  integer :: nfit !--dimension needed for fits
36  integer :: nfitmax  !--dimension of the fit arrays
37  integer,allocatable :: ifit(:)
38  logical,allocatable :: tafit(:)
39 
40  !--Fittedbonds variables
41  integer :: nbondex
42  integer :: ibn_tot,ibn_tot2,ibn_tots
43  integer,allocatable :: imat(:)
44  integer,dimension(:,:),allocatable :: ibnd_mat,ibmat,ibmat_large
45 
46  public ::             &
47    bond_atom_init,     &
48    bond_tafit_init,    &
49    bond_matrix_alloc,  &
50    bond_matrix_set,    &
51    bond_compute,       &
52    bond_fit_set,       &
53    bond_dealloc
54 
55 
56 contains

bond_lotf/bond_atom_init [ Functions ]

[ Top ] [ bond_lotf ] [ Functions ]

NAME

 bond_atom_init

FUNCTION

INPUTS

PARENTS

      m_lotf

CHILDREN

SOURCE

113  subroutine bond_atom_init(nneigx,nneig,neighl)
114 
115 
116 !This section has been created automatically by the script Abilint (TD).
117 !Do not modify the following lines by hand.
118 #undef ABI_FUNC
119 #define ABI_FUNC 'bond_atom_init'
120 !End of the abilint section
121 
122   implicit none
123 
124   !Arguments ------------------------
125   integer,intent(in) :: nneigx
126   integer,intent(in) :: nneig(:)
127   integer,intent(in) :: neighl(:,:)
128 
129   ! nneigx : max number of neighbours
130   ! tau0(3,natom) : atomic positions
131   ! neighl(nneigx,natom) : list of neighbours
132   ! nneig(natom) : number of neighbours
133   ! niter  : iteration number (itime)
134 
135   !Local ---------------------------
136   integer :: i,j,iat,jat, ibn
137   integer,allocatable :: ibnd_dum(:,:)
138   character(len=500) :: msg
139 
140 ! *************************************************************************
141 
142   !--Initialize nbondex
143    nbondex = ((nfitmax/2 * nneigx)+1)/2
144 
145   !--Now find initial numbers of active/border bonds :
146   !  (0)        CLEARS THE BOND(atom) MATRIX :
147    ABI_ALLOCATE(ibnd_mat,(2,nbondex))
148    ABI_ALLOCATE(ibnd_dum,(2,nfitmax*6))
149    ibnd_mat = 0
150    ibnd_dum = 0
151 
152    ibn_tot  = 0    !-- bonds between fitted atoms
153    ibn_tot2 = 0    !--existing but non optimized bonds with border atoms
154    ibn_tots = 0    !--total number of bonds in the fit + border zone
155 
156    do i =1,nfit
157      iat = ifit(i)
158      do j = 1,nneig(iat)
159        jat = neighl(j,iat)
160        if(tafit(jat)) then
161          if(jat > iat) then !--jat is a fitted atom
162            ibn_tot = ibn_tot + 1
163            ibnd_mat(:,ibn_tot) = (/iat,jat/)
164          end if
165        else    !--jat is a border atom
166          ibn_tot2 = ibn_tot2 + 1
167          ibnd_dum(:,ibn_tot2) = (/iat,jat/)
168        end if
169      end do
170    end do
171 
172    if(ibn_tot2 > (6*nfitmax)) then
173      write(msg,'(3a,i8,2a)')&
174 &     'ERROR: BOND_ATOM_INIT',ch10,&
175 &     'IBN_TOT2 =  ',ibn_tot2,ch10,&
176 &     ' ibnd_dum out of bounds, ibn_tot2 too large '
177      MSG_ERROR(msg)
178 
179    end if
180 
181   !--Reorder to keep 'variational' bonds first :
182    ibn_tots = ibn_tot + ibn_tot2
183    do ibn=ibn_tot+1,ibn_tots
184      ibnd_mat(:,ibn) = ibnd_dum(:,ibn-ibn_tot)
185    end do
186 
187 
188    ABI_DEALLOCATE(ibnd_dum)
189  end subroutine bond_atom_init

bond_lotf/bond_compute [ Functions ]

[ Top ] [ bond_lotf ] [ Functions ]

NAME

 bond_compute

FUNCTION

  Updates bond matrix, associates the bond to the atom neighlists

INPUTS

  nneig
  neighl

PARENTS

      m_lotf

CHILDREN

SOURCE

320  subroutine bond_compute(nneig,neighl)
321 
322 
323 !This section has been created automatically by the script Abilint (TD).
324 !Do not modify the following lines by hand.
325 #undef ABI_FUNC
326 #define ABI_FUNC 'bond_compute'
327 !End of the abilint section
328 
329   implicit none
330 
331   !Arguments ------------------------
332   integer,intent(in) :: nneig(:)
333   integer,intent(in) :: neighl(:,:)
334   !Local ---------------------------
335   integer :: ii,jj,iat,jat
336   integer, allocatable, dimension(:,:) :: ibnd_dum
337   character(len=500) :: msg
338 
339 ! *************************************************************************
340 
341    ABI_ALLOCATE(ibnd_dum,(2,nfitmax*6))
342    ibnd_dum(:,:) = 0
343    ibn_tot  = 0    ! bonds between the fitted atoms
344    ibn_tot2 = 0    ! existing but non optimized bonds with border atoms
345    ibn_tots = 0    ! total number of bonds
346 
347 
348    if(nfit <= 0) then
349      write(msg,'(a,i8)')' UPDLIS WARNING : nfit <= 0 = ', nfit
350      MSG_WARNING(msg)
351    end if
352 
353    do ii = 1,nfit
354      iat = ifit(ii)
355      do jj = 1,nneig(iat)
356        jat = neighl(jj,iat)
357        if(tafit(jat))then
358          if(jat > iat) then  ! jat is a fitted atom
359            ibn_tot = ibn_tot + 1
360            ibnd_mat(:,ibn_tot) = (/ iat, jat /)
361          end if
362        else  ! jat is a border atom
363          ibn_tot2 = ibn_tot2 + 1
364          ibnd_dum(:,ibn_tot2) = (/ iat, jat /)
365        end if
366      end do
367    end do
368 
369    if(ibn_tot2 > (6*nfitmax)) then
370      write(msg,'(3a,i8,2a)')&
371 &     'ERROR: BOND_ATOM_INIT',ch10,&
372 &     'IBN_TOT2 =  ',ibn_tot2,ch10,&
373 &     ' ibnd_dum out of bounds, ibn_tot2 too large '
374      MSG_ERROR(msg)
375    end if
376 
377   !--reorder to keep 'variational' bonds first :
378    ibn_tots = ibn_tot + ibn_tot2
379    do ii = ibn_tot+1,ibn_tots
380      ibnd_mat(:,ii) = ibnd_dum(:,ii-ibn_tot)
381    end do
382 
383    ABI_DEALLOCATE (ibnd_dum)
384  end subroutine bond_compute

bond_lotf/bond_dealloc [ Functions ]

[ Top ] [ bond_lotf ] [ Functions ]

NAME

 bond_dealloc

FUNCTION

  deallocate variables

INPUTS

PARENTS

      m_lotf

CHILDREN

SOURCE

458  subroutine  bond_dealloc()
459 
460 ! *************************************************************************
461 
462 !This section has been created automatically by the script Abilint (TD).
463 !Do not modify the following lines by hand.
464 #undef ABI_FUNC
465 #define ABI_FUNC 'bond_dealloc'
466 !End of the abilint section
467 
468    ABI_DEALLOCATE(ibnd_mat)
469    ABI_DEALLOCATE(tafit)
470    ABI_DEALLOCATE(ifit)
471    ABI_DEALLOCATE(ibmat)
472    ABI_DEALLOCATE(ibmat_large)
473    ABI_DEALLOCATE(imat)
474  end subroutine bond_dealloc
475 
476 end module bond_lotf

bond_lotf/bond_fit_set [ Functions ]

[ Top ] [ bond_lotf ] [ Functions ]

NAME

 bond_fit_set

FUNCTION

  set nfitmax (or control it), ifit,nfit

INPUTS

PARENTS

      m_lotf

CHILDREN

SOURCE

402  subroutine bond_fit_set(nax,nfitdum)
403 
404 
405 !This section has been created automatically by the script Abilint (TD).
406 !Do not modify the following lines by hand.
407 #undef ABI_FUNC
408 #define ABI_FUNC 'bond_fit_set'
409 !End of the abilint section
410 
411   implicit none
412 
413   !Arguments ------------------------
414   integer,intent(in) :: nax
415   integer,intent(in) :: nfitdum
416   !Local -----------------------------
417   integer  :: i
418   character(len=500) :: msg
419 
420 ! *************************************************************************
421 
422    if(.not. allocated(ifit)) then
423      nfitmax = nfitdum + 1000
424      ABI_ALLOCATE(ifit,(nfitmax))
425    elseif(nfitdum > nfitmax) then
426      write(msg,'(a)')' BOND_FIT_SET : PROBLEM OF dimensionS !! '
427      MSG_ERROR(msg)
428    end if
429 
430    ifit = 0
431    nfit = 0
432    do i = 1, nax
433      if(tafit(i)) then
434        nfit = nfit + 1
435        ifit(nfit) = i
436      end if
437    end do
438 
439  end subroutine bond_fit_set

bond_lotf/bond_matrix_alloc [ Functions ]

[ Top ] [ bond_lotf ] [ Functions ]

NAME

 bond_matrix_alloc

FUNCTION

  allocate imat,ibmat,ibmat_large

INPUTS

PARENTS

      m_lotf

CHILDREN

SOURCE

209  subroutine bond_matrix_alloc(nax,nneigx)
210 
211 
212 !This section has been created automatically by the script Abilint (TD).
213 !Do not modify the following lines by hand.
214 #undef ABI_FUNC
215 #define ABI_FUNC 'bond_matrix_alloc'
216 !End of the abilint section
217 
218   implicit none
219 
220   !Arguments ------------------------
221   integer,intent(in) :: nax
222   integer,intent(in) :: nneigx
223 
224 ! *************************************************************************
225 
226    ABI_ALLOCATE(ibmat,(nneigx,0:nfitmax))
227    ABI_ALLOCATE(ibmat_large,(nfitmax,nfitmax))
228    ABI_ALLOCATE(imat,(nax))
229  end subroutine bond_matrix_alloc

bond_lotf/bond_matrix_set [ Functions ]

[ Top ] [ bond_lotf ] [ Functions ]

NAME

 bond_matrix_set

FUNCTION

  Set or update bond matrix imat,ibmat,ibmat_large
  associates the bond to the atom neighlists

INPUTS

PARENTS

      m_lotf

CHILDREN

SOURCE

248  subroutine bond_matrix_set(nneig,neighl)
249 
250 
251 !This section has been created automatically by the script Abilint (TD).
252 !Do not modify the following lines by hand.
253 #undef ABI_FUNC
254 #define ABI_FUNC 'bond_matrix_set'
255 !End of the abilint section
256 
257   implicit none
258 
259   !Arguments ------------------------
260   integer,intent(in) :: nneig(:)
261   integer,intent(in) :: neighl(:,:)
262   !Local ---------------------------
263   integer :: ibn,ibn2,iat,jat,ii,jj,j2
264 
265 ! *************************************************************************
266 
267    ibmat(:,:) = 0
268    ibmat_large(:,:) = 0
269    imat(:) = 0
270 
271    ibn = 0
272    do ii =1,nfit
273      iat = ifit(ii)
274      do jj = 1,nneig(iat)
275        jat = neighl(jj,iat)
276        if(jat > iat.AND.tafit(jat)) then
277          ibn  = ibn  + 1
278          ibmat(jj,ii) = ibn
279        end if
280      end do
281      imat(iat) = ii
282    end do
283 
284   !--ibmat_large: useful for the glue potential
285    ibn2 = 0
286    do ii =1,nfit
287      iat = ifit(ii)
288      do jj = 1,nneig(iat)
289        jat = neighl(jj,iat)
290        if(jat > iat.AND.tafit(jat)) then
291          ibn2  = ibn2  + 1
292          j2 = imat(jat)
293          ibmat_large(j2,ii) = ibn2
294          ibmat_large(ii,j2) = ibn2
295        end if
296      end do
297    end do
298  end subroutine bond_matrix_set

bond_lotf/bond_tafit_init [ Functions ]

[ Top ] [ bond_lotf ] [ Functions ]

NAME

 bond_tafit_init

FUNCTION

INPUTS

PARENTS

      m_lotf

CHILDREN

SOURCE

74  subroutine bond_tafit_init(nax)
75 
76 
77 !This section has been created automatically by the script Abilint (TD).
78 !Do not modify the following lines by hand.
79 #undef ABI_FUNC
80 #define ABI_FUNC 'bond_tafit_init'
81 !End of the abilint section
82 
83   implicit none
84 
85   !Arguments ------------------------
86   integer,intent(in) :: nax
87 
88 ! *************************************************************************
89 
90    ABI_ALLOCATE(tafit,(nax))
91 
92   !--MMANCINI strange!!!!!
93   !  tafit(:nax) = tquant(:nax)
94    tafit(:nax) = .true.
95 
96  end subroutine bond_tafit_init