File: readdvs.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 (101 lines) | stat: -rw-r--r-- 2,810 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
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101

! Copyright (C) 2008 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 readdvs(iq,is,ia,ip)
use modmain
use modphonon
implicit none
! arguments
integer, intent(in) :: iq,is,ia,ip
! local variables
integer js,jas,iostat
integer version_(3),nspecies_
integer lmmaxo_,natoms_
integer nrmt_,ngridg_(3)
character(256) fext,fname
! allocatable arrays
complex(8), allocatable :: zfmt(:,:,:)
allocate(zfmt(lmmaxo,nrmtmax,natmtot))
call phfext(iq,is,ia,ip,fext)
fname='DVS'//trim(fext)
open(50,file=trim(fname),form='UNFORMATTED',status='OLD',iostat=iostat)
if (iostat.ne.0) then
  write(*,*)
  write(*,'("Error(readdvs): error opening ",A)') trim(fname)
  write(*,*)
  stop
end if
read(50) version_
if ((version(1).ne.version_(1)).or.(version(2).ne.version_(2)) &
 .or.(version(3).ne.version_(3))) then
  write(*,*)
  write(*,'("Warning(readdvs): different versions")')
  write(*,'(" current : ",I3.3,".",I3.3,".",I3.3)') version
  write(*,'(" file    : ",I3.3,".",I3.3,".",I3.3)') version_
  write(*,'(" in file ",A)') trim(fname)
end if
read(50) nspecies_
if (nspecies.ne.nspecies_) then
  write(*,*)
  write(*,'("Error(readdvs): differing nspecies")')
  write(*,'(" current : ",I4)') nspecies
  write(*,'(" file    : ",I4)') nspecies_
  write(*,'(" in file ",A)') trim(fname)
  write(*,*)
  stop
end if
read(50) lmmaxo_
if (lmmaxo.ne.lmmaxo_) then
  write(*,*)
  write(*,'("Error(readdvs): differing lmmaxo")')
  write(*,'(" current : ",I4)') lmmaxo
  write(*,'(" file    : ",I4)') lmmaxo_
  write(*,'(" in file ",A)') trim(fname)
  write(*,*)
  stop
end if
do js=1,nspecies
  read(50) natoms_
  if (natoms(js).ne.natoms_) then
    write(*,*)
    write(*,'("Error(readdvs): differing natoms for species ",I4)') js
    write(*,'(" current : ",I4)') natoms(js)
    write(*,'(" file    : ",I4)') natoms_
    write(*,'(" in file ",A)') trim(fname)
    write(*,*)
    stop
  end if
  read(50) nrmt_
  if (nrmt(js).ne.nrmt_) then
    write(*,*)
    write(*,'("Error(readdvs): differing nrmt for species ",I4)') js
    write(*,'(" current : ",I6)') nrmt(js)
    write(*,'(" file    : ",I6)') nrmt_
    write(*,'(" in file ",A)') trim(fname)
    write(*,*)
    stop
  end if
end do
read(50) ngridg_
if ((ngridg(1).ne.ngridg_(1)).or.(ngridg(2).ne.ngridg_(2)).or. &
 (ngridg(3).ne.ngridg_(3))) then
  write(*,*)
  write(*,'("Error(readdvs): differing ngridg")')
  write(*,'(" current : ",3I6)') ngridg
  write(*,'(" file    : ",3I6)') ngridg_
  write(*,'(" in file ",A)') trim(fname)
  write(*,*)
  stop
end if
read(50) zfmt,dvsir
do jas=1,natmtot
  js=idxis(jas)
  call zfmtpack(.true.,nrmt(js),nrmti(js),zfmt(:,:,jas),dvsmt(:,jas))
end do
close(50)
deallocate(zfmt)
return
end subroutine