[Wien] Plot bandstructure error!!!
Peter Blaha
pblaha at theochem.tuwien.ac.at
Fri Nov 11 16:49:41 CET 2005
I checked with your files and yes, there is a problem in spagh.f
It puts for each color definition the 0.0 on the stack (this 0.0 is not
used by setcolor) and thus for big plots (many k-points + many eigenvalues)
and depending on your systemmemory you will get this stackoverflow
(On my Linux it displays about half of the plot until the error comes out,
my printer can even print 2/3 of the plot).
Anyway, I modified spag.f and modules.f to fix this.
I also used the opportunity to reduce the filesize significantly by
removing unnecessary "translate statements", "plotting points outside
the displayed area" and removing the many unnecessary cXX color definitions.
I'll include spag.f and modules.f with this mail.
spaghetti will also be fixed with the newxt update.
Regards
> Thanks.
> My calculation contains 22 atoms with a P1 space group. After scf, I just
> followed the basic steps to calculate band. I inserted the EF and changed
> the E-window in case.insp file, then ran x spaghetti -c -so. I clicked the
> "plot band structure" button, and I can't see the fig, just a red cross
> image. I downloaded the hardcopy in postscript format then open it with
> acrobat distiller (I can't open it with gv either), the error information is
> this:
>
> %%[ Warning: Times-Roman not found, using Font Substitution. Font cannot be
> embedded. ]%%
> %%[ Error: stackoverflow; OffendingCommand: 0.2 ]%%
>
> Stack:
> [0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
> 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
> 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
> 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0
> 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0.....
> I can see the case.spaghetti_ps file. It is not empty and there is a
> "showpage" in the last line. But the last part of this file is all zero like
> this:
> " .........
> .........
> 0.0 c13
> 0.0 c14
> 0.0 c15
> 0.0 c16
> 0.0 c17
> 0.0 c18
> stroke
> showpage "
> btw. the error files are all empty.
> So what's the reason for this problem? Any suggestions will be appreciated.
> Best wished! Arin
>
> 2005/11/9, Peter Blaha <pblaha at theochem.tuwien.ac.at>:
> >
> > We can't hardly give any suggestions, when you do not provide details.
> >
> > > > > k-path with xcrysden, then calculated the band structure:
> > > > > x lapw1 -c -band
> > > > > x lapwso -c
> > > > > x spaghetti -so -c
> >
> > These steps look basically ok. If there was no error message, it should
> > have
> > produced a file case.spaghetti_ps
> > Have you looked into that file ? Is it empty ? does it have a "showpage"
> > in the
> > last line ?
> > Did you change anything in the case.insp file (besides putting the correct
> > EF)?
> > What means "cannot plot bandstructure" ?
> >
> > have you tried gv case.spaghetti_ps ?
> > Do you get any error message ?
> >
> > More than 850 other users seem to have no problems.
> >
> >
> >
> > > I didn't ran irrep.
> > > I just followed the steps: lapw1, lapwso, then edited the case.inspfile,
> > > then spaghetti.
> > > I asked someone else before. They seemed to have the same question,
> > > especially for complex calculation with a lot of atoms.
> > > How can I fix this problem? Any suggestions?thanks
> > >
> > >
> > > 2005/11/9, Peter Blaha <pblaha at theochem.tuwien.ac.at>:
> > > >
> > > > Could it be that you ran irrep before ?
> > > >
> > > > spaghetti will crash if you have an old case.irrep file (one which
> > does
> > > > not
> > > > fit to the current outputso file).
> > > >
> > > > remove any case.irrep* file and run spaghetti again.
> > > >
> > > > Any other message, when you run spaghetti ???
> > > >
> > > > > I am doing a complex calculation with SO. After calculation, I
> > selected
> > > > the
> > > > > k-path with xcrysden, then calculated the band structure:
> > > > > x lapw1 -c -band
> > > > > x lapwso -c
> > > > > x spaghetti -so -c
> > > > > Finally I plotted the band structure, but it could not be displayed.
> > I
> > > > > downloaded the .ps file and there was an error when I opened this
> > file.
> > > > I
> > > > > can see the eigenvalues in case.outputso file but can't plot the
> > band
> > > > > structure.
> > > > > It is the newest version of wien2k. I remember when I used an old
> > > > version,
> > > > > the band could not be plotted either, but sometimes I could open the
> > .ps
> > > > > file to plot the band. Now I don't know what is the reason. Can
> > anybody
> > > > give
> > > > > me a hint? Or this is a bug of the program? I appreciate your help!
> > > > > Best Regards
> > > > > Arin
> > > > >
> > > >
> > > >
> > > > P.Blaha
> > > >
> > --------------------------------------------------------------------------
> > > > Peter BLAHA, Inst.f. Materials Chemistry, TU Vienna, A-1060 Vienna
> > > > Phone: +43-1-58801-15671 FAX: +43-1-58801-15698
> > > > Email: blaha at theochem.tuwien.ac.at WWW:
> > http://info.tuwien.ac.at/theochem/
> > > >
> > --------------------------------------------------------------------------
> > > > _______________________________________________
> > > > Wien mailing list
> > > > Wien at zeus.theochem.tuwien.ac.at
> > > > http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien
> > > >
> > >
> >
> >
> > P.Blaha
> > --------------------------------------------------------------------------
> > Peter BLAHA, Inst.f. Materials Chemistry, TU Vienna, A-1060 Vienna
> > Phone: +43-1-58801-15671 FAX: +43-1-58801-15698
> > Email: blaha at theochem.tuwien.ac.at WWW: http://info.tuwien.ac.at/theochem/
> > --------------------------------------------------------------------------
> > _______________________________________________
> > Wien mailing list
> > Wien at zeus.theochem.tuwien.ac.at
> > http://zeus.theochem.tuwien.ac.at/mailman/listinfo/wien
> >
>
P.Blaha
--------------------------------------------------------------------------
Peter BLAHA, Inst.f. Materials Chemistry, TU Vienna, A-1060 Vienna
Phone: +43-1-58801-15671 FAX: +43-1-58801-15698
Email: blaha at theochem.tuwien.ac.at WWW: http://info.tuwien.ac.at/theochem/
--------------------------------------------------------------------------
-------------- next part --------------
!$hp9000_800 intrinsics on
program spag
!
!.....Generates postscript file showing the band structure of the
! Kohn-Sham eigenvalues obtained from the Wien2k program.
! Includes averaging over degenerate bands (SO) and multiple atoms
!
! last updated: Dec. 2001 Peter Blaha, Vienna Techn University
! f77->f90
! Mars 2002 Clas Persson, Uppsala University
! incorporated irreducible representations
! for drawing lines of the bands.
! Peter Blaha, testing
!
use reallocate
use irr_param
IMPLICIT REAL*8 (A-H,O-Z)
!
character aline*80,fname*80
CHARACTER*11 STATUS,FORM
character k_pat*7,ei_pat*28
character dummy*1,label1*12
character symbol*12,label*12
character title*80,lattice*4
logical flg1
!
integer,allocatable :: jatom_list(:) ! NATO)
real*8,pointer :: eigen(:,:),charr(:,:) ! NEVL,NKP)
integer,pointer :: n_ene(:),lines(:) ! NKP)
real*8,pointer :: vk(:,:) ! 3,NKP)
real*8,pointer :: xval(:) ! NKP)
logical,pointer :: break(:) ! nkp)
character*12,pointer :: k_name(:) ! nkp
dimension qtl(13),icomma(15)
data k_pat /' K='/
data ei_pat /'EIGENVALUES BELOW THE ENERGY'/
!-----------------------------------------------------------------------
!
!.....INITIALIZE VARS
n_kpt=0
NEVL=1000
NKP=1000
allocate ( eigen(nevl,nkp),charr(nevl,nkp) )
allocate ( n_ene(nkp),k_name(nkp),lines(nkp) )
allocate ( vk(3,NKP))
allocate ( xval(NKP))
allocate ( break(nkp))
!
allocate ( eneik(NEVL,NKP) )
allocate ( ngrp(4,NEVL,NKP), ngde(4,NEVL,NKP) ) !(4,NEVL,NKP)
allocate ( ndeg(NEVL,NKP) )
allocate ( neik(nkp),ikg(NKP),ntab4(NKP) )
allocate ( iscol(NEVL),iscol2(NEVL),iste(NEVL),ibnd(NEVL) ) !NEVL)
allocate ( lkg(48,48,NKP) )
allocate ( nlkg(48,NKP) )
allocate ( grpnam(NKP) )
allocate ( cnam(48,NKP) )
!
iarg=iargc()
if(iarg.ne.1) goto 900
call getarg(1,fname)
OPEN(1,FILE=fname,STATUS='OLD',ERR=8000)
8003 READ(1,*,END=8001) IUNIT,FNAME,STATUS,FORM,IRECL
OPEN(IUNIT,FILE=FNAME,STATUS=STATUS,FORM=FORM, &
ERR=8002)
GOTO 8003
8000 WRITE(*,*) ' ERROR IN OPENING SPAG.DEF !!!!'
STOP 'SPAG.DEF'
8002 WRITE(*,*) ' ERROR IN OPENING UNIT:',IUNIT
WRITE(*,*) ' FILENAME: ',FNAME,' STATUS: ',STATUS, &
' FORM:',FORM
STOP 'OPEN FAILED'
8001 CONTINUE
!.....write date
call wrtdate(6)
!
label=' '
label1=' '
do i=1,12
if(fname(i:i).ne.'.') then
label(i:i)=fname(i:i)
else
! label(i:i)=','
goto 100
endif
enddo
!
!.....READ K-VECTORS
!
100 CONTINUE
read(7,'(a80)',end=200) aline
if (aline(1:7).eq.k_pat) then
! THIS IS THE BEGINNING OF AN EIGENVALUE SECTION (output1)
n_kpt=n_kpt + 1
if(n_kpt.gt.nkp) then
nkp=nkp*2
call doreallocate(eigen,nevl,nkp)
call doreallocate(charr,nevl,nkp)
call doreallocate(n_ene,nkp)
call doreallocate(k_name,nkp)
call doreallocate(lines,nkp)
call doreallocate(vk,3,nkp)
call doreallocate(xval,nkp)
call doreallocate(break,nkp)
call doreallocate ( eneik,NEVL,NKP )
call doreallocate ( ngrp,4,NEVL,NKP)
call doreallocate ( ngde,4,NEVL,NKP)
call doreallocate ( ndeg,NEVL,NKP)
call doreallocate ( neik,nkp)
call doreallocate ( ikg,NKP)
call doreallocate ( ntab4,NKP)
call doreallocate ( lkg,48,48,NKP)
call doreallocate ( nlkg,48,NKP)
call doreallocate ( grpnam,NKP)
!
write(6,*) ' parameter NKP increased:',nkp
end if
n_ene(n_kpt)=0
call get_k (aline,vk(1,n_kpt),vk(2,n_kpt), &
vk(3,n_kpt),k_name(n_kpt))
read(7,'(a1)') dummy
read(7,'(a1)') dummy
! READ EIGENVALUES OF K-VECTOR
110 read(7,'(a80)') aline
if (aline(15:42).eq.ei_pat) then
goto 100
else
if(n_ene(n_kpt).ge.nevl-4) then
nevl=nevl*2
call doreallocate(eigen,nevl,nkp)
call doreallocate(charr,nevl,nkp)
call doreallocate ( eneik,NEVL,NKP )
call doreallocate ( ngrp,4,NEVL,NKP)
call doreallocate ( ngde,4,NEVL,NKP)
call doreallocate ( ndeg,NEVL,NKP)
call doreallocate ( iscol,NEVL)
call doreallocate ( iscol2,NEVL)
call doreallocate ( iste,NEVL)
call doreallocate ( ibnd,NEVL)
write(6,*) 'parameter nevl increased:',nevl
endif
call get_ei(aline,eigen(1,n_kpt), &
n_ene(n_kpt))
goto 110
endif
else
goto 100
endif
!
!.....ALL K-VECTORS HAVE BEEN READ; SEARCH FOR K-POINT WITH SMALLEST
! NUMBER OF EIGENVALUES
!
200 continue
write(*,*) 'number of k-points read in case.vector=',n_kpt
nu_min=999
do 205 j=1,n_kpt
if (n_ene(j).lt.nu_min) then
nu_min=n_ene(j)
k_min=j
endif
205 continue
write(6,*) 'smallest number eigenvalues at k=',k_min,' (', &
k_name(k_min),')'
write(6,*) ' =',nu_min
!
!.....READ IRREDUCIBLE REPRESENTATIONS OF THE EIGENVECTOR
nikir=0
read(30,fmt=511,end=444,err=444) tjunk
!
read(30,511) tjunk
write(6,*) 'using irred. representations'
440 read(30,fmt=502,end=445) ikir, skir(1), skir(2), skir(3)
nikir=nikir+1
read(30,504) grpnam(nikir),ikg(nikir),ntab4(nikir)
if(grpnam(nikir).eq.'Nop') then
write(*,*) 'k-point nr ',nikir,' not treated with irrep'
write(6,*) 'k-point nr ',nikir,' cannot be characterized according to any point group.'
write(6,*) 'here, we draw bands between closest energy values.'
write(6,*) 'thus, we treat the k-point as point group C1.'
grpnam(nikir)='C1 '
endif
read(30,*)
do 460 itab=1,ntab4(nikir)
read(30,507) cnam(itab,nikir)
!
!.........find number of group elements
flg1=.true.
nlkg(itab,nikir)=1
icnam=1
do while(flg1)
select case ( cnam(itab,nikir)(icnam:icnam) )
case ('E')
flg1=.false.
case ('I')
flg1=.false.
case ('C')
flg1=.false.
case ('2')
flg1=.false.
nlkg(itab,nikir)=2
case ('3')
flg1=.false.
nlkg(itab,nikir)=3
case ('4')
flg1=.false.
nlkg(itab,nikir)=4
case ('6')
flg1=.false.
nlkg(itab,nikir)=6
case ('8')
flg1=.false.
nlkg(itab,nikir)=8
end select
icnam=icnam+1
if(icnam.eq.6) then
write(*,*) 'ERROR spag: cannot find nr: ',cnam(itab,nikir)
write(6,*) 'ERROR spag: cannot find nr: ',itab,nikir,icnam
write(6,*) 'ERROR spag: cannot find nr: ',cnam(itab,nikir)
! call flush(6)
stop 'spag: cannot find nr of pntgrp '
endif
enddo
read(30,508) (lkg(itab,inge,nikir), inge=1,nlkg(itab,nikir))
460 continue
!
read(30,503) neik(nikir)
write(6,510) 'ik=',nikir,'; k=',skir(1), skir(2), skir(3), &
'point group= ', grpnam(nikir), &
'; no of symm.ops=', ikg(nikir), &
'; no of classes=', ntab4(nikir), &
'no of energies=', neik(nikir)
do itab=1,ntab4(nikir)
write(6,506) cnam(itab,nikir)
write(6,508) (lkg(itab,inge,nikir),inge=1,nlkg(itab,nikir))
enddo
do 480 ie=1,neik(nikir)
read(30,518) inum, ndeg(ie,nikir),eneik(ie,nikir), &
(ngrp(ipgr,ie,nikir),ngde(ipgr,ie,nikir),ipgr=1,4)
480 continue
! call flush(6)
goto 440
444 write(6,*) 'do not use irred. representations'
445 write(6,*) 'number of k-points read in case.irrep= ',nikir
!
499 continue
502 format(1i6,3f10.6)
503 format(1i10)
504 format(3x,a3,i3,i10)
506 format(x,a6,$)
507 format(x,a6)
508 format(48i3)
518 format(2i10,f10.6,4(2x,2i3))
510 format(a3,i6,a4,3f10.6,/,a13,a3,a17,i5,a16,i5,/,a15,i5)
511 format(a3)
!
!.....READ STRUCT FILE
!
read(20,'(a80)') title
read(20,'(a4,23x,i3)') lattice,natom
write(6,*) 'lattice is=',lattice,natom,' atoms'
allocate (jatom_list(natom+1))
read(20,*)
read(20,457) aa,bb,cc , alpha,beta,gamma
if(alpha.eq.0.d0) alpha=90.d0
if(beta.eq.0.d0) beta=90.d0
if(gamma.eq.0.d0) gamma=90.d0
457 format(6f10.7)
write(6,*) ' a=',aa
write(6,*) ' b=',bb
write(6,*) ' c=',cc
!
!.....TRANSFORM K-VECTORS TO CARTESIAN COORDINATES
!
call cartco (vk,n_kpt,lattice,aa,bb,cc,alpha,beta,gamma)
write(6,*) 'transforming k-points to cartesian coordinates...'
!
!.....SEARCH FOR THE ENDPOINTS OF THE STRAIGHT LINES THROUGH THE
! BRILLOUIN-ZONE; CALCULATE THE CUMULATIVE DISTANCES BETWEEN
! THE K-POINTS, WHICH ARE THE X-VALUES FOR THE BANDSTRUCTURE
! PLOT
!
write(6,*) 'checking for lines...'
call bz_lin(vk,n_kpt,lines,n_lin,xval,nbreak,break)
write(6,505) n_lin,nbreak
505 format('numbers of BRILLOUIN-zone lines found:',i3,/, &
'numbers of breaks',i3)
do 210 j=1,n_lin
write(6,*) 'line',j,' is: ',k_name(lines(j)),' to ', &
k_name(lines(j+1))
write(6,*) 'xmax=',xval(lines(j))
210 continue
!
!.....READ VIEWING PARAMETERS
!
call inview (ymin,ymax,xcm,ycm,xoffs,yoffs,height,efermi, &
lyflag,nb_min,nb_max, &
jatom,jatom_list,jtype,sizec,ticks,nticks)
xmin=0.d0
xmax=xval(n_kpt)
if(lyflag.eq.2) call con_ev(n_ene,n_kpt, &
efermi,eigen,nevl)
!
! ....READ QTL FILE
!
char0=0.01d0
do 208 i=1,n_kpt
do 208 ii=1,n_ene(i)
208 charr(ii,i)=char0
if(jatom.eq.0) goto 206
207 read(9,'(a80)',end=206) aline
if(aline(2:5).eq.'JATO') then
read(aline(8:9),'(i2)') iatom
if(iatom.eq.jatom_list(1)) then
i0=2
icomma(1)=31
do i=32,80
if(aline(i:i).eq.',') then
icomma(i0)=i
i0=i0+1
endif
enddo
icomma(i0)=icomma(i0-1)+2
label1=aline(icomma(jtype)+1:icomma(jtype+1)-1)
if(label1(1:1).eq.'0') label1(1:1)='s,'
if(label1(1:1).eq.'1') label1(1:1)='p,'
if(label1(1:1).eq.'2') label1(1:1)='d,'
if(label1(1:1).eq.'3') label1(1:1)='f,'
endif
endif
if (aline(2:5).ne.'BAND') goto 207
jband=0
211 jband=jband+1
if(jband.ge.nb_min.and.jband.le.nb_max) then
chmin=1.0d0
chmax=0.0d0
do 209 i=1,n_kpt
do 209 ii=1,natom+1
do jatom1=1,jatom
if(ii.ne.jatom_list(jatom1)) goto 5210
read(9,'(13X,F8.5,3X,12F8.5)',err=212,end=212) &
(qtl(j),j=1,13)
charr(jband,i)=charr(jband,i)+qtl(jtype)*sizec
! write(*,*) jband,i,charr(jband,i)
if(chmin.gt.charr(jband,i)) chmin=charr(jband,i)
if(chmax.lt.charr(jband,i)) chmax=charr(jband,i)
goto 209
5210 continue
enddo
read(9,*,err=212,end=212)
209 continue
write(6,*) ' band',jband,' char-min/max:',chmin,chmax
else
do 213 i=1,n_kpt
do 213 ii=1,natom+1
213 read(9,*,err=212,end=212)
end if
if(jband.eq.nu_min) goto 206
read(9,*,err=212,end=212)
goto 211
212 write(*,*) ' ERROR reading QTLs:'
write(*,*) ' band:',jband,' k-point:',i
write(*,*) ' execution continued ....'
206 continue
!.....averaging character of degenerate bands (max 6 fold)
do jk=1,n_kpt
je=1
555 if(eigen(je,jk)-eigen(je+1,jk).gt.-5.d-6) then
if(eigen(je,jk)-eigen(je+2,jk).gt.-5.d-6) then
if(eigen(je,jk)-eigen(je+3,jk).gt.-5.d-6) then
if(eigen(je,jk)-eigen(je+4,jk).gt.-5.d-6) then
if(eigen(je,jk)-eigen(je+5,jk).gt.-5.d-6) then
charr(je,jk)=(charr(je,jk)+charr(je+1,jk)+charr(je+2,jk)+ &
charr(je+3,jk)+charr(je+4,jk)+charr(je+5,jk))/6.d0
charr(je+1,jk)=charr(je,jk)
charr(je+2,jk)=charr(je,jk)
charr(je+3,jk)=charr(je,jk)
charr(je+4,jk)=charr(je,jk)
charr(je+5,jk)=charr(je,jk)
je=je+6
goto 556
endif
charr(je,jk)=(charr(je,jk)+charr(je+1,jk)+charr(je+2,jk)+ &
charr(je+3,jk)+charr(je+4,jk))/5.d0
charr(je+1,jk)=charr(je,jk)
charr(je+2,jk)=charr(je,jk)
charr(je+3,jk)=charr(je,jk)
charr(je+4,jk)=charr(je,jk)
je=je+5
goto 556
endif
charr(je,jk)=(charr(je,jk)+charr(je+1,jk)+charr(je+2,jk)+ &
charr(je+3,jk))/4.d0
charr(je+1,jk)=charr(je,jk)
charr(je+2,jk)=charr(je,jk)
charr(je+3,jk)=charr(je,jk)
je=je+4
goto 556
endif
charr(je,jk)=(charr(je,jk)+charr(je+1,jk)+charr(je+2,jk) &
)/3.d0
charr(je+1,jk)=charr(je,jk)
charr(je+2,jk)=charr(je,jk)
je=je+3
goto 556
endif
charr(je,jk)=(charr(je,jk)+charr(je+1,jk) &
)/2.d0
charr(je+1,jk)=charr(je,jk)
je=je+2
goto 556
endif
je=je+1
556 if(je.lt.n_ene(jk)) goto 555
write(6,557) jk,(je,eigen(je,jk),charr(je,jk),je=1,n_ene(jk))
enddo
557 format(i4,(i5,2f10.6))
!
!.....initialize ps
!
pi=acos(-1.d0)
write(6,*)xoffs,yoffs
call psinit (xoffs,yoffs,xcm,ycm,dlwid)
write(6,*) ' init done'
! labelling:
xval1=0.0d0
call movet (xval1,ycm+1.d0)
call writs (label,1.d0*height,'label ')
if(jatom.gt.0) then
write(label,878) (jatom_list(j),j=1,min(3,jatom))
878 format(' atom ',3i2)
call writs (label,1.d0*height,'label ')
call writs (label1,1.d0*height,'label ')
write(label,879) sizec
879 format(' size',f5.2)
call writs (label,1.d0*height,'label ')
endif
!
!.....draw vertical lines
!
write(6,*) n_lin
write(11,*) '%vertical lines'
do 220 j=1,n_lin+1
xval1=xval(lines(j))*xcm/xmax
write(11,140) dawid
call movet (xval1,0.d0)
call drawt (xval1,ycm)
220 continue
write(6,*) ' vertical lines done'
!
!.....draw horizontal lines and FERMI-energy
!
write(11,*) '%horizontal lines and FERMI-energy'
do 230 j=1,n_lin
if (.not.break(lines(j+1))) then
xval1=xval(lines(j))*xcm/xmax
xval2=xval(lines(j+1))*xcm/xmax
call movet (xval1,0.d0)
call drawt (xval2,0.d0)
call movet (xval1,ycm)
call drawt (xval2,ycm)
write(11,*) 'stroke'
write(11,140) dawid
!
if (iprtf.ne.0) then
!! if (efermi.ne.999.d0) then
yval1=ycm/(ymax-ymin)*(efermi-ymin)
write(11,'(A41)')
write(11,*) 'stroke'
write(11,140) dfwid
!
if(iprtf.eq.1) then
call movet (xval1,yval1)
call drawt (xval2,yval1)
elseif(iprtf.eq.2) then
xvaldel = (xval2-xval1)/18.5
xval2 = xval1+xvaldel/2.d0
do j1=1,19
call movet (xval1,yval1)
call drawt (xval2,yval1)
xval1 = xval1+xvaldel
xval2 = xval1+xvaldel/2.d0
enddo
elseif(iprtf.eq.3) then
xvaldel = (xval2-xval1)/41
xval2 = xvaldel
col='setgray'
call movet (xval1,yval1)
do j1=1,41
!c call movet (xval1,yval1)
call transt(xval1,yval1)
call pointi (0,0.01d0,iprto)
call transt(-xval1,-yval1)
xval1 = xval1+xvaldel
enddo
else
stop 'error in input, incorrect line for Fermi'
endif
!
write(11,*) 'stroke'
write(11,140) dawid
write(6,*) ' EF-line done'
endif
endif
230 continue
write(6,*) ' horizontal lines and FERMI-energy done'
write(11,*) 'stroke'
write(11,140) dlwid
!
!.....label x-axis with names of k-points
! write the names of the HIGH-SYMMETRY k-points
!
do 240 j=1,n_kpt
symbol=k_name(j)
if (symbol(1:1).eq.' ') goto 240
xval1=xval(j)*xcm/xmax
call movet(xval1-0.2d0,-0.5d0)
call writs (symbol,0.8d0*height,'labels')
240 continue
write(6,*) ' x-axis labeling done'
!
!.....label FERMI-energy
!
! if (efermi.ne.999.d0) then
if (iprtf.ne.0) then
yval1=ycm/(ymax-ymin)*(efermi-ymin)
call movet (xcm+0.1d0,yval1-0.2d0)
call writs ('E ',1.d0*height,'fermi ')
call movet (xcm+0.5d0,yval1-0.4d0)
call writs ('F ',0.8d0*height,'fermi ')
write(6,*) ' EF-labeling done'
endif
!
!.....label y-axis with 'Energy [Ry] '
!
yval1=ycm/2.d0-2.d0
call movet (-2.7d0,yval1)
write(11,*) '90. rotate'
if (lyflag.eq.1) then
call writs ('Energy (Ry) ',1.d0*height,'label ')
else if (lyflag.eq.2) then
call writs ('Energy (eV) ',1.d0*height,'label ')
endif
write(11,*) '-90. rotate'
write(11,*) 'stroke'
write(6,*) ' y-axis labeling done'
!
!.....label y-axis with ticks
!
if (lyflag.eq.1) then
ymin1=ymin
246 yval1=ycm/(ymax-ymin)*(ymin1-ymin)
call movet (0.0d0,yval1)
call drawt (-0.3d0,yval1)
call writz (ymin1,0.7*height,yval1,lyflag)
do 245 i=1,nticks
if ((ymin1+ticks*i/(nticks+1)).ge.ymax) goto 247
yval1=ycm/(ymax-ymin)*(ymin1+ticks*i/(nticks+1)-ymin)
call movet (0.0d0,yval1)
call drawt (-0.2d0,yval1)
245 continue
ymin1=ymin1+ticks
if (ymin1.le.ymax) goto 246
247 continue
else if (lyflag.eq.2) then
ymin1=efermi
256 yval1=ycm/(ymax-ymin)*(ymin1-ymin)
call movet (0.0d0,yval1)
call drawt (-0.3d0,yval1)
call movet (0.0d0,yval1)
call writz (ymin1,0.7d0*height,yval1,lyflag)
do 255 i=1,nticks
if ((ymin1+ticks*i/(nticks+1)).ge.ymax) goto 257
yval1=ycm/(ymax-ymin)*(ymin1+ticks*i/(nticks+1)-ymin)
call movet (0.0d0,yval1)
call drawt (-0.2d0,yval1)
255 continue
ymin1=ymin1+ticks
if (ymin1.le.ymax) goto 256
257 continue
ymin1=efermi
259 ymin1=ymin1-ticks
if (ymin1.lt.ymin) goto 258
yval1=ycm/(ymax-ymin)*(ymin1-ymin)
call movet (0.0d0,yval1)
call drawt (-0.3d0,yval1)
call movet (0.0d0,yval1)
call writz (ymin1,0.7d0*height,yval1,lyflag)
do 265 i=1,nticks
if ((ymin1+ticks*i/(nticks+1)).le.ymin) goto 258
yval1=ycm/(ymax-ymin)*(ymin1+ticks*i/(nticks+1)-ymin)
call movet (0.0d0,yval1)
call drawt (-0.2d0,yval1)
265 continue
goto 259
258 continue
endif
!
!.....set clipboard
!
write(11,*) 'stroke'
if (nbreak.eq.0) then
call clipin (0.d0,0.d0,xcm,ycm)
else
do 500 j=1,n_lin
if(break(lines(j+1))) then
xval1=xval(lines(j))*xcm/xmax
call clipin (0.d0,0.d0,xval1,ycm)
goto 501
endif
500 continue
501 continue
endif
!
!.....draw energies
!
write(6,*) 'drawing energy-bands...'
do 607 je2=1,NEVL
iscol(je2) =je2
iscol2(je2)=je2
607 continue
!
!.....for each k-point
do 250 jk=1,n_kpt
if (break(jk)) then
do 600 jjj=jk+1,n_kpt
if (break(jjj)) then
xval2=xval(jjj-1)*xcm/xmax
xval1=xval(jk)*xcm/xmax
write(11,*) 'stroke'
call clipin (xval1,0.d0,xval2,ycm)
goto 601
endif
600 continue
xval1=xval(jk)*xcm/xmax
write(11,*) 'stroke'
call clipin (xval1,0.d0,xcm,ycm)
601 continue
endif
!........order of the lines
do 608 je2=1,n_ene(jk)
608 iste(je2)=je2
!
if(nikir.ne.0.and.jk.lt.n_kpt) then
if(n_ene(jk ).ne.neik(jk )) stop 'no of egval 1'
if(n_ene(jk+1).ne.neik(jk+1)) stop 'no of egval 2'
if(k_name(jk )(1:1).eq.' '.or. &
k_name(jk+1)(1:1).eq.' ') call bndind(jk,nevl)
endif
!
!.....for each energy state
do 260 je=1,n_ene(jk)
!
!........define colors
col='0.0 setgray'
if(iprto(2).ne.0) then
if(iprto(2).eq.1) then
nprtc = 1
write(col,6001) nprtc
elseif(iprto(2).eq.2) then
nprtc = 2 + mod(iscol2(je),3)
write(col,6001) nprtc
elseif(iprto(2).eq.3) then
nprtc = 10 + mod(iscol2(je),12)
write(col,6010) nprtc
elseif(iprto(2).eq.4) then
nprtc = 10 + mod(iscol2(je),12)
write(col,6010) nprtc
else
stop 'color switch > 4 not implemented'
endif
if(nprtc.gt.22) stop 'spag: need more color definitions'
!
endif
! write(11,*) '0.0 ',col
!
6001 format('c0',I1,8X)
6010 format('c', I2,8X)
xval1=xval(jk)*xcm/xmax
yval1=ycm/(ymax-ymin)*(eigen(je,jk)-ymin)
if (yval1.gt.ycm) goto 260
if (yval1.lt.0.d0) goto 260
write(11,*) col
if(iprto(1).ne.1) then
call transt(xval1,yval1)
if (charr(je,jk).le.char0) then
call pointi (0,0.01d0,iprto)
else
call pointi (1,charr(je,jk),iprto)
endif
call transt(-xval1,-yval1)
endif
!........draw lines of the energy bands
if(k_name(jk )(1:1).eq.' '.or. &
k_name(jk+1)(1:1).eq.' ') then
if((iprto(1).ne.0).and.(jk.ne.n_kpt).and. &
(je.le.n_ene(jk+1)).and.iste(je).ne.0 ) then
xval11=xval(jk )*xcm/xmax
yval11=ycm/(ymax-ymin)*(eigen(je,jk)-ymin)
xval22=xval(jk+1)*xcm/xmax
yval22=ycm/(ymax-ymin)*(eigen(iste(je),jk+1)-ymin)
call writln(xval11,yval11,xval22,yval22)
endif
endif
260 continue
!
if(iprto(2).ne.4) then
do 609 je2=1,n_ene(jk)
609 iscol2(iste(je2))=iscol(je2)
do 610 je2=1,n_ene(jk)
610 iscol(je2)=iscol2(je2)
endif
!
250 continue
!
!..... write energies bandwise on file
!
write(6,*) 'writing energy-bands on file'
do 270 je=nb_min,min(nu_min,nb_max)
write(10,*) ' bandindex:',je
do 280 jk=1,n_kpt
280 write(10,227)vk(1,jk),vk(2,jk),vk(3,jk),xval(jk),eigen(je,jk)
! 280 write(10,227) xval(jk),eigen(je,jk)
270 continue
227 format(5f10.5)
!
call psend
write(6,*) ' plotcl done'
stop 'SPAGH END'
!
!.....ERROR HANDLING FOR FILES
!
! 7 write(6,*) 'cannot open LAPW1 output file'
! write(6,*) 'filename=',ofname
! stop 'end NOT ok'
!
! 20 write(6,*) 'cannot open STRUCT file'
! write(6,*) 'filename=',sfname
! stop 'end NOT ok'
!
! 8 write(6,*) 'cannot open VIEWING-parameter file'
! write(6,*) 'filename=',ifname
! stop 'end NOT ok'
140 format(f5.2,' setlinewidth')
900 STOP 'GTFNAM - Exactly one commandline argument has to be given.'
end
-------------- next part --------------
module irr_param
dimension skir(3),iprto(3)
real*8,pointer :: eneik(:,:) !NEVL,NKP)
integer,pointer :: ngrp(:,:,:), ngde(:,:,:) !(4,NEVL,NKP)
integer,pointer :: ndeg(:,:) !NEVL,NKP)
integer,pointer :: neik(:),ikg(:),ntab4(:) !NKP)
integer,pointer :: iscol(:),iscol2(:),iste(:),ibnd(:) !NEVL)
integer,pointer :: lkg(:,:,:) !48,48,NKP)
integer,pointer :: nlkg(:,:) !48,NKP)
character*12,pointer :: grpnam(:) !NKP)
character*12,pointer :: cnam(:,:) !48,NKP)
character*11 col
character*3 tjunk
integer iprtf,ifont
real*8 dlwid,dawid,dfwid
end module irr_param
More information about the Wien
mailing list