TABLE OF CONTENTS


ABINIT/metric [ Functions ]

[ Top ] [ Functions ]

FUNCTION

 Compute first dimensional primitive translation vectors in reciprocal space
 gprimd from rprimd, and eventually writes out.
 Then, computes metrics for real and recip space rmet and gmet using length
 dimensional primitive translation vectors in columns of rprimd(3,3) and gprimd(3,3).
  gprimd is the inverse transpose of rprimd.
  i.e. $ rmet_{i,j}= \sum_k ( rprimd_{k,i}*rprimd_{k,j} )  $
       $ gmet_{i,j}= \sum_k ( gprimd_{k,i}*gprimd_{k,j} )  $
 Also computes unit cell volume ucvol in $\textrm{bohr}^3$

COPYRIGHT

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

INPUTS

  rprimd(3,3)=dimensional primitive translations for real space (bohr)
  iout=unit number of output file.  If iout<0, do not write output.

OUTPUT

  gmet(3,3)=reciprocal space metric ($\textrm{bohr}^{-2}$).
  gprimd(3,3)=dimensional primitive translations for reciprocal space ($\textrm{bohr}^{-1}$)
  rmet(3,3)=real space metric ($\textrm{bohr}^{2}$).
  ucvol=unit cell volume ($\textrm{bohr}^{3}$).

PARENTS

      afterscfloop,bethe_salpeter,chkinp,clnup1,conducti_nc,conducti_paw
      conducti_paw_core,cut3d,d2frnl,dfpt_eltfrhar,dfpt_eltfrkin,dfpt_eltfrxc
      dfpt_looppert,dfpt_newvtr,dfpt_scfcv,dist2,elpolariz,emispec_paw,energy
      extrapwf,fftprof,finddistrproc,forces,forstr,get_npert_rbz,getkgrid
      ingeo,initaim,initberry,inkpts,inqpt,invacuum,invars2m,ks_ddiago
      linear_optics_paw,m_ab7_symmetry,m_crystal,m_cut3d,m_ddb,m_dens
      m_effective_potential,m_effective_potential_file,m_fft,m_fft_prof
      m_fit_polynomial_coeff,m_hamiltonian,m_io_kss,m_ioarr,m_mep,m_pawpwij
      m_screening,m_use_ga,m_vcoul,m_wfk,mag_constr,mag_constr_e,memory_eval
      mkcore_wvl,mlwfovlp_qp,moddiel,mpi_setup,mrgscr,newrho,newvtr,nres2vres
      odamix,optic,pawgrnl,prcref,prcref_PMA,pred_bfgs,pred_delocint
      pred_isothermal,pred_langevin,pred_lbfgs,pred_nose,pred_srkna14
      pred_verlet,prt_cif,prtefield,prtimg,psolver_rhohxc,rhohxc,scfcv
      screening,setup1,setup_bse,setup_screening,setup_sigma,sigma,smallprim
      stress,symmetrize_rprimd,testkgrid,thmeig,vdw_dftd2,vdw_dftd3
      wrt_moldyn_netcdf,wvl_initro,xchybrid_ncpp_cc,xfpack_vin2x,xfpack_x2vin

CHILDREN

      matr3inv,wrtout

SOURCE

 55 #if defined HAVE_CONFIG_H
 56 #include "config.h"
 57 #endif
 58 
 59 #include "abi_common.h"
 60 
 61 subroutine metric(gmet,gprimd,iout,rmet,rprimd,ucvol)
 62 
 63  use defs_basis
 64  use m_errors
 65  use m_profiling_abi
 66 
 67 !This section has been created automatically by the script Abilint (TD).
 68 !Do not modify the following lines by hand.
 69 #undef ABI_FUNC
 70 #define ABI_FUNC 'metric'
 71  use interfaces_14_hidewrite
 72  use interfaces_32_util
 73 !End of the abilint section
 74 
 75  implicit none
 76 
 77 !Arguments ------------------------------------
 78 !scalars
 79  integer,intent(in) :: iout
 80  real(dp),intent(out) :: ucvol
 81 !arrays
 82  real(dp),intent(in) :: rprimd(3,3)
 83  real(dp),intent(out) :: gmet(3,3),gprimd(3,3),rmet(3,3)
 84 
 85 !Local variables-------------------------------
 86 !scalars
 87  integer :: nu
 88  character(len=500) :: message
 89 !arrays
 90  real(dp) :: angle(3)
 91 
 92 ! *************************************************************************
 93 
 94 !Compute unit cell volume
 95  ucvol=rprimd(1,1)*(rprimd(2,2)*rprimd(3,3)-rprimd(3,2)*rprimd(2,3))+&
 96 & rprimd(2,1)*(rprimd(3,2)*rprimd(1,3)-rprimd(1,2)*rprimd(3,3))+&
 97 & rprimd(3,1)*(rprimd(1,2)*rprimd(2,3)-rprimd(2,2)*rprimd(1,3))
 98 
 99 !Check that the input primitive translations are not linearly dependent (and none is zero); i.e. ucvol~=0
100 !Also ask that the mixed product is positive.
101  if (abs(ucvol)<tol12) then 
102 !  write(std_out,*)"rprimd",rprimd,"ucvol",ucvol
103    write(message,'(5a)')&
104 &   'Input rprim and acell gives vanishing unit cell volume.',ch10,&
105 &   'This indicates linear dependency between primitive lattice vectors',ch10,&
106 &   'Action : correct either rprim or acell in input file.'
107    MSG_ERROR(message)
108  end if
109  if (ucvol<zero)then
110    write(message,'(2a,3(a,3es16.6,a),7a)')&
111 &   'Current rprimd gives negative (R1xR2).R3 . ',ch10,&
112 &   'Rprimd =',rprimd(:,1),ch10,&
113 &   '        ',rprimd(:,2),ch10,&
114 &   '        ',rprimd(:,3),ch10,&
115 &   'Action: if the cell size and shape are fixed (optcell==0),',ch10,&
116 &   '        exchange two of the input rprim vectors;',ch10,&
117 &   '        if you are optimizing the cell size and shape (optcell/=0),',ch10,&
118 &   '        maybe the move was too large, and you might try to decrease strprecon.'
119    MSG_ERROR(message)
120  end if
121 
122 !Generates gprimd
123  call matr3inv(rprimd,gprimd)
124 
125 !Write out rprimd, gprimd and ucvol
126  if (iout>=0) then
127    write(message,'(2a)')' Real(R)+Recip(G) ','space primitive vectors, cartesian coordinates (Bohr,Bohr^-1):'
128    call wrtout(iout,message,'COLL')
129    do nu=1,3
130      write(message, '(1x,a,i1,a,3f11.7,2x,a,i1,a,3f11.7)' ) &
131 &     'R(',nu,')=',rprimd(:,nu)+tol10,&
132 &     'G(',nu,')=',gprimd(:,nu)+tol10
133      call wrtout(iout,message,'COLL')
134    end do
135    write(message,'(a,1p,e15.7,a)') ' Unit cell volume ucvol=',ucvol+tol10,' bohr^3'
136    call wrtout(iout,message,'COLL')
137    call wrtout(std_out,message,'COLL')
138  end if
139 
140 !Compute real space metric.
141  rmet = MATMUL(TRANSPOSE(rprimd),rprimd)   
142 
143 !Compute reciprocal space metric.
144  gmet = MATMUL(TRANSPOSE(gprimd),gprimd)
145 
146 !Write out the angles
147  if (iout>=0) then
148    angle(1)=acos(rmet(2,3)/sqrt(rmet(2,2)*rmet(3,3)))/two_pi*360.0d0
149    angle(2)=acos(rmet(1,3)/sqrt(rmet(1,1)*rmet(3,3)))/two_pi*360.0d0
150    angle(3)=acos(rmet(1,2)/sqrt(rmet(1,1)*rmet(2,2)))/two_pi*360.0d0
151    write(message, '(a,3es16.8,a)' )' Angles (23,13,12)=',angle(1:3),' degrees'
152    call wrtout(iout,message,'COLL')
153    call wrtout(std_out,message,'COLL')
154  end if
155 
156 end subroutine metric