Branch data Line data Source code
1 : : !!****m* ABINIT/m_tdep_phdos
2 : : !! NAME
3 : : !! m_tdep_phdos
4 : : !!
5 : : !! FUNCTION
6 : : !! This module contains the routine computing the phonon density of states
7 : : !! and the phonon frequencies and eigenvectors.
8 : : !!
9 : : !! COPYRIGHT
10 : : !! Copyright (C) 2011-2026 ABINIT group (GA,FB,JB)
11 : : !! This file is distributed under the terms of the
12 : : !! GNU General Public License, see ~abinit/COPYING
13 : : !! or http://www.gnu.org/copyleft/gpl.txt .
14 : : !!
15 : : !! SOURCE
16 : :
17 : : #if defined HAVE_CONFIG_H
18 : : #include "config.h"
19 : : #endif
20 : :
21 : : #include "abi_common.h"
22 : :
23 : : module m_tdep_phdos
24 : :
25 : : use defs_basis
26 : : !FB use m_nctk
27 : : use m_errors
28 : : use m_abicore
29 : : use m_phonons
30 : : use m_xmpi
31 : : use m_io_tools
32 : : use m_ifc, only : ifc_type
33 : : use m_crystal, only : crystal_t
34 : : use m_tdep_phi2, only : Eigen_type, tdep_write_yaml, tdep_write_dij
35 : : use m_tdep_qpt, only : Qpoints_type, Qbz_type
36 : : use m_tdep_dataset, only : atdep_dataset_type, MPI_enreg_type
37 : :
38 : : implicit none
39 : : !!***
40 : :
41 : : public :: tdep_calc_phdos
42 : :
43 : : contains
44 : :
45 : : !=====================================================================================================
46 : :
47 : : !!****f* ABINIT/m_tdep_phdos/tdep_calc_phdos
48 : : !! NAME
49 : : !! tdep_calc_phdos
50 : : !!
51 : : !! FUNCTION
52 : : !! Compute the phonon density of states (on the q-point grid Invar%ngqpt2,
53 : : !! with a gaussian smearing) from the IFC, and write it in the _PHDOS file
54 : : !! and, normalized to one, in the _vdos.dat file. The frequency mesh is
55 : : !! enlarged until it contains all the phonon frequencies.
56 : : !! Then, compute the phonon frequencies, eigenvectors and dynamical matrices
57 : : !! along the q-point path (written in the _omega.dat, _eigenvectors.dat,
58 : : !! _dij.dat and YAML files) and on the irreducible q-points of Qbz.
59 : : !!
60 : : !! INPUTS
61 : : !! Invar = Input object containing the input variables.
62 : : !! MPIdata = Info on MPI parallelism.
63 : : !! Ifc = IFC object.
64 : : !! Qbz = Qbz object containing the q-points of the irreducible Brillouin zone.
65 : : !! Qpt = Qpoints object containing the q-point path.
66 : : !!
67 : : !! OUTPUT
68 : : !! PHdos = Phonon density of states (normalized to one).
69 : : !!
70 : : !! SIDE EFFECTS
71 : : !! Crystal = Crystal object describing the unitcell.
72 : : !! Eigen2nd_MP = Eigen object. On exit, contains the phonon frequencies,
73 : : !! eigenvectors and dynamical matrices on the irreducible
74 : : !! q-points of Qbz.
75 : : !! Eigen2nd_path = Eigen object. On exit, contains the phonon frequencies,
76 : : !! eigenvectors and dynamical matrices along the q-point path.
77 : : !!
78 : : !! SOURCE
79 : :
80 [ - + - + : 644 : subroutine tdep_calc_phdos(Invar,MPIdata,Crystal,Ifc,Qbz,Qpt,Eigen2nd_MP,Eigen2nd_path,PHdos)
- + - + -
+ - + - +
- + - + -
+ - + + +
+ + ]
81 : :
82 : : type(atdep_dataset_type),intent(in) :: Invar
83 : : type(phdos_t),intent(out) :: PHdos
84 : : type(ifc_type),intent(in) :: Ifc
85 : : type(crystal_t),intent(inout) :: Crystal
86 : : type(Qbz_type),intent(in) :: Qbz
87 : : type(Qpoints_type),intent(in) :: Qpt
88 : : type(MPI_enreg_type), intent(in) :: MPIdata
89 : : type(Eigen_type),intent(inout) :: Eigen2nd_path
90 : : type(Eigen_type),intent(inout) :: Eigen2nd_MP
91 : :
92 : : integer :: natom_unitcell
93 : : integer :: prtdos,iqpt,iq_ibz,iomega
94 : : integer :: dos_ngqpt(3)
95 : : integer :: count_wminmax(2)
96 : : character (len=fnlen):: phdos_fname
97 : : double precision :: dossmear,integ,domega
98 : : double precision :: dos_qshift(3)
99 [ - + ]: 46 : double precision, allocatable :: displ(:,:)
100 : : character(len=500) :: message
101 : : real(dp) :: wminmax(2)
102 : :
103 : 46 : write(ab_out,*)' '
104 : 46 : write(ab_out,*) '#############################################################################'
105 : 46 : write(ab_out,*) '############## Compute the phonon spectrum and DOS ##########################'
106 : 46 : write(ab_out,*) '#############################################################################'
107 : :
108 : 46 : natom_unitcell = Invar%natom_unitcell
109 : :
110 : : ! Compute the DOS
111 : : ! ===============
112 : 46 : prtdos=1 !Gaussian
113 : : ! prtdos=2 !Tetra
114 : 46 : dossmear=4.5d-6
115 : : ! dossmear=4.5d-5
116 [ + + ]: 184 : dos_qshift(:)=0.5d0
117 [ + + ]: 184 : dos_ngqpt(:)=Invar%ngqpt2(:)
118 : 46 : write(ab_out,'(a)') ' Compute the phonon DOS'
119 : : ! Only 1 shift in q-mesh
120 : 46 : wminmax = zero
121 : 46 : do
122 : : call PHdos%init(Crystal,Ifc,prtdos,Invar%dosdeltae,dossmear,dos_ngqpt,1,dos_qshift, &
123 : 92 : "freq_displ", wminmax, count_wminmax, XMPI_WORLD,prtout=.false.)
124 [ + + + + : 230 : if (all(count_wminmax == 0)) exit
+ + ]
125 : 46 : wminmax(1) = wminmax(1) - abs(wminmax(1)) * 0.05
126 : 46 : wminmax(2) = wminmax(2) + abs(wminmax(2)) * 0.05
127 : 46 : call phdos%free()
128 : 46 : write(message, "(a, 2f8.5)")" Initial frequency mesh not large enough. Recomputing PHDOS with wmin, wmax: ",wminmax
129 : 46 : call wrtout(std_out, message)
130 : : end do
131 : :
132 : : ! Print the DOS
133 : : ! =============
134 [ + - ]: 46 : phdos_fname = trim(Invar%output_prefix)//"_PHDOS"
135 [ + - ]: 46 : if (MPIdata%iam_master) then
136 : 46 : write(ab_out,'(a)') ' Writing phonon DOS in file PHDOS'
137 : 46 : call phdos%print(phdos_fname)
138 : : end if
139 : 46 : domega=(Invar%dosdeltae*Ha_meV)
140 : 46 : integ=0.d0
141 [ + + ]: 84990 : do iomega=1,PHdos%nomega
142 : 84990 : integ=integ + domega*PHdos%phdos(iomega)
143 : : end do
144 [ + + ]: 84990 : PHdos%phdos(:)=PHdos%phdos(:)/integ
145 [ + - ]: 46 : if (MPIdata%iam_master) then
146 : 46 : write(ab_out,'(a)') ' Writing vibrational DOS in file vdos.dat'
147 : 46 : open(unit=56,file=trim(Invar%output_prefix)//'_vdos.dat')
148 [ + + ]: 84990 : do iomega=1,PHdos%nomega
149 [ + + ]: 84944 : if (Invar%enunit.eq.0) write(56,'(2(f18.6,1x))') PHdos%omega(iomega)*Ha_eV*1000,PHdos%phdos(iomega)
150 [ + + ]: 84944 : if (Invar%enunit.eq.1) write(56,'(2(f18.6,1x))') PHdos%omega(iomega)*Ha_cmm1 ,PHdos%phdos(iomega)
151 [ + + ]: 84944 : if (Invar%enunit.eq.2) write(56,'(2(f18.6,1x))') PHdos%omega(iomega)*1000 ,PHdos%phdos(iomega)
152 [ + + ]: 84990 : if (Invar%enunit.eq.3) write(56,'(2(f18.6,1x))') PHdos%omega(iomega)*Ha_THz ,PHdos%phdos(iomega)
153 : : end do
154 : 46 : close(56)
155 : : end if
156 : :
157 : : ! Compute the frequencies (2)
158 : : ! =======================
159 [ + - - + : 1372398 : ABI_MALLOC(displ,(2*3*natom_unitcell*3*natom_unitcell,Qpt%nqpt)); displ(:,:)=zero
+ - - + +
- + - + -
- + - + +
+ + + ]
160 [ + + + + : 4783736 : Eigen2nd_path%eigenval=zero ; Eigen2nd_path%eigenvec=zero ; Eigen2nd_path%dynmat=zero
+ + + + +
+ + + + +
+ + + + +
+ + + + +
+ + + + ]
161 [ + + ]: 9336 : do iqpt=1,Qpt%nqpt
162 : : call Ifc%fourq(Crystal,Qpt%qpt_red(:,iqpt),Eigen2nd_path%eigenval(:,iqpt),displ(:,iqpt),&
163 : : & out_eigvec=Eigen2nd_path%eigenvec(:,:,:,:,:,iqpt),&
164 : 9336 : & out_d2cart=Eigen2nd_path%dynmat (:,:,:,:,:,iqpt))
165 : : end do
166 [ - + ]: 46 : ABI_FREE(displ)
167 : :
168 : : ! Write the Dij, eigenvalues and eigenvectors (in ASCII and YAML)
169 : : ! ======================================================
170 [ + - ]: 46 : if (MPIdata%iam_master) then
171 : 46 : write(ab_out,'(a)') ' Writing dynamical matrix in file dij.dat'
172 : 46 : write(ab_out,'(a)') ' Writing phonon frequencies in file omega.dat'
173 : 46 : write(ab_out,'(a)') ' Writing phonon eigenvectors in file eigenvectors.dat'
174 : 46 : open(unit=51,file=trim(Invar%output_prefix)//'_eigenvectors.dat')
175 : 46 : open(unit=52,file=trim(Invar%output_prefix)//'_dij.dat')
176 : 46 : open(unit=53,file=trim(Invar%output_prefix)//'_omega.dat')
177 [ + + ]: 46 : if (Invar%enunit.eq.0) write(53,'(a)') '# Phonon frequencies in meV'
178 [ + + ]: 46 : if (Invar%enunit.eq.1) write(53,'(a)') '# Phonon frequencies in cm-1'
179 [ + + ]: 46 : if (Invar%enunit.eq.2) write(53,'(a)') '# Phonon frequencies in mHa'
180 [ + + ]: 46 : if (Invar%enunit.eq.3) write(53,'(a)') '# Phonon frequencies in THz'
181 [ + + ]: 9336 : do iqpt=1,Qpt%nqpt
182 : 9336 : call tdep_write_dij(Eigen2nd_path,iqpt,Invar,Qpt%qpt_red(:,iqpt))
183 : : end do
184 : 46 : close(51)
185 : 46 : close(52)
186 : 46 : close(53)
187 : 46 : call tdep_write_yaml(Eigen2nd_path,Qpt,Invar%output_prefix)
188 : : end if
189 : :
190 : : ! Compute the frequencies (2)
191 : : ! =======================
192 [ - + + - : 2121588 : ABI_MALLOC(displ,(2*3*natom_unitcell*3*natom_unitcell,Qbz%nqibz)); displ(:,:)=zero
- + + - +
- + - - +
- + + + +
+ ]
193 [ + + + + : 14564642 : Eigen2nd_MP%eigenval=zero ; Eigen2nd_MP%eigenvec=zero ; Eigen2nd_MP%dynmat=zero
+ + + + +
+ + + + +
+ + + + +
+ + + + +
+ + + + ]
194 [ + + ]: 7242 : do iq_ibz=1,Qbz%nqibz
195 : : call Ifc%fourq(Crystal,Qbz%qibz(:,iq_ibz),Eigen2nd_MP%eigenval(:,iq_ibz),displ(:,iq_ibz),&
196 : : & out_eigvec=Eigen2nd_MP%eigenvec(:,:,:,:,:,iq_ibz),&
197 : 7242 : & out_d2cart=Eigen2nd_MP%dynmat (:,:,:,:,:,iq_ibz))
198 : : end do
199 [ - + ]: 46 : ABI_FREE(displ)
200 : :
201 : 46 : end subroutine tdep_calc_phdos
202 : : !!***
203 : :
204 : : !=====================================================================================================
205 : :
206 : : end module m_tdep_phdos
207 : : !!***
|
