TABLE OF CONTENTS

ABINIT/dist2 [ Functions ]

[ Top ] [ Functions ]

NAME

  dist2

FUNCTION

  dist2(v1,v2,rprimd,option) calculates the distance of v1 and v2 in a crystal by
  repeating the unit cell

COPYRIGHT

  Copyright (C) 1998-2017 ABINIT group (DJA)
  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

  v1,v2
  rprimd: dimensions of the unit cell. if not given 1,0,0/0,1,0/0,0,1 is assumed
  option: 0 v1, v2 given in cartesian coordinates (default) / 1 v1,v2 given in reduced coordinates

OUTPUT

  dist2

PARENTS

  ioniondist

CHILDREN

SOURCE

 33 #if defined HAVE_CONFIG_H
 34 #include "config.h"
 35 #endif
 36 
 37 #include "abi_common.h"
 38 
 39 
 40 function dist2(v1,v2,rprimd,option)
 41 
 42  use defs_basis
 43 
 44  use m_numeric_tools,   only : wrap2_pmhalf
 45 
 46 !This section has been created automatically by the script Abilint (TD).
 47 !Do not modify the following lines by hand.
 48 #undef ABI_FUNC
 49 #define ABI_FUNC 'dist2'
 50  use interfaces_41_geometry, except_this_one => dist2
 51 !End of the abilint section
 52 
 53  implicit none
 54 
 55 !Arguments ------------------------------------
 56 !scalars
 57  integer,intent(in),optional                      :: option
 58  real(dp)                                         :: dist2
 59 !arrays
 60  real(dp),intent(in),optional      :: rprimd(3,3)
 61  real(dp),intent(in)                  :: v1(3),v2(3)
 62 
 63 !Local variables-------------------------------
 64 !scalars
 65  integer :: i1,i2,i3,opt,s1,s2,s3
 66  real(dp):: min2,norm2,ucvol
 67 !arrays
 68  integer :: limits(3)
 69  real(dp) :: corner(3),dred(3),dtot(3),dv(3),dwrap(3),sh(3)
 70  real(dp) :: gmet(3,3),gprimd(3,3),rmet(3,3)
 71  real(dp) :: vprimd(3,3)
 72 
 73 ! *************************************************************************
 74 
 75  if (.not.PRESENT(rprimd)) then
 76    vprimd=reshape((/1,0,0,  0,1,0,  0,0,1/),(/3,3/))
 77  else
 78    vprimd=rprimd
 79  end if
 80 
 81  call metric(gmet,gprimd,-1,rmet,vprimd,ucvol)
 82 
 83  dv(:)=v2(:)-v1(:)
 84 
 85 !If in cartesian coordinates, need to be transformed to reduced coordinates.
 86  opt=0
 87  if(present(option))then
 88    opt=option
 89  end if
 90  if(opt==0)then
 91    dred(:)=gprimd(1,:)*dv(1)+gprimd(2,:)*dv(2)+gprimd(3,:)*dv(3)
 92  else
 93    dred(:)=dv(:)
 94  end if
 95 
 96 !Wrap in the ]-1/2,1/2] interval
 97  call wrap2_pmhalf(dred(1),dwrap(1),sh(1))
 98  call wrap2_pmhalf(dred(2),dwrap(2),sh(2))
 99  call wrap2_pmhalf(dred(3),dwrap(3),sh(3))
100 
101 !Compute the limits of the parallelipipedic box that contains the Wigner-Seitz cell
102 !The reduced coordinates of the corners of the Wigner-Seitz cell are computed (multiplied by two)
103 !Then, the maximal values of these reduced coordinates are stored.
104  limits(:)=0
105  do s1=-1,1,2
106    do s2=-1,1,2
107      do s3=-1,1,2
108        corner(:)=gmet(:,1)*s1*rmet(1,1)+gmet(:,2)*s2*rmet(2,2)+gmet(:,3)*s3*rmet(3,3) 
109        limits(1)=max(limits(1),ceiling(abs(corner(1))+tol14))
110        limits(2)=max(limits(2),ceiling(abs(corner(2))+tol14))
111        limits(3)=max(limits(3),ceiling(abs(corner(3))+tol14))
112      end do
113    end do
114  end do
115 
116 !Use all relevant primitive real space lattice vectors to find the minimal difference vector
117  min2=huge(zero)
118  do i1=-limits(1),limits(1)
119    do i2=-limits(2),limits(2)
120      do i3=-limits(3),limits(3)
121        dtot(1)=dwrap(1)+i1
122        dtot(2)=dwrap(2)+i2
123        dtot(3)=dwrap(3)+i3
124        norm2=dtot(1)*rmet(1,1)*dtot(1)+dtot(2)*rmet(2,2)*dtot(2)+dtot(3)*rmet(3,3)*dtot(3)+&
125 &       2*(dtot(1)*rmet(1,2)*dtot(2)+dtot(2)*rmet(2,3)*dtot(3)+dtot(3)*rmet(3,1)*dtot(1))
126        min2=min(norm2,min2)
127      end do
128    end do
129  end do
130  dist2=sqrt(min2)
131 
132  end function dist2