TABLE OF CONTENTS
ABINIT/holocell [ Functions ]
NAME
holocell
FUNCTION
Examine whether the trial conventional cell described by cell_base is coherent with the required holohedral group. Possibly enforce the holohedry and modify the basis vectors. Note : for iholohedry=4, the tetragonal axis is not required to be along the C axis.
COPYRIGHT
Copyright (C) 2000-2018 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 . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
enforce= if 0, only check; if =1, enforce exactly the holohedry iholohedry=required holohegral group iholohedry=1 triclinic 1bar iholohedry=2 monoclinic 2/m iholohedry=3 orthorhombic mmm iholohedry=4 tetragonal 4/mmm iholohedry=5 trigonal 3bar m iholohedry=6 hexagonal 6/mmm iholohedry=7 cubic m3bar m tolsym=tolerance for the symmetry operations
OUTPUT
foundc=1 if the basis vectors supports the required holohedry ; =0 otherwise
SIDE EFFECTS
cell_base(3,3)=basis vectors of the conventional cell (changed if enforce==1, otherwise unchanged)
PARENTS
symlatt,symmetrize_rprimd
CHILDREN
SOURCE
45 #if defined HAVE_CONFIG_H 46 #include "config.h" 47 #endif 48 49 #include "abi_common.h" 50 51 52 subroutine holocell(cell_base,enforce,foundc,iholohedry,tolsym) 53 54 use defs_basis 55 use m_errors 56 use m_profiling_abi 57 58 !This section has been created automatically by the script Abilint (TD). 59 !Do not modify the following lines by hand. 60 #undef ABI_FUNC 61 #define ABI_FUNC 'holocell' 62 !End of the abilint section 63 64 implicit none 65 66 !Arguments ------------------------------------ 67 !scalars 68 integer,intent(in) :: enforce,iholohedry 69 integer,intent(out) :: foundc 70 real(dp),intent(in) :: tolsym 71 !arrays 72 real(dp),intent(inout) :: cell_base(3,3) 73 74 !Local variables ------------------------------ 75 !scalars 76 integer :: allequal,ii,jj,orth 77 real(dp):: aa,reldiff,scprod1 78 character(len=500) :: msg 79 !arrays 80 integer :: ang90(3),equal(3) 81 real(dp) :: length(3),metric(3,3),norm(3),rbasis(3,3),rconv(3,3),rconv_new(3,3) 82 real(dp) :: rnormalized(3,3),symmetrized_length(3) 83 84 !************************************************************************** 85 86 do ii=1,3 87 metric(:,ii)=cell_base(1,:)*cell_base(1,ii)+& 88 & cell_base(2,:)*cell_base(2,ii)+& 89 & cell_base(3,:)*cell_base(3,ii) 90 end do 91 92 !Examine the angles and vector lengths 93 ang90(:)=0 94 if(metric(1,2)**2<tolsym**2*metric(1,1)*metric(2,2))ang90(3)=1 95 if(metric(1,3)**2<tolsym**2*metric(1,1)*metric(3,3))ang90(2)=1 96 if(metric(2,3)**2<tolsym**2*metric(2,2)*metric(3,3))ang90(1)=1 97 orth=0 98 if(ang90(1)==1 .and. ang90(2)==1 .and. ang90(3)==1) orth=1 99 equal(:)=0 100 if(abs(metric(1,1)-metric(2,2))<tolsym*half*(metric(1,1)+metric(2,2)))equal(3)=1 101 if(abs(metric(1,1)-metric(3,3))<tolsym*half*(metric(1,1)+metric(3,3)))equal(2)=1 102 if(abs(metric(2,2)-metric(3,3))<tolsym*half*(metric(2,2)+metric(3,3)))equal(1)=1 103 allequal=0 104 if(equal(1)==1 .and. equal(2)==1 .and. equal(3)==1) allequal=1 105 106 foundc=0 107 if(iholohedry==1) foundc=1 108 if(iholohedry==2 .and. ang90(1)+ang90(3)==2 ) foundc=1 109 if(iholohedry==3 .and. orth==1) foundc=1 110 if(iholohedry==4 .and. orth==1 .and. & 111 & (equal(3)==1 .or. equal(2)==1 .or. equal(1)==1) ) foundc=1 112 if(iholohedry==5 .and. allequal==1 .and. & 113 & (abs(metric(1,2)-metric(2,3))<tolsym*metric(2,2)) .and. & 114 & (abs(metric(1,2)-metric(1,3))<tolsym*metric(1,1)) ) foundc=1 115 if(iholohedry==6 .and. equal(3)==1 .and. & 116 & ang90(1)==1 .and. ang90(2)==1 .and. & 117 & (2*metric(1,2)-metric(1,1))<tolsym*metric(1,1) ) foundc=1 118 if(iholohedry==7 .and. orth==1 .and. allequal==1) foundc=1 119 120 !------------------------------------------------------------------------------------- 121 !Possibly enforce the holohedry (if it is to be enforced !) 122 123 if(foundc==1.and.enforce==1.and.iholohedry/=1)then 124 125 ! Copy the cell_base vectors, and possibly fix the tetragonal axis to be the c-axis 126 if(iholohedry==4.and.equal(1)==1)then 127 rconv(:,3)=cell_base(:,1) ; rconv(:,1)=cell_base(:,2) ; rconv(:,2)=cell_base(:,3) 128 else if (iholohedry==4.and.equal(2)==1)then 129 rconv(:,3)=cell_base(:,2) ; rconv(:,2)=cell_base(:,1) ; rconv(:,1)=cell_base(:,3) 130 else 131 rconv(:,:)=cell_base(:,:) 132 end if 133 134 ! Compute the length of the three conventional vectors 135 length(1)=sqrt(sum(rconv(:,1)**2)) 136 length(2)=sqrt(sum(rconv(:,2)**2)) 137 length(3)=sqrt(sum(rconv(:,3)**2)) 138 139 ! Take care of the first conventional vector aligned with rbasis(:,3) (or aligned with the trigonal axis if rhombohedral) 140 ! and choice of the first normalized direction 141 if(iholohedry==5)then 142 rbasis(:,3)=third*(rconv(:,1)+rconv(:,2)+rconv(:,3)) 143 else 144 rbasis(:,3)=rconv(:,3) 145 end if 146 norm(3)=sqrt(sum(rbasis(:,3)**2)) 147 rnormalized(:,3)=rbasis(:,3)/norm(3) 148 149 ! Projection of the first conventional vector perpendicular to rbasis(:,3) 150 ! and choice of the first normalized direction 151 scprod1=sum(rnormalized(:,3)*cell_base(:,1)) 152 rbasis(:,1)=rconv(:,1)-rnormalized(:,3)*scprod1 153 norm(1)=sqrt(sum(rbasis(:,1)**2)) 154 rnormalized(:,1)=rbasis(:,1)/norm(1) 155 156 ! Generation of the second vector, perpendicular to the third and first 157 rnormalized(1,2)=rnormalized(2,3)*rnormalized(3,1)-rnormalized(3,3)*rnormalized(2,1) 158 rnormalized(2,2)=rnormalized(3,3)*rnormalized(1,1)-rnormalized(1,3)*rnormalized(3,1) 159 rnormalized(3,2)=rnormalized(1,3)*rnormalized(2,1)-rnormalized(2,3)*rnormalized(1,1) 160 161 ! Compute the vectors of the conventional cell, on the basis of iholohedry 162 if(iholohedry==2)then 163 rconv_new(:,3)=rconv(:,3) 164 rconv_new(:,1)=rconv(:,1) 165 rconv_new(:,2)=rnormalized(:,2)*length(2) ! Now, the y axis is perpendicular to the two others, that have not been changed 166 else if(iholohedry==3.or.iholohedry==4.or.iholohedry==7)then 167 if(iholohedry==7)then 168 symmetrized_length(1:3)=sum(length(:))*third 169 else if(iholohedry==4)then 170 symmetrized_length(3)=length(3) 171 symmetrized_length(1:2)=half*(length(1)+length(2)) 172 else if(iholohedry==3)then 173 symmetrized_length(:)=length(:) 174 end if 175 do ii=1,3 176 rconv_new(:,ii)=rnormalized(:,ii)*symmetrized_length(ii) 177 end do 178 else if(iholohedry==5)then 179 ! In the normalized basis, they have coordinates (a,0,c), and (-a/2,+-sqrt(3)/2*a,c) 180 ! c is known, but a is computed from the knowledge of the average length of the initial vectors 181 aa=sqrt(sum(length(:)**2)*third-norm(3)**2) 182 rconv_new(:,1)=aa*rnormalized(:,1)+rbasis(:,3) 183 rconv_new(:,2)=aa*half*(-rnormalized(:,1)+sqrt(three)*rnormalized(:,2))+rbasis(:,3) 184 rconv_new(:,3)=aa*half*(-rnormalized(:,1)-sqrt(three)*rnormalized(:,2))+rbasis(:,3) 185 else if(iholohedry==6)then 186 187 ! In the normalized basis, they have coordinates (a,0,0), (-a/2,+-sqrt(3)/2*a,0), and (0,0,c) 188 ! c is known, but a is computed from the knowledge of the average length of the initial vectors 189 aa=half*(length(1)+length(2)) 190 rconv_new(:,1)=aa*rnormalized(:,1) 191 rconv_new(:,2)=aa*half*(-rnormalized(:,1)+sqrt(three)*rnormalized(:,2)) 192 rconv_new(:,3)=rconv(:,3) 193 end if 194 195 ! Check whether the modification make sense 196 do ii=1,3 197 do jj=1,3 198 reldiff=(rconv_new(ii,jj)-rconv(ii,jj))/length(jj) 199 ! Allow for twice tolsym 200 if(abs(reldiff)>two*tolsym)then 201 write(msg,'(a,6(2a,3es14.6))')& 202 & 'Failed rectification of lattice vectors to comply with Bravais lattice identification, modifs are too large',ch10,& 203 & ' rconv =',rconv(:,1),ch10,& 204 & ' ',rconv(:,2),ch10,& 205 & ' ',rconv(:,3),ch10,& 206 & ' rconv_new=',rconv_new(:,1),ch10,& 207 & ' ',rconv_new(:,2),ch10,& 208 & ' ',rconv_new(:,3) 209 MSG_ERROR_CLASS(msg, "TolSymError") 210 end if 211 end do 212 end do 213 214 ! Copy back the cell_base vectors 215 if(iholohedry==4.and.equal(1)==1)then 216 cell_base(:,3)=rconv_new(:,2) ; cell_base(:,2)=rconv_new(:,1) ; cell_base(:,1)=rconv_new(:,3) 217 else if (iholohedry==4.and.equal(2)==1)then 218 cell_base(:,3)=rconv_new(:,1) ; cell_base(:,1)=rconv_new(:,2) ; cell_base(:,2)=rconv_new(:,3) 219 else 220 cell_base(:,:)=rconv_new(:,:) 221 end if 222 223 end if 224 225 end subroutine holocell