TABLE OF CONTENTS


ABINIT/cgpr [ Functions ]

[ Top ] [ Functions ]

NAME

 cgpr

FUNCTION

 perform Polak-Ribiere conjugate gradient on a function f
 implementation based on the cg recipe of "numerical recipe"

COPYRIGHT

 Copyright (C) 1998-2018 ABINIT group (DCA, XG, MT)
 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

 dp_dum_vdp: function  to be minimized (return a dp from a vector of dp)
 vdp_dum_vdp: derivative of f
 dtol: precision precision required for the minimization
 itmax: number of iterations allowed (each linmin will be done with at max 10 times
 this number

OUTPUT

 fmin: value of f at the minimum
 lastdelta: absolute value of the last delta between steps

SIDE EFFECTS

 v: vector on which minimization is to be performed, starting point
 and resulting min

PARENTS

      prcrskerker1,prcrskerker2

CHILDREN

      linmin

SOURCE

 39 #if defined HAVE_CONFIG_H
 40 #include "config.h"
 41 #endif
 42 
 43 #include "abi_common.h"
 44 
 45 
 46 subroutine cgpr(nv1,nv2,dp_dum_v2dp,v2dp_dum_v2dp,sub_dum_dp_v2dp_v2dp,dtol,itmax,v,fmin,delta)
 47 
 48  use defs_basis
 49  use m_profiling_abi
 50 
 51 !This section has been created automatically by the script Abilint (TD).
 52 !Do not modify the following lines by hand.
 53 #undef ABI_FUNC
 54 #define ABI_FUNC 'cgpr'
 55  use interfaces_62_cg_noabirule, except_this_one => cgpr
 56 !End of the abilint section
 57 
 58  implicit none
 59 
 60 !Arguments ------------------------------------
 61  include "dummy_functions.inc"
 62 !scalars
 63  integer,intent(in) :: itmax,nv1,nv2
 64  real(dp),intent(in) :: dtol
 65  real(dp),intent(out) :: delta,fmin
 66 !arrays
 67  real(dp),intent(inout) :: v(nv1,nv2)
 68 
 69 !Local variables-------------------------------
 70 !scalars
 71  integer :: iiter
 72  real(dp) :: fv,gam,gscal,gscal2,sto
 73 !arrays
 74  real(dp) :: grad0(nv1,nv2),grad1(nv1,nv2),grad2(nv1,nv2),grad3(nv1,nv2)
 75 !no_abirules
 76 
 77 !************************************************************************
 78  fv = dp_dum_v2dp(nv1,nv2,v(:,:))
 79  grad0(:,:) = -v2dp_dum_v2dp(nv1,nv2,v(:,:))
 80  grad1(:,:) = grad0(:,:)
 81  grad2(:,:) = grad0(:,:)
 82  do iiter=1,itmax
 83   call linmin(nv1,nv2,dp_dum_v2dp,v2dp_dum_v2dp,sub_dum_dp_v2dp_v2dp,v,grad0,fmin)
 84 ! return if the min is reached
 85   sto=dtol*(abs(fmin)+abs(fv)+tol14)
 86   delta=abs(fv-fmin)
 87   delta=abs(delta)
 88   if((delta.lt.sto).or.(iiter==itmax)) then
 89 !  DEBUG
 90 !  write(std_out,*) 'cgpr (01cg) : stop cond for cgpr:',sto,'delta:',delta,'fv:',fv,'fmin:',fmin
 91 !  ENDDEBUG
 92    return
 93   end if
 94 ! a new step
 95   fv=fmin
 96   grad0(:,:)=v2dp_dum_v2dp(nv1,nv2,v(:,:))
 97   gscal=dotproduct(nv1,nv2,grad1(:,:),grad1(:,:))
 98   grad3(:,:)=grad0(:,:)+grad1(:,:)
 99   gscal2=dotproduct(nv1,nv2,grad3(:,:),grad0(:,:))
100   gam=gscal2/gscal
101   grad1(:,:)=-grad0(:,:)
102   grad2(:,:)=grad1(:,:)+gam*grad2(:,:)
103   grad0(:,:)=grad2(:,:)
104 ! DEBUG
105 ! write(std_out,*) 'cgpr (01cg) :================================================================================='
106 ! write(std_out,*) 'cgpr (01cg) : step',iiter,'delta:',delta ,'fv',fv,'fmin',fmin
107 ! write(std_out,*) 'cgpr (01cg) :================================================================================='
108 ! ENDDEBUG
109  end do
110 end subroutine cgpr