TABLE OF CONTENTS


ABINIT/convert_notation [ Functions ]

[ Top ] [ Functions ]

NAME

  convert_notation

FUNCTION

  notation: Convert notation index of the elastic tensor to use the  
            second derivative of gylm
            for example :  to calculate d2(gylm)/d(eps_32)
                          - 
                          - voigt notation       => 32   
                          - normal notation      => 3 3 2 2
                          - notation for gylmgr2 => 32 32 32 32 => 4 4 4   
            The calculation of the second derivative of gylm wrt strains requires  four
            derivative (see the formula)

COPYRIGHT

  Copyright (C) 2013-2018 ABINIT group (AM)
  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

  eps_alpha 
  eps_beta  
  eps_delta 
  eps_gamma 

OUTPUT

  mu4(4) = array with index for the second derivative of gylm

SIDE EFFECTS

  return the four index in mu4 for the calculation of the second derivative of gylm

NOTES

PARENTS

      pawgrnl

CHILDREN

SOURCE

 45 #if defined HAVE_CONFIG_H
 46 #include "config.h"
 47 #endif
 48 
 49 #include "abi_common.h"
 50 
 51 subroutine convert_notation(mu4,eps_alpha,eps_beta,eps_gamma,eps_delta)
 52 
 53  use defs_basis
 54  use m_profiling_abi
 55  use m_errors
 56 
 57 !This section has been created automatically by the script Abilint (TD).
 58 !Do not modify the following lines by hand.
 59 #undef ABI_FUNC
 60 #define ABI_FUNC 'convert_notation'
 61 !End of the abilint section
 62 
 63  implicit none
 64 
 65 !Arguments ------------------------------------
 66  !scalar
 67  integer,intent(in)  :: eps_alpha,eps_beta
 68  integer,optional,intent(in)  :: eps_gamma,eps_delta
 69  !array
 70  integer,intent(inout) :: mu4(4)
 71 
 72 !Local variables-------------------------------
 73  integer :: eps1,eps2,i,j,k      
 74  integer,allocatable :: mu_temp(:)
 75 
 76 ! *************************************************************************
 77  
 78 !DBG_ENTER("COLL")
 79 
 80  ABI_ALLOCATE(mu_temp,(4))
 81  if (present(eps_gamma).and.present(eps_delta)) then
 82    mu_temp(1)=eps_alpha
 83    mu_temp(2)=eps_beta
 84    mu_temp(3)=eps_gamma
 85    mu_temp(4)=eps_delta
 86  else
 87    mu_temp(1)=eps_alpha
 88    mu_temp(2)=eps_beta
 89    mu_temp(3)= 0
 90    mu_temp(4)= 0
 91  end if
 92  k=1
 93  do i=1,2
 94    eps1=mu_temp(i)
 95    do j=1,2
 96      eps2=mu_temp(2+j)
 97      if(eps1==eps2) then
 98        if(eps1==1) mu4(k)=1;
 99        if(eps1==2) mu4(k)=2;
100        if(eps1==3) mu4(k)=3;
101      else
102        if((eps1==3.and.eps2==2).or.(eps1==2.and.eps2==3)) mu4(k)=4;
103        if((eps1==3.and.eps2==1).or.(eps1==1.and.eps2==3)) mu4(k)=5;
104        if((eps1==1.and.eps2==2).or.(eps1==2.and.eps2==1)) mu4(k)=6;
105      end if
106      k=k+1
107    end do
108  end do
109 
110  ABI_DEALLOCATE(mu_temp)
111 
112 !DBG_EXIT("COLL")
113 
114 end subroutine convert_notation