TABLE OF CONTENTS
ABINIT/chkdilatmx [ Functions ]
NAME
chkdilatmx
FUNCTION
Check whether the new rprimd does not give a too large number of plane waves, compared to the one booked for rprimd, taking into account the maximal dilatation dilatmx. Actually check whether the new Fermi sphere is inside the old one, dilated.
COPYRIGHT
Copyright (C) 2003-2017 ABINIT group (XG) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
chkdilatmx_ = if 1, will prevent to have any vector outside the Fermi sphere, possibly by rescaling (three times at most), and then stopping the execution if 0, simply send a warning, but continues execution dilatmx = maximal dilatation factor (usually the input variable) rprimd = new primitive vectors rprimd_orig = original primitive vectors (usually the input variable)
OUTPUT
dilatmx_errmsg=Emptry string if calculation can continue. If the calculation cannot continue, dilatmx_errmsg will contain the message that should be reported in the output file. Client code should handle a possible problem with the following test: if (LEN_TRIM(dilatmx_errmsg) then dump dilatmx_errmsg to the main output file. handle_error end if
PARENTS
driver,mover
CHILDREN
matr3eigval,matr3inv
SOURCE
47 #if defined HAVE_CONFIG_H 48 #include "config.h" 49 #endif 50 51 #include "abi_common.h" 52 53 54 subroutine chkdilatmx(chkdilatmx_,dilatmx,rprimd,rprimd_orig,dilatmx_errmsg) 55 56 use defs_basis 57 use m_errors 58 use m_profiling_abi 59 60 !This section has been created automatically by the script Abilint (TD). 61 !Do not modify the following lines by hand. 62 #undef ABI_FUNC 63 #define ABI_FUNC 'chkdilatmx' 64 use interfaces_32_util 65 !End of the abilint section 66 67 implicit none 68 69 !Arguments ------------------------------------ 70 !scalars 71 integer,intent(in) :: chkdilatmx_ 72 real(dp),intent(in) :: dilatmx 73 character(len=500),intent(out) :: dilatmx_errmsg 74 !arrays 75 real(dp),intent(inout) :: rprimd(3,3) 76 real(dp),intent(in) :: rprimd_orig(3,3) 77 78 !Local variables------------------------------- 79 !scalars 80 integer :: ii,jj,mu 81 real(dp) :: dilatmx_new 82 !arrays 83 real(dp) :: eigval(3),gprimd_orig(3,3),met(3,3),old_to_new(3,3) 84 real(dp) :: eigval_orig(3), alpha 85 character(len=500) :: message 86 87 ! ************************************************************************* 88 89 !Generates gprimd 90 call matr3inv(rprimd_orig,gprimd_orig) 91 92 !Find the matrix that transform an original xcart to xred, then to the new xcart 93 do mu=1,3 94 old_to_new(mu,:)=rprimd(mu,1)*gprimd_orig(:,1)+& 95 & rprimd(mu,2)*gprimd_orig(:,2)+& 96 & rprimd(mu,3)*gprimd_orig(:,3) 97 end do 98 99 !The largest increase in length will be obtained thanks 100 !to the diagonalization of the corresponding metric matrix : 101 !it is the square root of its largest eigenvalue. 102 do ii=1,3 103 do jj=1,3 104 met(ii,jj)=old_to_new(1,ii)*old_to_new(1,jj)+& 105 & old_to_new(2,ii)*old_to_new(2,jj)+& 106 & old_to_new(3,ii)*old_to_new(3,jj) 107 end do 108 end do 109 110 call matr3eigval(eigval,met) 111 112 dilatmx_new=sqrt(maxval(eigval(:))) 113 114 dilatmx_errmsg = "" 115 if(dilatmx_new>dilatmx+tol6)then 116 117 ! MJV 2014 07 22: correct rprim to maximum jump allowed by dilatmx 118 ! XG 20171011 : eigenvalues of "old_to_old" tensor are of course the unity ! 119 120 if(chkdilatmx_/=0)then 121 alpha = (dilatmx - one) / (dilatmx_new - one) 122 ! for safety, only 90 percent of max jump 123 alpha = 0.9_dp * alpha 124 125 rprimd = alpha * rprimd + (one - alpha) * rprimd_orig 126 127 write(dilatmx_errmsg,'(3a,es16.6,4a,es16.6,2a,es16.6,a)')& 128 & 'The new primitive vectors rprimd (an evolving quantity)',ch10,& 129 & 'are too large with respect to the old rprimd and the accompanying dilatmx:',dilatmx,ch10,& 130 & 'This large change of unit cell parameters is not allowed by the present value of dilatmx.',ch10,& 131 & 'An adequate value would have been dilatmx_new=',dilatmx_new,ch10,& 132 & 'Calculation continues with limited jump, by rescaling the projected move by the factor',alpha,'.' 133 else 134 write(message, '(3a,es16.6,2a,es16.6,2a)' )& 135 & 'The new primitive vectors rprimd (an evolving quantity)',ch10,& 136 & 'are too large, given the initial rprimd and the accompanying dilatmx:',dilatmx,ch10,& 137 & 'An adequate value would have been dilatmx_new=',dilatmx_new,ch10,& 138 & 'As chkdilatmx=0, assume experienced user. Execution will continue.' 139 MSG_WARNING(message) 140 endif 141 142 end if 143 144 end subroutine chkdilatmx