TABLE OF CONTENTS
ABINIT/findminscf [ Functions ]
NAME
findminscf
FUNCTION
Compute the minimum of a function whose value and derivative are known at two points, using different algorithms. Also deduce different quantities at this predicted point, and at the two other points
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (XG, GMR) 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
choice=1,uses a linear interpolation of the derivatives =2,uses a quadratic interpolation based on the values of the function, and the second derivative at mid-point etotal_1=first value of the function etotal_2=second value of the function dedv_1=first value of the derivative dedv_2=second value of the derivative lambda_1=first value of the argument lambda_2=second value of the argument
OUTPUT
dedv_predict=predicted value of the derivative (usually zero, except if choice=4, if it happens that a minimum cannot be located, and a trial step is taken) d2edv2_predict=predicted value of the second derivative (not if choice=4) d2edv2_1=first value of the second derivative (not if choice=4) d2edv2_2=second value of the second derivative (not if choice=4) etotal_predict=predicted value of the function lambda_predict=predicted value of the argument status= 0 if everything went normally ; 1 if negative second derivative 2 if some other problem
PARENTS
scfcge
CHILDREN
wrtout
SOURCE
54 #if defined HAVE_CONFIG_H 55 #include "config.h" 56 #endif 57 58 #include "abi_common.h" 59 60 61 subroutine findminscf(choice,dedv_1,dedv_2,dedv_predict,& 62 & d2edv2_1,d2edv2_2,d2edv2_predict,& 63 & etotal_1,etotal_2,etotal_predict,& 64 & lambda_1,lambda_2,lambda_predict,errid,errmess) 65 66 use defs_basis 67 use m_profiling_abi 68 use m_errors 69 70 !This section has been created automatically by the script Abilint (TD). 71 !Do not modify the following lines by hand. 72 #undef ABI_FUNC 73 #define ABI_FUNC 'findminscf' 74 use interfaces_14_hidewrite 75 !End of the abilint section 76 77 implicit none 78 79 !Arguments ------------------------------------ 80 !scalars 81 integer,intent(in) :: choice 82 integer,intent(out) :: errid 83 character(len=500), intent(out) :: errmess 84 real(dp),intent(in) :: dedv_1,dedv_2,etotal_1,etotal_2,lambda_1,lambda_2 85 real(dp),intent(out) :: d2edv2_1,d2edv2_2,d2edv2_predict,dedv_predict 86 real(dp),intent(out) :: etotal_predict,lambda_predict 87 88 !Local variables------------------------------- 89 !scalars 90 real(dp) :: cc,d2edv2_mid,d_lambda,dedv_2bis 91 real(dp) :: dedv_mid2,etotal_2bis 92 character(len=500) :: message 93 94 ! ************************************************************************* 95 96 !DEBUG 97 !write(std_out,*)' findmin : enter' 98 !write(std_out,*)' choice,lambda_1,lambda_2=',choice,lambda_1,lambda_2 99 !ENDDEBUG 100 101 errid = AB7_NO_ERROR 102 d_lambda=lambda_1-lambda_2 103 104 if(choice==1) then 105 106 ! Use the derivative information to predict lambda 107 d2edv2_mid=(dedv_1-dedv_2)/d_lambda 108 lambda_predict=lambda_2-dedv_2/d2edv2_mid 109 dedv_predict=dedv_2+(lambda_predict-lambda_2)*d2edv2_mid 110 d2edv2_1=d2edv2_mid 111 d2edv2_2=d2edv2_mid 112 d2edv2_predict=d2edv2_mid 113 ! also use the first energy to predict new energy 114 etotal_predict=etotal_1+dedv_1*(lambda_predict-lambda_1)& 115 & +0.5_dp*d2edv2_1*(lambda_predict-lambda_1)**2 116 etotal_2bis=etotal_1+dedv_1*(lambda_2-lambda_1)& 117 & +0.5_dp*d2edv2_1*(lambda_2-lambda_1)**2 118 119 if(d2edv2_mid<0.0_dp)then 120 errid = AB7_ERROR_MIXING_INTERNAL 121 write(errmess,'(a,es18.10,a)')'The second derivative is negative, equal to ',d2edv2_mid,'.' 122 MSG_WARNING(errmess) 123 end if 124 125 else if(choice==2) then 126 127 ! Use energies and first derivative information 128 ! etotal = aa + bb * lambda + cc * lambda**2 129 dedv_mid2=(etotal_1-etotal_2)/d_lambda 130 cc=(dedv_1-dedv_mid2)/d_lambda 131 lambda_predict=lambda_1-0.5_dp*dedv_1/cc 132 d2edv2_1=2*cc 133 d2edv2_2=d2edv2_1 134 d2edv2_predict=d2edv2_1 135 if(d2edv2_predict<0.0_dp)then 136 errid = AB7_ERROR_MIXING_INTERNAL 137 write(errmess, '(a,es18.10,a,a,a)' )& 138 & 'The second derivative is negative, equal to',d2edv2_predict,'.',ch10,& 139 & '=> Pivoting ' 140 MSG_WARNING(errmess) 141 if(etotal_2 < etotal_1)then 142 lambda_predict=lambda_2-0.5_dp*(lambda_1-lambda_2) 143 else 144 lambda_predict=lambda_1-0.5_dp*(lambda_2-lambda_1) 145 end if 146 end if 147 dedv_predict=dedv_1+(lambda_predict-lambda_1)*d2edv2_1 148 dedv_2bis=dedv_1+(lambda_2-lambda_1)*d2edv2_1 149 etotal_predict=etotal_1+dedv_1*(lambda_predict-lambda_1)& 150 & +0.5_dp*d2edv2_1*(lambda_predict-lambda_1)**2 151 152 end if 153 154 write(message, '(a,es12.4,a,es18.10)' ) & 155 & ' findmin : lambda_predict ',lambda_predict,' etotal_predict ',etotal_predict 156 call wrtout(std_out,message,'COLL') 157 158 end subroutine findminscf