TABLE OF CONTENTS


ABINIT/fred2fcart [ Functions ]

[ Top ] [ Functions ]

NAME

 fred2fcart

FUNCTION

 Convert reduced forces into cartesian forces

COPYRIGHT

 Copyright (C) 1998-2017 ABINIT group (DCA, XG, GMR, FJ, 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/doc/developers/contributors.txt.

INPUTS

  fred(3,natom)=symmetrized grtn = d(etotal)/d(xred)
  natom=Number of atoms in the unitary cell
  Favgz_null=TRUE if the average cartesian force has to be set to zero
             FALSE if it is set to zero only in x,y directions (not z)
  gprimd(3,3)=dimensional primitive translations for reciprocal space(bohr^-1)

OUTPUT

  fcart(3,natom)=forces in cartesian coordinates (Ha/Bohr)

SIDE EFFECTS

NOTES

    Unlike fred, fcart has been corrected by enforcing
    the translational symmetry, namely that the sum of force
    on all atoms is zero (except is a slab is used)

PARENTS

      forces,m_mep

CHILDREN

SOURCE

41 #if defined HAVE_CONFIG_H
42 #include "config.h"
43 #endif
44 
45 #include "abi_common.h"
46 
47 
48 subroutine fred2fcart(favg,Favgz_null,fcart,fred,gprimd,natom)
49 
50  use defs_basis
51  use m_profiling_abi
52 
53 !This section has been created automatically by the script Abilint (TD).
54 !Do not modify the following lines by hand.
55 #undef ABI_FUNC
56 #define ABI_FUNC 'fred2fcart'
57 !End of the abilint section
58 
59  implicit none
60 
61 !Arguments ------------------------------------
62 !scalars
63  integer,intent(in) :: natom
64  logical :: Favgz_null
65 !arrays
66  real(dp),intent(out) :: fcart(3,natom)
67  real(dp),intent(in) :: fred(3,natom)
68  real(dp),intent(in) :: gprimd(3,3)
69  real(dp),intent(out) :: favg(3)
70 
71 !Local variables-------------------------------
72 !scalars
73  integer :: iatom,mu
74 !arrays
75 
76 ! *************************************************************************
77 
78 !Note conversion to cartesian coordinates (bohr) AND
79 !negation to make a force out of a gradient
80  favg(:)=zero
81  do iatom=1,natom
82    do mu=1,3
83      fcart(mu,iatom)= - (gprimd(mu,1)*fred(1,iatom)+&
84 &     gprimd(mu,2)*fred(2,iatom)+&
85 &     gprimd(mu,3)*fred(3,iatom))
86      favg(mu)=favg(mu)+fcart(mu,iatom)
87    end do
88  end do
89 
90 !Subtract off average force from each force component
91 !to avoid spurious drifting of atoms across cell.
92  favg(:)=favg(:)/dble(natom)
93  if(.not.Favgz_null) favg(3)=zero
94  do iatom=1,natom
95    fcart(:,iatom)=fcart(:,iatom)-favg(:)
96  end do
97 
98 end subroutine fred2fcart