The ABINIT Users Mailing List ( CLOSED )

["Matteo Giantomassi" <Matteo.Giantomassi@uclouvain.be>]

>> You have to do this by hand for the moment. The idea is that for each
>> mode
>> you calculate the contribution in the phonon eigenvectors (normalized)
>> which
>> is due to a given atom, e.g.
>>
>> weight(lambda, q, tau) = sum_i  |u_lambda,q (i, tau)|^2
>>
>> for direction i, atom tau, phonon mode lambda and qpoint q. The sum of
>> these
>> over all tau gives 1 (for a given lambda,q), so you have a weight.
>> Multiply
>> this by the phonon dos, and you get a PDOS. This also works fine for
>> crystals and any other system.
>>
>> PDOS_tau (omega) = sum_q,lambda  weight(lambda, q, tau)
>> delta(omega-omega_q,lambda)
>>
>> It should be fairly easy to add this to the anaddb output at some point
>> -
>> someone may even have already done it! I notice Matteo Giantomassi has a
>> routine mkphdos with similar functionality: Matteo does this work yet?
>
>  Hi Matthieu,
>
>  Well, it should work but some hacking is needed to enter this part of the
> code (see below)
>  I wrote the routine some time ago and I've used it for one of my
>  applications. It provides the PHDOS calculated either with tetrahedrons
> or
>  standard gaussian smearing. Type-projected contributions to the total
> PHDOS
>  are also available.
>  Atom-projected contributions are not directly available because the
> different
>  contributions arising from atoms of the same type are summed up.
>  Anyway calculating atom-projected contributions is just a matter of
> changing a bit the routine (17ddb/mkphdos.F90).

  I had a look at my implementation and it seems that I took into
  account the possibility of calculating atom-projected PHDOS
   The array pjdos is indeed defined as

 pjdos(nomega,3,natom)

 so it contains the partial projected PHDOS for each atom and for
 each cartesian direction. One can easily extract from this array
 the contribution to the PHDOS given by a particular atom.


>
>  No example and reference output files are provided yet (I know, I know,
> shame on me!).
>  Digging into my files I found an anaddb input that can be used to
> calculate the PHDOS,
>  The flags that have to be added in the input file are:
>
> ifcflag  1  #IFC's are needed
>
> #DOS related variables
> prtdos  1              # 1 for gaussian, 2 for tetrahedrons
> ng2qpt  35 35 7        # Dense q-mesh for Fourier interpolation
> q2shft  0.0 0.0 0.0    # Shift for the dense mesh
> dosdeltae  9.112670E-7 # 0.1 cm-1  frequency step for PHDOS
> dosmear    1.366900E-5 # 3 cm-1    broadening (only for Gaussian method)
>
> At the end of the calculation the results are stored in the PHDOS file
> There are several columns:
>
>  # omega     PHDOS    IPHDOS   PJDOS[1]  IPJDOS[1] ..
>
>  PHDOS refers to the total phonon DOS, IPHDOS is the integrated DOS, then
> we have the
>  contribution due to all the atoms of the first type as well as the
> corresponding integrated dos.
>
> 3 Comments:
>
> 1) To call mkphdos you have to modify a bit 17ddb/invars9
>   commenting out the call to leave_new in the following piece of code
>
>  if(anaddb_dtset%prtdos/=0 .and. anaddb_dtset%ifcflag/=1) then
>   write(message, '(5a)' )&
> &  ' invars9 : ERROR -',ch10,&
> &  '  ifcflag must be 1 when the calculation of the phonon DOS is required
> ',ch10,&
> &  '  Action : correct ifcflag in your input file.'
>   call wrtout(6,message,'COLL')
>   call leave_new('COLL')
>  end if
>
>
> 2) If I remember well in old versions of abinit there's a small bug (just
> a missing factor)
>    in the calculation of the integrated type-projected PHDOS. So don't
> rely on the columns
>    reporting the type-projected IPHDOS. If you need, you can easily
> integrate the curve using
>    a plotting software such as xmgrace.
>
> 3) Be careful when using the tetrahedron method. The CPU time doesn't
> scale well wrt
>    to the number of Fourier interpolated q-points. This is a general
> problem I found also
>    in case of electronic DOS calculations. One should do some serious
> profiling of the tetrahedron
>    routines and the fortran loops to individuate the bottleneck.
>
> Hope it helps.
> Best Regards
> Matteo Giantomassi
>
>