TABLE OF CONTENTS


ABINIT/prcrskerker1 [ Functions ]

[ Top ] [ Functions ]

NAME

 prcrskerker1

FUNCTION

 preconditionning by a real-space conjugate gradient on residual
 using a model dielectric function in real space

COPYRIGHT

 Copyright (C) 1998-2018 ABINIT group (PMA)
 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/Infos/contributors .

INPUTS

  nfft=number of fft grid points
  nspden=number of spin-density components
  ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft
  dielar(7)=input parameters for dielectric matrix:
                diecut,dielng,diemac,diemix,diegap,dielam,diemixmag.
  gprimd(3,3)=dimensional primitive translations in fourier space (bohr**-1)
  rprimd(3,3)=dimensional primitive translations in real space (bohr)
  vresid(nfft,nspden)=residual potential
  base(nfft) = real space function used as a basis to guess a fine dielectric funtion
  see the calling routine to know the content

OUTPUT

  vrespc(nfft,nspden)=preconditioned residual of the potential

SIDE EFFECTS

WARNINGS

 This is experimental code : input, ouptput, results and any other feature may vary greatly.

NOTES

  needs severe cleaning and this is abuse of modules as common blocks...

PARENTS

      prcref,prcref_PMA

CHILDREN

      cgpr,frskerker1__end,frskerker1__init,laplacian,prc_mem_init

SOURCE

 48 #if defined HAVE_CONFIG_H
 49 #include "config.h"
 50 #endif
 51 
 52 #include "abi_common.h"
 53 
 54 
 55 subroutine prcrskerker1(dtset,mpi_enreg,nfft,nspden,ngfft,dielar,etotal,gprimd,vresid,vrespc,base)
 56 
 57 
 58  use defs_basis
 59  use defs_abitypes
 60  use frskerker1
 61  use mod_prc_memory
 62  use m_profiling_abi
 63 
 64 !This section has been created automatically by the script Abilint (TD).
 65 !Do not modify the following lines by hand.
 66 #undef ABI_FUNC
 67 #define ABI_FUNC 'prcrskerker1'
 68  use interfaces_56_recipspace
 69  use interfaces_62_cg_noabirule
 70 !End of the abilint section
 71 
 72  implicit none
 73 
 74 !Arguments ------------------------------------
 75 !scalars
 76  integer,intent(in) :: nfft,nspden
 77  real(dp) :: etotal
 78  type(MPI_type),intent(in) :: mpi_enreg
 79  type(dataset_type),intent(in) :: dtset
 80 !arrays
 81  integer,intent(in) :: ngfft(18)
 82  real(dp),intent(in) :: base(nfft),dielar(7),gprimd(3,3)
 83  real(dp),intent(in) :: vresid(nfft,nspden)
 84  real(dp),intent(out) :: vrespc(nfft,nspden)
 85 
 86 !Local variables-------------------------------
 87 !scalars
 88  integer ::  ifft,ispden,n1,n2,n3
 89  real(dp) :: base_delta,base_max,base_min,dielng,diemac,diemix
 90  real(dp) :: diemixmag
 91  real(dp) :: rdummy1,rdummy2
 92  logical ::  new_prc_func
 93 !arrays
 94  real(dp) :: deltaW(nfft,nspden)
 95  real(dp) :: g2cart(nfft)
 96  real(dp) :: mat(nfft,nspden)
 97 
 98 ! *************************************************************************
 99 
100 !DEBUG
101 !write(std_out,*)' prckerker1 : enter '
102 !ENDDEBUG
103 !if(cycle==0) then
104  call prc_mem_init(nfft)
105 
106  if(cycle==0) then
107    new_prc_func=.TRUE.
108    energy_min=etotal
109  else if(etotal < energy_min) then
110    new_prc_func=.TRUE.
111    energy_min=etotal
112  else
113    new_prc_func=.FALSE.
114  end if
115 
116 
117  dielng=dielar(2)
118  diemac=dielar(3)
119  diemix=dielar(4)
120  diemixmag=dielar(7)
121 !******************************************************************
122 !compute the diemac(r)                                          **
123 !******************************************************************
124 !this task will be devoted to a general function later
125  n1=ngfft(1)
126  n2=ngfft(2)
127  n3=ngfft(3)
128 !base_cp=base
129  if(new_prc_func) then
130    base_min=base(1)
131    base_max=base(1)
132    do ifft=1,nfft
133      base_min = min(base_min,base(ifft))
134      base_max = max(base_max,base(ifft))
135    end do
136    base_delta = base_max - base_min
137 !  if(cycle.lt.2) then
138    rdiemac(:) = (((base(:)-base_min) / (base_delta) ) *(diemac-one) + one)
139 !  else
140 !  rdiemac(:) = rdiemac(:)*0.5_dp+0.5_dp*(((base(:)-base_min) / (base_delta) ) *(diemac-one) + one)
141 !  end if
142 !  if(cycle==0) rdiemac(:) = (((base(:)-base_min) / (base_delta) ) *(diemac-one) + one)
143 !  rdiemac(:) = exp(((base(:)-base_min) / (base_delta) *log(diemac)))
144  end if
145  cycle=cycle+1
146 !if(cycle==5) cycle=0
147 !end if
148 !******************************************************************
149 !compute deltaW                                                 **
150 !******************************************************************
151  vrespc=vresid !starting point
152  call laplacian(gprimd,mpi_enreg,nfft,nspden,ngfft,dtset%paral_kgb,rdfuncr=vrespc,laplacerdfuncr=deltaW,g2cart_out=g2cart) ! put the laplacian of the residuals into deltaW
153 !call laplacian(vrespc,buffer,ngfft,gprimd) ! put the laplacian of the residuals into deltaW
154 !do ifft=1,nfft
155 !if (buffer(ifft,1)/=deltaW(ifft,1)) then
156 !stop
157 !end if
158 !end do
159  deltaW(:,1)= diemix*(((one/rdiemac(:))*vresid(:,1))-(((dielng)**2)*deltaW(:,1)))
160  if (nspden>1.and.(diemixmag>=zero)) then
161    do ispden=2,nspden
162      deltaW(:,ispden)= abs(diemixmag)*(((one/rdiemac(:))*vresid(:,ispden))-(((dielng)**2)*deltaW(:,ispden)))
163    end do
164  end if
165 !call random_number(deltaW)
166 !call random_number(vrespc)
167 !******************************************************************
168 !Finding the preconditionned residuals which minimizes          **
169 !half*(vrespc*(1-dielng2/4pi2 nabla2) vrespc) - vrespc * deltaW **
170 !***********************************************************************
171  vrespc(:,1)=diemix*vrespc(:,1)
172  if (nspden>1) vrespc(:,2:nspden)=abs(diemixmag)*vrespc(:,2:nspden)
173 !buffer=vrespc
174 
175 
176 !==============================================================================
177 !==============================================================================
178 !! Original loop
179 !==============================================================================
180 !==============================================================================
181 
182  call frskerker1__init(dtset,mpi_enreg,nfft,ngfft,nspden,dielng,deltaW,gprimd,mat,g2cart)
183 
184 !call cgpr(pf_rscgres,dpf_rscgres,newvres,real(1e-40,dp),700,vrespc,rdummy1,rdummy2)
185 !rdummy1 = pf_rscgres(nfft,nspden,vrespc)
186  call cgpr(nfft,nspden,frskerker1__pf,frskerker1__dpf,frskerker1__newvres,&
187 & real(1e-10,dp),700,vrespc,rdummy1,rdummy2)
188  call frskerker1__end()
189 
190 !==============================================================================
191 !==============================================================================
192 !! Original loop end
193 !==============================================================================
194 !==============================================================================
195 
196 
197 !cplex=1
198 !qphon(:)=zero
199 !call moddiel(cplex,dielar,nfft,ngfft,nspden,1,0,qphon,rprimd,vresid,buffer)
200 !c1=0
201 !do ifft=1,nfft,1
202 !if((abs(buffer(ifft,1)-vrespc(ifft,1))/(abs(buffer(ifft,1)+vrespc(ifft,1))*half)) > 5e-3) then
203 !c1=c1+1
204 !end if
205 !end do
206 !call laplacian(vrespc,buffer,ngfft,gprimd)
207 !buffer=vrespc(:,:)-buffer(:,:)*dielng**2
208 !c2=0
209 !do ifft=1,nfft,1
210 !if((abs(buffer(ifft,1)-deltaW(ifft,1))/(abs(buffer(ifft,1)+deltaW(ifft,1))*half)) > 5e-3) then
211 !c2=c2+1
212 !end if
213 !end do
214 !!!  !stop
215 !call laplacian(gprimd,mpi_enreg,nfft,nspden,ngfft,&
216 !& g2cart_out=g2cart)
217 
218 !vrespc=vresid
219 !do ispden=1,nspden
220 !call fourdp(1, gvrespc(:,:,ispden), vrespc(:,ispden),-1,mpi_enreg,nfft,ngfft,0)
221 !end do
222 !filtering
223 !do ispden=1,nspden
224 !do ifft=1,nfft
225 !!    gvrespc(:,ifft,ispden)=(one-exp(-g2cart(ifft)*15.0_dp))*gvrespc(:,ifft,ispden)
226 !!      gvrespc(:,ifft,ispden)=(exp(-g2cart(ifft)*10.0_dp))*gvrespc(:,ifft,ispden)
227 !!      gvrespc(:,ifft,ispden)=(one-one/(exp(-0.002_dp/g2cart(ifft)**2)+one))*gvrespc(:,ifft,ispden)
228 !gvrespc(:,ifft,ispden)=(two-2_dp/(exp(-0.008_dp/(g2cart(ifft)+0.0012_dp))+one))*gvrespc(:,ifft,ispden)
229 !gvrespc(:,ifft,ispden)=min(one,(sqrt(g2cart(ifft)/0.006_dp))**(one))*gvrespc(:,ifft,ispden)
230 !end do
231 !end do
232 !change resulting potential to real space
233 !do ispden=1,nspden
234 !call fourdp(1,gvrespc(:,:,ispden),vrespc(:,ispden),1,mpi_enreg,nfft,ngfft,0)
235 !end do
236 !vrespc=vrespc*diemix
237 !maxg2=g2cart(1)
238 !ming2=g2cart(5)
239 !do ifft=1,nfft
240 !maxg2=max(g2cart(ifft),maxg2)
241 !if(g2cart(ifft) .gt. zero) ming2=min(g2cart(ifft),ming2)
242 !end do
243 !stop
244 
245 !DEBUG
246 !write(std_out,*)' prckerker1 : exit '
247 !ENDDEBUG
248 
249 end subroutine prcrskerker1