TABLE OF CONTENTS


ABINIT/density_rec [ Functions ]

[ Top ] [ Functions ]

NAME

 density_rec

FUNCTION

 This routine computes the density using  the coefficients corresponding to
 continued fraction at a point from a fixed potential. 

COPYRIGHT

 Copyright (C) 2008-2018 ABINIT group (SLeroux,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 .

INPUTS

  coordx, coordy, coordz=coordonnees of the computed point
  an, bn2 : coefficient given by density_rec. Input if get_rec_coef=0, output else
  nrec=order of density_rec
  fermie=fermi energy (Hartree)
  tsmear=temperature (Hartree)
  rtrotter=real trotter parameter
  tol=tolerance criteria for stopping density_rec
  inf_ucvol=infinitesimal unit cell volume
  dim_trott = max(0,2*trotter-1)

OUTPUT

  rho_out=result of the continued fraction multiplied by a multiplicator

SIDE EFFECTS

PARENTS

      fermisolverec

CHILDREN

      timab,trottersum

NOTES

  at this time :
       - exppot should be replaced by ?
       - coord should be replaced by ?
       - need a rectangular box (rmet diagonal matrix)

SOURCE

 47 #if defined HAVE_CONFIG_H
 48 #include "config.h"
 49 #endif
 50 
 51 #include "abi_common.h"
 52 
 53 
 54 subroutine density_rec(an,bn2,rho_out,nrec, &
 55 &                     fermie,tsmear,rtrotter, &
 56 &                     dim_trott,tol,inf_ucvol)
 57 
 58  use m_profiling_abi
 59  
 60  use defs_basis
 61  use m_rec_tools,only :       trottersum
 62 
 63 !This section has been created automatically by the script Abilint (TD).
 64 !Do not modify the following lines by hand.
 65 #undef ABI_FUNC
 66 #define ABI_FUNC 'density_rec'
 67  use interfaces_18_timing
 68 !End of the abilint section
 69 
 70  implicit none
 71 
 72 !Arguments -------------------------------
 73 !scalars
 74  integer,intent(in) :: nrec
 75  integer,intent(in) :: dim_trott
 76  real(dp),intent(in) :: fermie,tol,tsmear,inf_ucvol,rtrotter
 77  real(dp), intent(out) :: rho_out
 78 !arrays
 79  real(dp),intent(in) :: an(0:nrec),bn2(0:nrec)
 80 !Local variables-------------------------------
 81 !not used, debugging purpose only
 82 !for debugging purpose, detailled printing only once for density and ekin
 83 !scalars
 84  integer, parameter :: minrec = 3
 85  integer  :: irec
 86  real(dp) :: beta,mult,prod_b2,error,errold
 87  real(dp) :: pi_on_rtrotter,twortrotter,exp1,exp2
 88  complex(dpc) :: cinv2rtrotter,coeef_mu,facrec0
 89 ! character(len=500) :: msg
 90 !arrays
 91  real(dp) :: tsec(2)
 92  complex(dpc) :: acc_rho(0:nrec)
 93  complex(dpc) :: D(0:dim_trott),Dold(0:dim_trott)
 94  complex(dpc) :: N(0:dim_trott),Nold(0:dim_trott)
 95 !**************************************************************************
 96 
 97  call timab(605,1,tsec)
 98  
 99 !##############################################################
100 !--Initialisation of metrics 
101  mult = two/inf_ucvol   !non-spined system
102  beta = one/tsmear
103 
104 !--Variables for optimisation
105  pi_on_rtrotter = pi/rtrotter
106  twortrotter = two*rtrotter
107  exp1 = exp((beta*fermie)/(rtrotter))
108  exp2 = exp(beta*fermie/(twortrotter))
109  cinv2rtrotter = cmplx(one/twortrotter,zero,dp)
110  coeef_mu = cmplx(one/exp2,zero,dp)
111 
112  N = czero;  D = cone
113  facrec0 = cone
114  Nold = czero; Dold = czero
115 !--Initialisation of accumulated density 
116  acc_rho = czero  
117 !--Initialisation of estimated error
118  prod_b2 = twortrotter/exp1
119  errold = zero
120  
121 
122 !##############################################################
123 !--Main loop
124  maindo : do irec = 0, nrec
125    
126 !  ######################################################
127 !  --Density computation
128 !  !--using the property that: sum_i(bi*c)^2|(z-ai*c)=1/c*sum_i(bi)^2|(z/c-ai)
129 !  !and for c =exp(-beta*fermie/(two*rtrotter)
130 
131    call trottersum(dim_trott,error,prod_b2,pi_on_rtrotter,&
132 &   facrec0,coeef_mu,exp1,&
133 &   an(irec),bn2(irec),&
134 &   N,D,Nold,Dold)
135    
136    if(irec/=nrec .and. irec>=minrec)then
137      if((bn2(irec+1)<tol14).or.(mult*error<tol.and.errold<tol)) exit maindo
138    end if
139    errold = mult*error
140  end do maindo
141 !--Accumulated density
142  rho_out = mult*real(cone-sum(N/D)*cinv2rtrotter,dp)
143  
144  call timab(605,2,tsec)
145  
146  end subroutine density_rec