TABLE OF CONTENTS
- ABINIT/m_anharmonics_terms
- defs_abitypes/anharmonics_terms_type
- m_anharmonics_terms/anharmonics_terms_evaluateIFCStrainCoupling
- m_anharmonics_terms/anharmonics_terms_init
- m_effective_potential/anharmonics_terms_evaluateElastic
ABINIT/m_anharmonics_terms [ Functions ]
NAME
m_anharmonics_term
FUNCTION
Module with datatype and tools for the anharmonics terms
COPYRIGHT
Copyright (C) 2010-2024 ABINIT group (AM) This file is distributed under the terms of the GNU General Public Licence, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
SOURCE
19 #if defined HAVE_CONFIG_H 20 #include "config.h" 21 #endif 22 23 #include "abi_common.h" 24 25 module m_anharmonics_terms 26 27 use defs_basis 28 use m_errors 29 use m_abicore 30 use m_polynomial_coeff 31 use m_ifc, only : ifc_type 32 33 implicit none 34 35 public :: anharmonics_terms_init 36 public :: anharmonics_terms_free 37 public :: anharmonics_terms_freeCoeffs 38 public :: anharmonics_terms_evaluateElastic 39 public :: anharmonics_terms_evaluateIFCStrainCoupling 40 public :: anharmonics_terms_setCoeffs 41 public :: anharmonics_terms_setElastic3rd 42 public :: anharmonics_terms_setElastic4th 43 public :: anharmonics_terms_setElasticDispCoupling 44 public :: anharmonics_terms_setStrainPhononCoupling
defs_abitypes/anharmonics_terms_type [ Types ]
[ Top ] [ defs_abitypes ] [ Types ]
NAME
anharmonics_terms_type
FUNCTION
datatype for a effective potential constructed.
SOURCE
57 type, public :: anharmonics_terms_type 58 59 integer :: ncoeff = 0 60 ! nterm store the number of coefficients 61 62 logical :: has_elastic3rd 63 ! Flag to know if the 3rd derivatives with respect to strain is present 64 65 logical :: has_elastic4th 66 ! Flag to know if the 3rd derivatives with respect to strain is present 67 68 logical :: has_strain_coupling 69 ! Flag to know if the 3rd derivatives with respect to strain and 2 atom disp is present 70 71 logical :: has_elastic_displ 72 ! Flag to know if the 3rd derivatives with respect to 2 strain and 3 atom disp is present 73 74 logical :: bounded 75 ! True : the model is bounded 76 77 type(polynomial_coeff_type),dimension(:),allocatable :: coefficients 78 ! array with all the coefficients from polynomial coefficients 79 80 real(dp) :: elastic3rd(6,6,6) 81 ! elastic_constant(6,6,6) 82 ! Elastic tensor Hartree 83 84 real(dp) :: elastic4th(6,6,6,6) 85 ! elastic_constant(6,6,6) 86 ! Elastic tensor Hartree 87 88 real(dp), allocatable :: elastic_displacement(:,:,:,:) 89 ! elastic_displacement(6,6,3,natom) 90 ! internal strain tensor 91 92 type(ifc_type),dimension(:),allocatable :: phonon_strain 93 ! Array of ifc with phonon_strain coupling for each strain 94 95 end type anharmonics_terms_type
m_anharmonics_terms/anharmonics_terms_evaluateIFCStrainCoupling [ Functions ]
[ Top ] [ m_anharmonics_terms ] [ Functions ]
NAME
anharmonics_terms_evaluateIFCStrainCoupling
FUNCTION
This fonction compute the harmonic part of the energy of the supercell in the eff_pot
INPUTS
strain_phonon(6)<type(ifc_type) = strain-phonon coupling disp(3,natom_sc) = atomics displacement between configuration and the reference natom = number of atoms in the supercell natom_uc = number of atoms in the unit cell sc_size(3) = size of the supercell cells(ncell) = number of the cells into the supercell (1,2,3,4,5) ncell = total number of cell to treat index_cells(3,ncell) = indexes of the cells into supercell (-1 -1 -1 ,...,1 1 1) comm=MPI communicator
OUTPUT
energy = contribution of the ifc to the energy fcart(3,natom) = contribution of the ifc to the forces strten(6) = contribution to the stress tensor PARENT effective_potential_evaluate
SOURCE
723 subroutine anharmonics_terms_evaluateIFCStrainCoupling(phonon_strain,disp,energy,fcart,natom,natom_uc,& 724 & sc_size,strain,strten,cells,ncell,& 725 & index_cells,comm) 726 727 !Arguments ------------------------------- 728 ! scalars 729 real(dp),intent(out) :: energy 730 integer,intent(in) :: natom,natom_uc,ncell 731 integer,intent(in) :: comm 732 ! array 733 integer,intent(in) :: cells(ncell),index_cells(ncell,3) 734 integer,intent(in) :: sc_size(3) 735 type(ifc_type),intent(in) :: phonon_strain(6) 736 real(dp),intent(in) :: disp(3,natom) 737 real(dp),intent(out) :: fcart(3,natom) 738 real(dp),intent(out) :: strten(6) 739 real(dp),intent(in) :: strain(6) 740 !Local variables------------------------------- 741 ! scalar 742 integer :: alpha 743 integer :: i1,i2,i3,ia,ib,icell,ii 744 integer :: irpt,jj,kk,ll,mu,nu 745 integer :: ierr 746 real(dp):: ifc 747 ! array 748 integer :: cell_atom2(3) 749 character(500) :: msg 750 751 ! ************************************************************************* 752 753 if (any(sc_size <= 0)) then 754 write(msg,'(a,a)')' sc_size can not be inferior or equal to zero' 755 ABI_ERROR(msg) 756 end if 757 758 ! Initialisation of variables 759 energy = zero 760 fcart(:,:) = zero 761 strten(:) = zero 762 763 do icell = 1,ncell 764 ii = (cells(icell)-1)*natom_uc 765 i1=index_cells(icell,1); i2=index_cells(icell,2); i3=index_cells(icell,3) 766 do alpha=1,6 767 do irpt = 1,phonon_strain(alpha)%nrpt 768 ! get the cell of atom2 (0 0 0, 0 0 1...) 769 cell_atom2(1) = i1 + phonon_strain(alpha)%cell(1,irpt) 770 cell_atom2(2) = i2 + phonon_strain(alpha)%cell(2,irpt) 771 cell_atom2(3) = i3 + phonon_strain(alpha)%cell(3,irpt) 772 call getPBCIndexes_supercell(cell_atom2(1:3),sc_size(1:3)) 773 ! index of the second atom in the displacement array 774 jj = (cell_atom2(1)-1)*sc_size(2)*sc_size(3)*natom_uc+& 775 & (cell_atom2(2)-1)*sc_size(3)*natom_uc+& 776 & (cell_atom2(3)-1)*natom_uc 777 do ib = 1, natom_uc 778 ll = jj + ib 779 do nu=1,3 780 do ia = 1, natom_uc 781 kk = ii + ia 782 do mu=1,3 783 ifc = phonon_strain(alpha)%atmfrc(mu,ia,nu,ib,irpt) 784 ! accumule energy 785 energy = energy + sixth*strain(alpha)*disp(mu,kk)*disp(nu,ll)*ifc 786 ! accumule forces 787 fcart(mu,kk) = fcart(mu,kk) + half*strain(alpha)*disp(nu,ll)*ifc 788 ! accumule stresses 789 strten(alpha) = strten(alpha) + half*disp(mu,kk)*disp(nu,ll)*ifc 790 end do 791 end do 792 end do 793 end do 794 end do 795 end do 796 end do 797 798 ! MPI_SUM 799 call xmpi_sum(energy, comm, ierr) 800 call xmpi_sum(fcart , comm, ierr) 801 call xmpi_sum(strten, comm, ierr) 802 803 end subroutine anharmonics_terms_evaluateIFCStrainCoupling
m_anharmonics_terms/anharmonics_terms_init [ Functions ]
[ Top ] [ m_anharmonics_terms ] [ Functions ]
NAME
anharmonics_terms_init
FUNCTION
Initialize anharmonics_terms datatype
INPUTS
natom = number of atoms in primitive cell ncoeff = number of coefficient for the fited polynome bounded = optional, flag to now if the model in bounded elastic3rd(6,6,6) = optional,3rd order of the elastic constants elastic4th(6,6,6,6) = optional,4st order of the elastic constants elastic_displacement(6,6,3,natom) = optional,elastic constant - force coupling phonon_strain<type(ifc_type)>(6) = optional,phonon strain couling coeffs<type(polynomial_coeff_type)>(ncoeff) = optional,datatype with polynomial coefficients
OUTPUT
anharmonics_terms<type(anharmonics_terms_type)> = anharmonics_terms datatype to be initialized
SOURCE
124 subroutine anharmonics_terms_init(anharmonics_terms,natom,ncoeff,& 125 & bounded,elastic3rd,elastic4th,elastic_displacement,& 126 & phonon_strain,coeffs) 127 128 !Arguments ------------------------------------ 129 !scalars 130 integer, intent(in) :: natom,ncoeff 131 type(anharmonics_terms_type), intent(out) :: anharmonics_terms 132 real(dp),optional,intent(in) :: elastic_displacement(6,6,3,natom) 133 real(dp),optional,intent(in) :: elastic3rd(6,6,6),elastic4th(6,6,6,6) 134 type(polynomial_coeff_type),optional :: coeffs(ncoeff) 135 type(ifc_type),optional,intent(in) :: phonon_strain(6) 136 logical,optional,intent(in) :: bounded 137 !arrays 138 !Local variables------------------------------- 139 !scalar 140 !arrays 141 character(len=500) :: msg 142 143 ! ************************************************************************* 144 145 call anharmonics_terms_free(anharmonics_terms) 146 147 ! Check the number of atoms 148 if (natom < 1) then 149 write(msg, '(a,a,a,i10,a)' )& 150 & 'The cell must have at least one atom.',ch10,& 151 & 'The number of atom is ',natom,'.' 152 ABI_BUG(msg) 153 end if 154 155 !Allocation of phonon strain coupling array (3rd order) 156 if(present(phonon_strain)) then 157 call anharmonics_terms_setStrainPhononCoupling(anharmonics_terms,natom,phonon_strain) 158 end if 159 160 !Set the 3rd order elastic tensor 161 anharmonics_terms%elastic3rd = zero 162 anharmonics_terms%has_elastic3rd = .FALSE. 163 if(present(elastic3rd))then 164 call anharmonics_terms_setElastic3rd(anharmonics_terms,elastic3rd) 165 end if 166 167 !Set the 3rd order elastic tensor 168 anharmonics_terms%elastic4th = zero 169 anharmonics_terms%has_elastic4th = .FALSE. 170 if(present(elastic4th))then 171 call anharmonics_terms_setElastic4th(anharmonics_terms,elastic4th) 172 end if 173 174 !Allocation of 3rd order with respecto to 2 strain and 1 atomic displacement 175 if(present(elastic_displacement))then 176 call anharmonics_terms_setElasticDispCoupling(anharmonics_terms,natom,elastic_displacement) 177 end if 178 179 anharmonics_terms%ncoeff = 0 180 181 !Allocation of the coefficient 182 if(present(coeffs))then 183 if(ncoeff /= size(coeffs))then 184 write(msg, '(a)' )& 185 & ' ncoeff has not the same size than coeffs array, ' 186 ABI_BUG(msg) 187 end if 188 call anharmonics_terms_setCoeffs(coeffs,anharmonics_terms,ncoeff) 189 end if 190 191 !Set the flag bounded 192 if(present(bounded))then 193 anharmonics_terms%bounded = bounded 194 else 195 anharmonics_terms%bounded = .FALSE. 196 end if 197 198 end subroutine anharmonics_terms_init
m_effective_potential/anharmonics_terms_evaluateElastic [ Functions ]
[ Top ] [ m_effective_potential ] [ Functions ]
NAME
anharmonics_terms_evaluateElastic
FUNCTION
Compute the energy, stresses and forces related to the application of strain
INPUTS
disp(3,natom_sc) = atomics displacement between configuration and the reference natom = number of atom in the supercell natom_uc = number of atom in the unit cell ncell = number of cell strain(6) = strain to apply elastic3rd(6,6,6) = 3 order derivatives with respect to to 3 strain elastic4th(6,6,66,) = 4 order derivatives with respect to to 4 strain elastic_displacement(6,6,3,natom) = 3 order derivatives with respect to 2 strain and 1 Atom disp
OUTPUT
energy = contribution of the ifc to the energy fcart(3,natom) = contribution of the ifc to the forces strten(6) = contribution to the stress tensor
SOURCE
609 subroutine anharmonics_terms_evaluateElastic(disp,energy,fcart,natom,natom_uc,ncell,strten,strain,& 610 & elastic3rd,elastic4th,elastic_displacement) 611 612 real(dp),intent(out):: energy 613 integer, intent(in) :: natom,natom_uc,ncell 614 ! array 615 real(dp),optional,intent(in) :: elastic3rd(6,6,6),elastic4th(6,6,6,6) 616 real(dp),optional,intent(in) :: elastic_displacement(6,6,3,natom) 617 real(dp),intent(out):: strten(6) 618 real(dp),intent(out):: fcart(3,natom) 619 real(dp),intent(in) :: disp(3,natom) 620 real(dp),intent(in) :: strain(6) 621 622 !Local variables------------------------------- 623 ! scalar 624 integer :: ia,ii,mu,alpha,beta,gamma,delta,d1,d2 625 real(dp):: cijk 626 logical :: has_elastic3rd,has_elastic4th,has_elastic_displ 627 ! array 628 ! ************************************************************************* 629 630 !Reset output and flags 631 energy = zero 632 fcart = zero 633 strten = zero 634 has_elastic3rd = .FALSE. 635 has_elastic4th = .FALSE. 636 has_elastic_displ = .FALSE. 637 d1=0;d2=0 638 639 !Set the flags 640 if(present(elastic3rd)) has_elastic3rd = .TRUE. 641 if(present(elastic4th)) then 642 has_elastic4th = .TRUE. 643 d1=1;d2=6 644 end if 645 if(present(elastic_displacement)) has_elastic_displ = .TRUE. 646 647 !1-Treat 3rd order elastic constants 648 if (has_elastic3rd.or.has_elastic4th) then 649 do alpha=1,6 650 do beta=1,6 651 do gamma=1,6 652 cijk = ncell*elastic3rd(alpha,beta,gamma) 653 ! Accumulate energy 654 energy = energy + sixth*cijk*strain(alpha)*strain(beta)*strain(gamma) 655 ! Accumulate stresses contributions 656 strten(alpha)=strten(alpha)+ half*cijk*strain(beta)*strain(gamma) 657 do delta=d1,d2 658 cijk = ncell*elastic4th(alpha,beta,gamma,delta) 659 ! Accumulate energy 660 energy = energy + (1/24.)*cijk*strain(alpha)*strain(beta)*& 661 & strain(gamma)*strain(delta) 662 ! Accumulate stresses contributions 663 strten(alpha)=strten(alpha)+ sixth*cijk*strain(beta)*strain(gamma)*& 664 & strain(delta) 665 end do 666 end do 667 end do 668 end do 669 end if 670 671 !2-Part due to the internat strain 672 if(has_elastic_displ)then 673 ii = 1 674 do ia = 1,natom 675 do mu = 1,3 676 do beta=1,6 677 do alpha=1,6 678 cijk = elastic_displacement(alpha,beta,mu,ii) 679 ! Accumulte for this atom 680 energy = energy + sixth*cijk*strain(alpha)*strain(beta)*disp(mu,ia) 681 fcart(mu,ia) = fcart(mu,ia) + half*cijk*strain(alpha)*strain(beta) 682 strten(alpha) = strten(alpha) + half*cijk*strain(beta)*disp(mu,ia) 683 end do 684 end do 685 end do 686 ii = ii +1 687 ! Reset to 1 if the number of atoms is superior than in the initial cell 688 if(ii==natom_uc+1) ii = 1 689 end do 690 end if 691 692 end subroutine anharmonics_terms_evaluateElastic