TABLE OF CONTENTS


ABINIT/m_rec_tools [ Modules ]

[ Top ] [ Modules ]

NAME

  m_rec_tools

FUNCTION

  This module provides some functions more or less generic used
  in the Recursion Mathod

COPYRIGHT

 Copyright (C) 2002-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 .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt.

NOTES

PARENTS

CHILDREN

SOURCE

25 #if defined HAVE_CONFIG_H
26 #include "config.h"
27 #endif
28 
29 #include "abi_common.h"
30 
31 module m_rec_tools
32 
33  use defs_basis
34  use m_abicore
35 
36  use defs_rectypes, only : recparall_type
37 
38  implicit none
39 
40  private
41 
42  public ::            &
43    get_pt0_pt1,       &  !--To get pt0 pt1 from inf,sup
44    reshape_pot,       &  !--To rescale the potential between 2 grids
45    trottersum            !--To calculate the trotter sum
46 
47 CONTAINS  !===========================================================

m_rec_tools/get_pt0_pt1 [ Functions ]

[ Top ] [ m_rec_tools ] [ Functions ]

NAME

 get_pt0_pt1

FUNCTION

 utility function to get pt0 and pt1 for given inf and sup

INPUTS

 ngfft(3) = fine grid (corresponds to dtset%ngfft(1:3))
 inf = inferior point in-line coordinate on the coarse grid
 sup = superior point in-line coordinate on the coarse grid

OUTPUT

 recpar%pt0<type(vec_int)>=Intial point for this proc in x,y,z
 recpar%pt1<type(vec_int)>=Final point for this proc in x,y,z
 recpar%min_pt=inferior point in-line coordinate on the fine grid
 recpar%max_pt=superior point in-line coordinate on the dine grid

SIDE EFFECTS

PARENTS

      first_rec,m_rec

CHILDREN

SOURCE

 77 subroutine get_pt0_pt1(ngfft,gratio,inf,sup,recpar)
 78 
 79 
 80 !This section has been created automatically by the script Abilint (TD).
 81 !Do not modify the following lines by hand.
 82 #undef ABI_FUNC
 83 #define ABI_FUNC 'get_pt0_pt1'
 84 !End of the abilint section
 85 
 86  implicit none
 87 
 88 !Arguments ------------------------------------
 89   integer,intent(in)        :: gratio
 90   integer,intent(in)        :: inf,sup
 91   type(recparall_type),intent(inout) :: recpar
 92   integer,intent(in)        :: ngfft(3)
 93 !Local ---------------------------
 94   integer :: x,y,z,count
 95   integer :: boz,boy,pt
 96 ! *********************************************************************
 97   count = 0
 98   do z = 0,ngfft(3)-1,gratio
 99     boz = z*ngfft(2)
100     do y = 0,ngfft(2)-1,gratio
101       boy = (boz+y)*ngfft(1)
102       do x = 0,ngfft(1)-1,gratio
103         pt = boy+x
104         if(count >= inf) then
105           if(count == inf) then
106             recpar%pt0%x = x; recpar%pt0%y = y; recpar%pt0%z = z
107             recpar%min_pt = pt
108           end if
109           recpar%pt1%x = x; recpar%pt1%y = y; recpar%pt1%z = z
110           recpar%max_pt = pt
111         endif
112         count = count+1
113         if(count == sup ) return
114       end do
115     end do
116   end do
117 end subroutine get_pt0_pt1

m_rec_tools/reshape_pot [ Functions ]

[ Top ] [ m_rec_tools ] [ Functions ]

NAME

 reshape_pot

FUNCTION

 Reshape array on

INPUTS

 nfft=size of the big grid
 nfftrec=size of the cut grid
 trasl(3) center point
 ngfft(3) dimensions of the big grid
 ngfftrec(3) dimensions of the cut grid
 pot(0:nfft-1) 3d-array on the big grid

OUTPUT

 potloc is a cut of pot around trasl

PARENTS

      first_rec,nlenergyrec,vtorhorec

CHILDREN

SOURCE

145 subroutine reshape_pot(trasl,nfft,nfftrec,ngfft,ngfftrec,pot,potloc)
146 
147 
148 !This section has been created automatically by the script Abilint (TD).
149 !Do not modify the following lines by hand.
150 #undef ABI_FUNC
151 #define ABI_FUNC 'reshape_pot'
152 !End of the abilint section
153 
154  implicit none
155 
156  !Arguments ------------------------------------
157  integer,  intent(in) :: nfft,nfftrec
158  integer,  intent(in) :: trasl(3)
159  integer,  intent(in) :: ngfftrec(3),ngfft(3)
160  real(dp), intent(in) :: pot(0:nfft-1)
161  real(dp), intent(out):: potloc(0:nfftrec-1)
162  !Local ----------------------------------------
163  ! scalars
164  integer :: zz,yy,xx,parz,pary
165  integer :: modi,modj,modk
166  !character(len=500) :: msg
167  ! *********************************************************************
168  do zz = 0,ngfftrec(3)-1
169    modk = modulo(zz-trasl(3),ngfft(3))*ngfft(2)
170    parz = ngfftrec(2)*zz
171    do yy = 0,ngfftrec(2)-1
172      modj = (modulo(yy-trasl(2),ngfft(2))+modk)*ngfft(1)
173      pary = ngfftrec(1)*(yy+parz)
174      do xx = 0,ngfftrec(1)-1
175        modi = modulo(xx-trasl(1),ngfft(1))+modj
176        potloc(xx+pary) = pot(modi)
177      end do
178    end do
179  end do
180 
181 end subroutine reshape_pot

m_rec_tools/trottersum [ Functions ]

[ Top ] [ m_rec_tools ] [ Functions ]

NAME

 trottersum

FUNCTION

 Calculate the contribution to the partial fraction decomposition
 due to a recursion step.

INPUTS

  dim_trott=dimension of the partial fraction decomposition (PFD)
  pi_on_rtrotter=parameter pi/rtrotter
  an,bn2=recursion coefficients at irec
  exp1=numerical factor  (see density_rec)
  coeef_mu=numerical factor

OUTPUT

 SIZE EFFECTS
  D,N=denominator and numerator accumalator of PFD
  Dold,Nold=denominator and numerator of PFD (old values)
  facrec0=used to select irec=0
  error=estimated error of recursion at this step
  prod_b2=numerical factor

PARENTS

      density_rec,recursion,recursion_nl

CHILDREN

SOURCE

215 subroutine trottersum(dim_trott,error,&
216      &                prod_b2,pi_on_rtrotter,&
217      &                facrec0,coeef_mu,exp1,&
218      &                an,bn2,&
219      &                N,D,Nold,Dold)
220 
221 
222 !This section has been created automatically by the script Abilint (TD).
223 !Do not modify the following lines by hand.
224 #undef ABI_FUNC
225 #define ABI_FUNC 'trottersum'
226 !End of the abilint section
227 
228  implicit none
229 
230  !Arguments ------------------------------------
231  !scalars
232  integer,  intent(in) :: dim_trott
233  real(dp), intent(in) :: an,bn2,exp1,pi_on_rtrotter
234  real(dp), intent(inout) :: error,prod_b2
235  complex(dp), intent(in) :: coeef_mu
236  complex(dp), intent(inout) :: facrec0
237  !arrays
238  complex(dpc),intent(inout) :: D(0:dim_trott),Dold(0:dim_trott)
239  complex(dpc),intent(inout) :: N(0:dim_trott),Nold(0:dim_trott)
240  !Local ----------------------------------------
241  ! scalars
242  integer :: itrot
243  real(dp) :: arg
244  complex(dpc) :: Dnew,Nnew,zj
245  !character(len=500) :: msg
246  ! *********************************************************************
247 
248  error = zero
249  prod_b2 = prod_b2 * exp1 * bn2
250 
251  do itrot=0,dim_trott
252    arg = pi_on_rtrotter*(real( itrot,dp) + half )
253    zj = cmplx(cos(arg),sin(arg),dp)*coeef_mu
254 
255    Nnew = zj*facrec0 + &
256         & (zj - cmplx(an  ,zero,dp))*N(itrot) - &
257         &       cmplx(bn2 ,zero,dp)*Nold(itrot)
258    Dnew = (zj - cmplx(an  ,zero,dp))*D(itrot) - &
259         &       cmplx(bn2 ,zero,dp)*Dold(itrot)
260    Nold(itrot) = N(itrot)
261    Dold(itrot) = D(itrot)
262    N(itrot) = Nnew
263    D(itrot) = Dnew
264 
265    !--Error estimator
266    error = error + abs(prod_b2/(D(itrot)*Dold(itrot)))
267  end do
268  facrec0 = czero
269 
270 end subroutine trottersum