TABLE OF CONTENTS


ABINIT/dbeta [ Functions ]

[ Top ] [ Functions ]

NAME

 dbeta

FUNCTION

  Calculate the rotation matrix d^l_{m{\prim}m}(beta) using Eq. 4.14 of
  M.E. Rose, Elementary Theory of Angular Momentum,
             John Wiley & Sons, New-York, 1957

COPYRIGHT

 Copyright (C) 1998-2017 ABINIT group (FJ, MT, NH)
 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

  cosbeta= cosinus of beta (=Euler angle)
  ll= index l
  mm= index m
  mp= index m_prime

OUTPUT

  dbeta= rotation matrix

NOTES

  - This file comes from the file crystal_symmetry.f
    by N.A.W. Holzwarth and A. Tackett for the code pwpaw
  - Assume l relatively small so that factorials do not cause
    roundoff error

PARENTS

     setsymrhoij

CHILDREN

SOURCE

 40 #if defined HAVE_CONFIG_H
 41 #include "config.h"
 42 #endif
 43 
 44 #include "abi_common.h"
 45 
 46 
 47 function dbeta(cosbeta,ll,mp,mm)
 48 
 49  use defs_basis
 50 
 51  use m_special_funcs, only : factorial
 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 'dbeta'
 57 !End of the abilint section
 58 
 59  implicit none
 60 
 61 !Arguments ---------------------------------------------
 62 !scalars
 63  integer,intent(in) :: ll,mm,mp
 64  real(dp) :: dbeta
 65  real(dp),intent(in) :: cosbeta
 66 
 67 !Local variables ------------------------------
 68 !scalars
 69  integer,parameter :: mxterms=200
 70  integer :: ii,ina,inb,inc,ml,ms
 71  real(dp) :: arg,cosbetab2,pref,sinbetab2,sum,tt
 72 !************************************************************************
 73  dbeta=zero
 74 
 75 !Special cases
 76  if (abs(cosbeta-1._dp).lt.tol10) then
 77    if (mp.eq.mm) dbeta=1
 78  else if (abs(cosbeta+1._dp).lt.tol10) then
 79    if (mp.eq.-mm) dbeta=(-1)**(ll+mm)
 80  else
 81 
 82 !  General case
 83    cosbetab2=sqrt((1+cosbeta)*0.5_dp)
 84    sinbetab2=sqrt((1-cosbeta)*0.5_dp)
 85    ml=max(mp,mm)
 86    ms=min(mp,mm)
 87    if (ml.ne.mp) sinbetab2=-sinbetab2
 88    tt=-(sinbetab2/cosbetab2)**2
 89    pref=sqrt((factorial(ll-ms)*factorial(ll+ml))&
 90 &   /(factorial(ll+ms)*factorial(ll-ml)))&
 91 &   /factorial(ml-ms)*(cosbetab2**(2*ll+ms-ml))&
 92 &   *((-sinbetab2)**(ml-ms))
 93    sum=1._dp
 94    arg=1._dp
 95    ina=ml-ll
 96    inb=-ms-ll
 97    inc=ml-ms+1
 98    do ii=1,mxterms
 99      if (ina.eq.0.or.inb.eq.0) exit
100      arg=(arg*ina*inb*tt)/(ii*inc)
101      sum=sum+arg
102      ina=ina+1
103      inb=inb+1
104      inc=inc+1
105    end do
106    dbeta=pref*sum
107  end if
108 
109 end function dbeta