File: genscph.f90

package info (click to toggle)
elkcode 5.4.24-2
  • links: PTS, VCS
  • area: main
  • in suites: buster, sid
  • size: 12,840 kB
  • sloc: f90: 48,415; fortran: 22,457; perl: 965; makefile: 384; sh: 369; python: 105; ansic: 67
file content (63 lines) | stat: -rw-r--r-- 1,731 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63

! Copyright (C) 2002-2005 J. K. Dewhurst, S. Sharma and C. Ambrosch-Draxl.
! This file is distributed under the terms of the GNU General Public License.
! See the file COPYING for license details.

subroutine genscph(p,dph)
use modmain
use modphonon
use modstore
implicit none
! arguments
integer, intent(in) :: p
real(8), intent(in) :: dph
! local variables
integer is,ia,na,i
real(8) vc(3),t1
if ((p.ne.0).and.(p.ne.1)) then
  write(*,*)
  write(*,'("Error(genscph): phase (p) should be 0 or 1 : ",I8)') p
  write(*,*)
  stop
end if
! find the smallest supercell which contains the q-vector
call findscq(iqph,avec_,nscph,vscph)
! construct supercell atomic positions and magnetic fields
do is=1,nspecies
  na=0
  do ia=1,natoms_(is)
    do i=1,nscph
      na=na+1
      if (na.gt.maxatoms) then
        write(*,*)
        write(*,'("Error(genscph): too many atoms in supercell : ",I8)') na
        write(*,'(" for species ",I4)') is
        write(*,'("Adjust maxatoms in modmain and recompile code")')
        write(*,*)
        stop
      end if
      vc(:)=vscph(:,i)+atposc_(:,ia,is)
! add small periodic displacement
      if ((isph.eq.is).and.(iaph.eq.ia)) then
        t1=dot_product(vqc(:,iqph),vscph(:,i))
        if (p.eq.0) then
          vc(ipph)=vc(ipph)+dph*cos(t1)
        else
          vc(ipph)=vc(ipph)+dph*sin(t1)
        end if
      end if
! convert to new lattice coordinates
      call r3mv(ainv,vc,atposl(:,na,is))
      call r3frac(epslat,atposl(:,na,is))
! set muffin-tin fields and fixed spin moments if required
      if (spinpol) then
        bfcmt0(:,na,is)=bfcmt0_(:,ia,is)
        mommtfix(:,na,is)=mommtfix_(:,ia,is)
      end if
    end do
  end do
  natoms(is)=na
end do
return
end subroutine