c version of 23.5.05
c program to calculate results of Fremlin scheme 2
program prrep
dimension store(100000), partyn(300), partyv(300), text(1000),
. arrayp(660,20), arrayv(660,20),
. partyq(300), elect0(660), elect2(660), list(300),
. jlist(300), kist(300), prior(300), ncand(660),
. constn(660), parmp0(300), parmp2(300), textt(1000)
integer partyn, partyv, arrayp, arrayv, partyq, elect0, elect2,
. parmp0, parmp2, constn
character*1 store, text, textt
nstorx = 100000
npartx = 300
ntx = 1000
ncandx = 20
nstor = 0
nparty = 0
do 12 nc = 1,660
do 12 j = 1,20
arrayp(nc,j) = 0
arrayv(nc,j) = 0
12 continue
call prrd1(store,nstor,nstorx,partyn,nparty,npartx,
. nconst,arrayp,arrayv,ncandx,constn)
nvtot = 0
do 22 np = 1,nparty
22 partyv(np) = 0
do 25 nc = 1,nconst
do 24 j = 1,ncandx
np = arrayp(nc,j)
if (np.eq.0) goto 24
ncand(nc) = j
nv = arrayv(nc,j)
partyv(np) = partyv(np) + nv
nvtot = nvtot + nv
24 continue
25 continue
c do 31 np = 1,nparty
c strpr = float(partyv(np))*float(nconst)/float(nvtot)
c if (strpr.lt.1.) goto 31
c nt = 0
c call strget(store,nstorx,text,ntx,partyn(np),nt)
c write (*,302) partyv(np), strpr, (text(i),i=1,nt)
c302 format(i10,f7.1,2x,1000a1)
c31 continue
c pause
call scheme0(arrayv,ncand,nconst,elect0)
call scheme2(arrayv,arrayp,elect2,nconst,nparty,npartx,
. ncandx,partyv,nvtot,ncand,partyn,store,nstorx)
c statistics
4 do 42 np = 1,nparty
parmp0(np) = 0
42 parmp2(np) = 0
do 45 nc = 1,nconst
j0 = elect0(nc)
np0 = arrayp(nc,j0)
parmp0(np0) = parmp0(np0) + 1
j2 = elect2(nc)
np2 = arrayp(nc,j2)
parmp2(np2) = parmp2(np2) + 1
45 continue
5 call sorti(jlist,partyv,nparty,kist,npartx,npartx)
do 51 i51 = 1,nparty
np = jlist(i51)
mp0 = parmp0(np)
mp2 = parmp2(np)
strpr = float(partyv(np))*float(nconst)/float(nvtot)
if ((mp0.eq.0).and.(mp2.eq.0).and.(strpr.lt.1.)) goto 51
nt = 0
call strget(store,nstorx,text,ntx,partyn(np),nt)
k = 15 - nt
if (k.gt.0) call strbl1(text,nt,ntx,k)
write (*,502) (text(i),i=1,15), partyv(np), strpr, mp0, mp2
502 format(2x,15a1,i10,f7.1,2i5,2x)
51 continue
c satisfaction
nsatv0 = 0
nsatv2 = 0
do 55 nc = 1,nconst
do 54 j = 1,ncand(nc)
nv = arrayv(nc,j)
if (elect0(nc).eq.j) nsatv0 = nsatv0 + nv
if (elect2(nc).eq.j) nsatv2 = nsatv2 + nv
54 continue
55 continue
satv0 = 100.*float(nsatv0)/float(nvtot)
satv2 = 100.*float(nsatv2)/float(nvtot)
write (*,56) satv0, satv2
56 format(1x,'Satisfaction percentages',2f5.1)
c output of results
call output(elect0,elect2,ncand,arrayv)
call out2(elect0,elect2,arrayp,store,nstorx,partyn,npartx,
. constn, nconst)
end
subroutine prrd1(store,nstor,nstorx,partyn,nparty,npartx,
. nconst,arrayp,arrayv,ncandx,constn)
dimension store(nstorx), partyn(npartx),
. text(1000), text1(1000), textp(1000), textc(1000), textc0(1000),
. arrayp(660,20), arrayv(660,20), constn(660), tindep(11)
integer partyn, arrayp, arrayv, constn
character*1 store, text, t, text1, textp, textc, textc0, tindep
data tindep/'I','n','d','e','p','e','n','d','e','n','t'/
ntx = 1000
ntc0 = 0
nc = 0
open (unit=3,file='prrep.in',status='old')
2 read (3,21,end=8) (text(i),i=1,ntx)
21 format(1000a1)
call stripc(text,ntx,nt)
if (nt.eq.0) goto 2
c write (*,21987) (text(i),i=1,nt)
c21987 format(1000a1)
call strsg1(text,nt,ntx,textc,ntc,ntx)
c Constituency_name
if (ntc.ne.ntc0) goto 225
do 222 i = 1,ntc
if (textc(i).ne.textc0(i)) goto 225
222 continue
goto 23
225 nc = nc + 1
do 226 i = 1,ntc
226 textc0(i) = textc(i)
ntc0 = ntc
j = 0
call strstr(textc,ntc,ntx,store,nstor,nstorx,constn(nc),iflag)
23 if (j.ge.ncandx) goto 93
j = j + 1
call strsg1(text,nt,ntx,text1,nt1,ntx)
c Region
call strsg1(text,nt,ntx,text1,nt1,ntx)
c Electorate
call strsg1(text,nt,ntx,text1,nt1,ntx)
c Last_name
call strsg1(text,nt,ntx,textp,ntp,ntx)
c Party_name
if (ntp.ne.11) goto 24
do 235 i = 1,11
if (textp(i).ne.tindep(i)) goto 24
235 continue
goto 242
24 call strsrh(textp,ntp,ntx,store,nstorx,partyn,npartx,np)
if (np.ne.0) goto 25
242 call strstr(textp,ntp,ntx,store,nstor,nstorx,n,iflag)
if (n.eq.-1) goto 91
if (nparty.eq.npartx) goto 92
nparty = nparty + 1
partyn(nparty) = n
np = nparty
25 arrayp(nc,j) = np
call strsg1(text,nt,ntx,text1,nt1,ntx)
c 2005_votes
call strti1(text1,nt1,0,nt11,nv)
c write (*,28987) nc, np, nv, (textp(i),i=1,ntp),
c . (textc0(i),i=1,ntc0)
c28987 format(i4,i4,i10,2x,1000a1)
arrayv(nc,j) = nv
goto 2
8 continue
nconst = nc
goto 99
91 write (*,911) nstorx, (text1(i),i=1,nt1)
911 format(/,'[PRRD1] nstorx =',i8,' too small;',
. /,'party name',2x,1000a1)
goto 99
92 write (*,921) npartx, (text1(i),i=1,nt1)
921 format(/,'[PRRD1] npartx =',i8,' too small;',
. /,'party name',2x,1000a1)
goto 99
93 write (*,931) ncandx, (text(i),i=1,nt1)
931 format(/,'[PRRD1] ncandx =',i8,' too small;',
. /,'text is',2x,1000a1)
goto 99
99 continue
close (3)
return
end
subroutine scheme2(arrayv,arrayp,elect,nconst,nparty,npartx,
. ncandx,partyv,nvtot,ncand,partyn,store,nstorx)
dimension arrayp(660,20), arrayv(660,20), elect(660),
. partyv(npartx), partyq(300), parmp(300), ncand(660), list(300),
. jlist(300), kist(300), arrayx(660,20), partyn(npartx),
. store(nstorx)
integer elect, arrayp, arrayv, partyv, partyq, parmp, partyn
character*1 store
mdc = -987654321
do 22 np = 1,nparty
parmp(np) = 0
partyq(np) = float(nconst)*float(partyv(np))/float(nvtot)
22 continue
c check for eligibility of candidates on 20% rule
c elect1(nc) = j if elected candidate in constituency nc
c is from party arrayp(nc,j) with vote arrayv(nc.j)
c parmp(np) = number of candidates so far elected from party np
do 39 nc = 1,nconst
nvc = 0
do 31 j = 1,ncand(nc)
nv = arrayv(nc,j)
list(j) = nv
nvc = nvc + nv
31 continue
call sorti(jlist,list,ncand(nc),kist,npartx,npartx)
j1 = jlist(1)
j2 = jlist(2)
nv1 = arrayv(nc,j1)
nv2 = arrayv(nc,j2)
np1 = arrayp(nc,j1)
iflelg = 0
do 35 j = 1,ncand(nc)
nv = arrayv(nc,j)
if (nv.ge.(nvc/5)) goto 33
arrayx(nc,j) = mdc
goto 35
33 continue
c nv = nv - (nv1+nv2)/2
if (j.eq.j1) nv = nv - nv2
if (j.ne.j1) nv = nv - nv1
iflelg = iflelg + 1
34 arrayx(nc,j) = nv
35 continue
if (iflelg.ge.2) goto 39
elect(nc) = j1
parmp(np1) = parmp(np1) + 1
39 continue
c basic algorithm begins
call state(nparty,parmp,npartx,store,nstorx,partyn)
5 p0 = 0.
np0 = 0
do 52 np = 1,nparty
if (parmp(np).ge.partyq(np)) p = 0.
if ((parmp(np).eq.0).and.(partyq(np).gt.0))
. p = (float(nconst)+1.)*float(partyv(np))
if ((parmp(np).gt.0).and.(partyq(np).gt.parmp(np)))
. p = float(partyv(np))/float(parmp(np))
if (p.le.p0) goto 52
p0 = p
np0 = np
52 continue
if (p0.le.0.) goto 7
c party np0 is the party next due for an MP
nc0 = 0
j0 = 0
q0 = mdc
do 55 nc = 1,nconst
if (elect(nc).ne.0) goto 55
do 54 j = 1,ncand(nc)
if (arrayp(nc,j).ne.np0) goto 54
if (arrayx(nc,j).le.q0) goto 54
nc0 = nc
q0 = arrayx(nc,j)
j0 = j
54 continue
55 continue
if (q0.gt.mdc) goto 6
c party np0 has no electable candidates; reduce its quota
nt = 0
partyq(np0) = parmp(np0)
goto 5
6 elect(nc0) = j0
parmp(np0) = parmp(np0) + 1
goto 5
7 call state(nparty,parmp,npartx,store,nstorx,partyn)
continue
c end of phase 2
c now fill in with first-past the post
do 75 nc = 1,nconst
if (elect(nc).ne.0) goto 75
nv0 = 0
j0 = 0
do 73 j = 1,ncand(nc)
nv = arrayv(nc,j)
if (nv.le.nv0) goto 73
nv0 = nv
j0 = j
73 continue
elect(nc) = j0
75 continue
return
end
subroutine state(nparty,parmp,npartx,store,nstorx,partyn)
dimension parmp(npartx), partyn(npartx), text(1000)
integer parmp, partyn
character*1 store, text
ntx = 1000
write (*,15)
15 format(1x,'State of parties')
do 3 np = 1,nparty
if (parmp(np).eq.0) goto 3
nt = 0
call strget(store,nstorx,text,ntx,partyn(np),nt)
write (*,22) parmp(np), (text(i),i=1,nt)
22 format(i6,2x,1000a1)
3 continue
return
end
subroutine scheme0(arrayv,ncand,nconst,elect)
dimension arrayv(660,20), ncand(660), elect(660)
integer arrayv, elect
do 3 nc = 1,nconst
j0 = 0
nv0 = 0
do 25 j = 1,ncand(nc)
nv = arrayv(nc,j)
if (nv.le.nv0) goto 25
j0 = j
nv0 = nv
25 continue
elect(nc) = j0
3 continue
return
end
c output in html of list of successful candidates
subroutine output(elect0,elect2,ncand,arrayv)
dimension elect0(660), elect2(660), text(1000), text1(1000),
. textp(1000), textc(1000), textc0(1000), textsn(1000),
. text0(1000), text2(1000), ncand(660), arrayv(660,20),
. cellbd(9)
integer elect0, elect2, arrayv
character*1 text, text1, textp, textc0, textsn, tfi, tmi,
. text0, text2, textc, cellbd
data cellbd/'<','/','t','d','>','<','t','d','>'/
ntx = 1000
ntc0 = 0
nc = 0
open (unit=3,file='prrep.in',status='old')
open (unit=7,file='prrep.out.htm')
read (3,21) (text(i),i=1,ntx)
call strip(text,ntx,nt)
write (7,18) (text(i),i=1,nt)
c write (*,18) (text(i),i=1,nt)
18 format(1000a1)
2 read (3,21,end=8) (text(i),i=1,ntx)
21 format(1000a1)
call stripc(text,ntx,nt)
if (nt.eq.0) goto 2
call strsg1(text,nt,ntx,textc,ntc,ntx)
c Constituency_name
if (ntc.ne.ntc0) goto 225
do 222 i = 1,ntc
if (textc(i).ne.textc0(i)) goto 225
222 continue
goto 23
225 nc = nc + 1
do 226 i = 1,ntc
226 textc0(i) = textc(i)
ntc0 = ntc
j = 0
nt0 = 0
nt2 = 0
23 j = j + 1
if ((elect0(nc).ne.j).and.(elect2(nc).ne.j)) goto 2
call strsg1(text,nt,ntx,text1,nt1,ntx)
c Region
call strsg1(text,nt,ntx,text1,nt1,ntx)
c Electorate
call strsg1(text,nt,ntx,textsn,ntsn,ntx)
c Candidate name
call strjoi(textsn,ntsn,ntx,cellbd,9,9,iflag)
if (elect0(nc).eq.j)
. call strjoi(text0,nt0,ntx,textsn,ntsn,ntx,iflag)
if (elect2(nc).eq.j)
. call strjoi(text2,nt2,ntx,textsn,ntsn,ntx,iflag)
call strsg1(text,nt,ntx,textp,ntp,ntx)
c Party_name
call strjoi(textp,ntp,ntx,cellbd,9,9,iflag)
if (elect0(nc).eq.j)
. call strjoi(text0,nt0,ntx,textp,ntp,ntx,iflag)
if (elect2(nc).eq.j)
. call strjoi(text2,nt2,ntx,textp,ntp,ntx,iflag)
call strsg1(text,nt,ntx,text1,nt1,ntx)
c 2005_votes
call strjoi(text1,nt1,ntx,cellbd,9,9,iflag)
if (elect0(nc).eq.j)
. call strjoi(text0,nt0,ntx,text1,nt1,ntx,iflag)
if (elect2(nc).eq.j)
. call strjoi(text2,nt2,ntx,text1,nt1,ntx,iflag)
26 if ((j.lt.elect0(nc)).or.(j.lt.elect2(nc))) goto 2
if (elect0(nc).eq.elect2(nc)) goto 28
call strjoi(text0,nt0,ntx,text2,nt2,ntx,iflag)
c calculate rank
j2 = elect2(nc)
nv2 = arrayv(nc,j2)
nr = 1
do 24 i = 1,ncand(nc)
if (arrayv(nc,i).gt.nv2) nr = nr + 1
24 continue
call strit1(nr,0,nfill,text0,ntx,nt0)
28 do 2802 i = 1,ntc0
if (textc0(i).eq.'_') textc0(i) = ' '
2802 continue
do 2803 i = 1,nt0
if (text0(i).eq.'_') text0(i) = ' '
2803 continue
write (7,29) (textc0(i),i=1,ntc0),(cellbd(i),i=1,9),
. (text0(i),i=1,nt0),
. '<','/','t','d','>','<','/','t','r','>'
c write (*,29) (textc0(i),i=1,ntc0),' ',' ',(text0(i),i=1,nt0),
c . '<','/','t','d','>','<','/','t','r','>'
29 format('
',1000a1)
goto 2
8 continue
close (3)
close (7)
return
end
c output of list of successful parties
subroutine out2(elect0,elect2,arrayp,store,nstorx,partyn,npartx,
. constn, nconst)
dimension elect0(660), elect2(660), arrayp(660,20),
. store(nstorx), partyn(npartx), constn(660), text(1000)
integer elect0, elect2, arrayp, partyn, constn
character*1 store, text
ntx = 1000
open (unit=7,file='prrep.out2')
open (unit=3,file='prrep.in')
read (3,21) (text(i),i=1,ntx)
21 format(1000a1)
call strip(text,ntx,nt)
write (7,21) (text(i),i=1,nt)
do 5 nc = 1, nconst
nt = 0
call strget(store,nstorx,text,ntx,constn(nc),nt)
call strbl1(text,nt,ntx,2)
j = elect0(nc)
np = arrayp(nc,j)
call strget(store,nstorx,text,ntx,partyn(np),nt)
call strbl1(text,nt,ntx,2)
j = elect2(nc)
np = arrayp(nc,j)
call strget(store,nstorx,text,ntx,partyn(np),nt)
c write (*,34) (text(i),i=1,nt)
write (7,34) (text(i),i=1,nt)
34 format(1000a1)
5 continue
close (3)
close (7)
return
end
c subroutine to delete blanks from the beginning of a text string text(1)-
c text(ntx); returns nt as position of last non-blank character
c (nt=0 if the whole string blank).
c The ASCII character 13 is ignored if found at the end of a line,
c as it is presumed to have been generated by an MSWord editor
subroutine strip(text,ntx,nt)
dimension text(ntx)
character*1 text, t, tx
tx = char(13)
nt = 0
do 1 k=1,ntx
if (text(k).ne.' ') goto 2
1 continue
goto 9
2 k=k-1
do 3 i=1,ntx-k
t = text(i+k)
text(i) = t
if ((t.ne.' ').and.(t.ne.tx)) nt = i
3 continue
if (k.eq.0) goto 9
do 4 i = ntx-k+1,ntx
4 text(i) = ' '
9 return
end
c subroutine to delete blanks from the beginning of a text string text(1)-
c text(ntx); stops at %;
c returns nt as position of last non-blank character
c (nt=0 if the whole string blank).
subroutine stripc(text,ntx,nt)
dimension text(ntx)
character*1 text, t
nt = 0
do 1 k=1,ntx
t = text(k)
if (t.eq.'%') goto 9
if (t.ne.' ') goto 2
1 continue
goto 9
2 k=k-1
do 3 i=1,ntx-k
t = text(i+k)
if (t.eq.'%') goto 9
text(i) = t
if (t.ne.' ') nt = i
3 continue
9 return
end
c subroutine to read a text string as a string of integers.
subroutine strtis(text,ntxtx,list,nlistx,nl)
dimension text(ntxtx), list(nlistx)
character*1 text
nl = 0
call strlen(text,ntxtx,nt)
if (nt.eq.0) goto 9
nt0 = 0
2 if (nl.ge.nlistx) goto 8
nl = nl+1
call strti1(text,ntxtx,nt0,nt1,list(nl))
if (nt1.ge.nt) goto 9
nt0 = nt1
goto 2
8 write (6,81) nlistx, (text(i),i=1,nt)
81 format(/,1x,'[STRTIS] nlistx = ',i7,' which is too small;',
1 /,'text string is ',10000a1)
goto 9
9 continue
return
end
c subroutine to read a text string as an integer. It
c starts at text(nt0+1) and keeps going until it thinks
c the integer has finished. Initial blanks are passed over.
c nt1 is the location in text() of the end of the integer.
subroutine strti1(text,ntxtx,nt0,nt1,i)
dimension text(ntxtx)
character*1 text, t
c write (6,2987) (text(i),i=1,ntxtx)
c2987 format(2x,'[STRTI1]',2x,50a1)
i = 0
ksgn = 1
nt1 = nt0
if (nt0.ge.ntxtx) goto 40
n0 = nt0
if (n0.lt.0) n0 = 0
mode = 1
do 20 j = n0+1,ntxtx
t = text(j)
if (t.eq.' ') goto (20,30) mode
if (t.eq.'-') goto (5,25) mode
m = ichar(t)
if ((m.lt.48).or.(m.gt.57)) goto 25
goto (10,15) mode
5 ksgn = -1
mode = 2
goto 20
10 mode = 2
15 i = 10*i + m - 48
nt1 = j
20 continue
goto (25,30) mode
25 write (6,27) (text(k),k=n0+1,j)
27 format(/,1x,'[STRTI1] query: string reads',/,5x,1000a1)
30 i = i *ksgn
40 return
end
subroutine strstr(text,ntxt,ntxtx,store,nstor,nstorx,
1 n,iflag)
c this stores text(1)-text(ntxt) in store(nstor+1)-store(nstor+ntxt);
c returning nstor augmented by ntxt , and
c n = nstor + ntxt*nstorx . If the length nstorx of store() is
c insufficient it does nothing and returns iflag.and.1 = 1 and
c n = -1 .
dimension store(nstorx), text(ntxtx)
character*1 store, text
if (ntxt.eq.0) goto 13
if ((nstor+ntxt).gt.nstorx) goto 15
do 5 i = 1,ntxt
5 store(nstor+i) = text(i)
10 n = nstor + nstorx*ntxt
nstor = nstor + ntxt
goto 25
13 n = 0
goto 25
15 n = -1
iflby2 = iflag/2
iflp1 = iflag - 2*iflby2
if (iflp1.ne.0) goto 25
iflag = iflag + 1
write (6,17) nstorx
17 format(/,1x,'[STRSTR] nstorx =',i6,' too small; string is')
write (5,18) (text(i),i=1,ntxt)
18 format(1x,10000a1)
25 return
end
subroutine strsrh(text,ntxt,ntxtx,store,nstorx,index,nindx,n)
c this subroutine looks for the string text(1)-text(ntxt) in the
c store() array, restricting itself to regions indexed
c by index() . It returns n where
c the string appears in the region coded by index(n) ; if there
c isn't one, it returns n = 0 .
dimension text(ntxtx), store(nstorx), index(nindx)
character*1 text, store
n = 0
do 10 l = 1,nindx
il = index(l)
if (il.lt.0) goto 10
k = il/nstorx
if (k.ne.ntxt) goto 10
if (k.eq.0) goto 15
i0 = il - k*nstorx
if ((i0+k).gt.nstorx) goto 10
do 5 i = 1,k
if (store(i0+i).ne.text(i)) goto 10
5 continue
goto 15
10 continue
goto 20
15 n = l
20 return
end
subroutine strjoi(text1,nt1,ntxtx1,text2,nt2,ntxtx2,iflag)
c this subroutine adds the string text2(1)-text2(nt2) onto the
c end of the string text1(1)-text1(nt1) . If ntxtx1 is too
c small it returns iflag.and.2 = 2 and just adds the new string
c up to the end of text1() .
dimension text1(ntxtx1), text2(ntxtx2)
character*1 text1, text2
m = ntxtx1 - nt1
if (m.gt.nt2) m = nt2
if (m.eq.0) goto 15
do 10 i = 1,m
10 text1(nt1+i) = text2(i)
15 continue
iflby2 = iflag/2
iflby4 = iflby2/2
iflp2 = 2*iflby2 - 4*iflby4
if ((m.eq.nt2).or.(iflp2.eq.2)) goto 25
iflag = iflag + 2
write (6,17) ntxtx1
17 format(/,1x,'[STRJOI] ntxtx1 =',i6,' too small; strings are')
write (6,19) (text1(i),i=1,nt1)
write (6,19) (text2(i),i=1,nt2)
19 format(1x,10000a1)
25 nt1 = nt1 + m
return
end
subroutine strget(store,nstorx,text,ntxtx,n,ntxt)
c this subroutine looks for the string in store() with code n and
c puts it into text(ntxt+1)-... . If the string is too long it
c gets truncated. If there is no string with this code then
c ?? is returned.
dimension store(nstorx), text(ntxtx)
character*1 store, text
m = n/nstorx
i0 = n - nstorx*m
if ((n.lt.0).or.((i0+m).gt.nstorx)) goto 15
if (m.le.(ntxtx-ntxt)) goto 8
write (6,4) ntxtx,(store(i0+i),i=1,m)
4 format(/,1x,'[STRGET] ntxtx =',i6,' too small; string is',
1 /,1x,10000a1)
m = ntxtx - ntxt
8 if (m.eq.0) goto 20
do 10 i = 1,m
10 text(i+ntxt) = store(i0 + i)
ntxt = ntxt + m
goto 20
15 do 17 i = 1,2
17 call strch1(text,ntxt,ntxtx,'?')
20 return
end
subroutine strlen(text,ntxtx,ntxt)
c this subroutine looks for the last non-blank character in
c text()
dimension text(ntxtx)
character*1 text
do 5 ntxt = ntxtx,1,-1
5 if (text(ntxt).ne.' ') goto 10
ntxt = 0
10 return
end
subroutine strsg1(text1,nt1,nt1x,text2,nt2,nt2x)
c this subroutine takes the string text1(1)-text1(nt1) and
c looks for (1) the first "word" (i.e. string without blanks),
c fitting this into text2(1)-text2(nt2) (2) the rest of the
c string (ignoring blanks at beginning), fitting this
c into text1(1)-text1(nt1). If text2() isn't long enough
c then the first word is truncated
dimension text1(nt1x), text2(nt2x)
character*1 t, text1, text2
mode = 1
nt2 = 0
nt = 0
do 20 i = 1,nt1
t = text1(i)
if (t.eq.' ') goto (20,6,20,10) mode
goto (2,4,8,10) mode
2 mode = 2
4 nt2 = nt2 + 1
if (nt2.le.nt2x) text2(nt2) = t
goto 20
6 mode = 3
goto 20
8 mode = 4
10 nt = nt + 1
text1(nt) = t
20 continue
if (nt2.le.nt2x) goto 30
write (6,22) nt2x,(text1(i),i=1,nt1)
22 format(/,1x,'[STRSG1] nt2x =',i6,' too small; string is',
1 /,1x,10000a1)
nt2 = nt2x
30 nt1 = nt
return
end
subroutine strch1(text,ntxt,ntxtx,t)
c this subroutine adds a single term to text()
dimension text(ntxtx)
character*1 t, text
if (ntxt.lt.ntxtx) goto 10
write (6,5) ntxtx, (text(i),i=1,ntxtx),t
5 format(/,1x,'[STRCH1] ntxtx =',i6,' too small; string is',
1 /,1x,10000a1)
goto 20
10 ntxt = ntxt + 1
text(ntxt) = t
20 return
end
subroutine strwr1(text,ntxt,ntxtx,nch)
c this subroutine sends a line to text to unit nch ; if nch = -1
c it does nothing; if nch = 0 it types the line.
dimension text(ntxtx)
character*1 text
if (ntxt.eq.0) goto 15
if (nch.gt.0) write(nch,10) (text(i),i=1,ntxt)
if (nch.eq.0) write (6,11) (text(i),i=1,ntxt)
10 format(10000a1)
11 format(1x,1000a1)
goto 25
15 if (nch.gt.0) write(nch,17)
if (nch.eq.0) write (6,17)
17 format(1x)
25 return
end
subroutine strbl1(text,ntxt,ntxtx,k)
c this subroutine adds k blanks to the end of the string
c text(1)-text(ntxt)
dimension text(ntxtx)
character*1 text
c write (*,2) ntxtx, ntxt, (text(i), i=1,ntxt)
c2 format('[STRBL1] input:',2i6,2x,255a1)
nt = ntxt + k
if (nt.gt.ntxtx) nt = ntxtx
if (nt.le.ntxt) goto 10
do 5 i = ntxt+1,nt
5 text(i) = ' '
ntxt = nt
10 continue
c write (*,11) ntxtx, ntxt, (text(i), i=1,ntxt)
c11 format('[STRBL1] output:',2i6,2x,255a1)
return
end
C THIS SUBROUTINE GIVES A NON-STANDARD FORTRAN OUTPUT FOR DECIMAL
C NUMBERS. IN THE INPUT
C strft1(X0,NPLACE,NFIELD,NUSED,text,ntxtx,ntxt)
C X0 REPRESENTS THE NUMBER TO BE WRITTEN, NPLACE GIVES THE DEGREE OF
C ACCURACY, NFIELD THE OUTPUT FIELD WIDTH, & NUSED REPORTS BACK TO
C THE CALLING PROGRAM. The output is located in text(ntxt+1)-
c with ntxt augmented at the end; it's truncated if text() isn't
c long enough.
C IF NFIELD > 0 THE OUTPUT IS RIGHT-JUSTIFIED IN A FIELD OF WIDTH
C NFIELD; IF NFIELD < 0 IT IS LEFT-JUSTIFIED IN A FIELD OF WIDTH
C -NFIELD; IF NFIELD = 0 THE OUTPUT IS GIVEN IN WHATEVER FIELD WIDTH
C IS NECESSARY.
C IF NPLACE > 1000 THE OUTPUT IS GIVEN TO NPLACE-1000 SIGNINFICANT
C FIGURES; IF NPLACE <= 1000 IT IS GIVEN TO NPLACE DECIMAL PLACES
C (-VE NPLACE PERMITTED). IF THE GIVEN FIELD WIDTH IS TOO SMALL THE
C PROGRAM DOES ITS BEST, REVERTING TO ** ONLY IN EXTREMIS. IT WILL
C USE AN E FORMAT IF THIS SAVES SPACE.
C IF NUSED IS RETURNED > 0, THIS IS THE NUMBER OF SPACES SPENT IN
C OUTPUTTING THE NUMBER X0. IF NUSED IS RETURNED -1, THEN THERE WAS
C NOT ENOUGH SPACE BUT SOMETHIN WAS TRIED; -2 MEANS THAT ** WAS
C RESORTED TO.
SUBROUTINE strft1(X0,NPLACE,NFIELD,NUSED,text,ntxtx,ntxt)
DIMENSION LIST(100), text(ntxtx)
character*1 list, text
DOUBLE PRECISION X,X1,X2
NTRY = 1
3 X1 = X0
NPL = NPLACE
NFD = NFIELD
IF (NFD.LT.0) NFD = -NFD
IF (NFD.EQ.0) NFD = 100
IF (NFD.LT.2) GOTO 220
NUSED = NFD
I1 = 0
IF (X0) 15,5,20
5 I0 = 2
LIST(1) = '0'
LIST(2) = '.'
IF (NPLACE.GT.1000) NPL = NPL-1001
IF (NPL.GT.(NFD-I0)) NPL = NFD - I0
NSIG = 0
IF (NPL.LE.0) GOTO 240
DO 10 I = 1,NPL
I0 = I0 + 1
10 LIST(I0) = '0'
GOTO 240
15 I1 = 1
LIST(1) = '-'
X1 = -X0
20 M1 = DLOG10(X1) + 1000.
M1 = M1 - 1000
X2 = X1/(10.**M1)
IF (X2.GE.1.) GOTO 22
M1 = M1 - 1
X2 = X1/(10.**M1)
22 KFLAG = 0
MX = M1
X = X2
I0 = I1
NSIG = MX + NPL + 1
IF (NPLACE.GT.1000) NSIG = NPL - 1000
IF (NSIG.GT.8) NSIG = 8
IF (NSIG) 5,24,26
24 IF (X.LE.5.) GOTO 5
MX = MX + 1
X = X/10.
NSIG = 1
26 X = X + 5./(10.**NSIG)
IF (X.LT.10.) GOTO 30
X = X/10.
MX = MX + 1
IF ((NPLACE.LE.1000).OR.((NTRY.GE.3).AND.((MX.EQ.-9).
1 OR.((MX.GT.-5).AND.(MX.LT.0))))) NSIG = NSIG + 1
30 IF (MX.LE.-5) GOTO 33
IF ((MX-NSIG).LT.2) GOTO 60
IF ((NTRY.EQ.1).AND.((MX-NSIG).EQ.2)) GOTO 55
33 DO 35 J = 1,NSIG
L = X
X = (X-DFLOAT(L))*10.
I0 = I0 + 1
LIST(I0) = char(l+48)
IF (J.NE.1) GOTO 35
I0 = I0 + 1
LIST(I0) = '.'
35 CONTINUE
I0 = I0 + 1
LIST(I0) = 'E'
IF (MX.GT.0) GOTO 40
I0 = I0 + 1
LIST(I0) = '-'
MX = -MX
40 IF (MX.LE.9) GOTO 45
L = MX/10
MX = MX - 10*L
I0 = I0 + 1
LIST(I0) = char(l+48)
45 I0 = I0 + 1
LIST(I0) = char(mx+48)
GOTO 200
55 KFLAG = 1
60 IF (MX.GE.0) GOTO 70
IF (NTRY.GT.1) GOTO 62
KFLAG = 1
I0 = I0 + 1
LIST(I0) = '0'
62 I0 = I0 + 1
LIST(I0) = '.'
IF (MX.EQ.-1) GOTO 70
DO 65 J = MX,-2
I0 = I0 + 1
65 LIST(I0) = '0'
70 DO 75 J = 1,NSIG
L = X
X = (X-DFLOAT(L))*10.
I0 = I0 + 1
LIST(I0) = char(l+48)
IF (MX.NE.0) GOTO 75
I0 = I0 + 1
LIST(I0) = '.'
75 MX = MX - 1
IF (MX.LT.0) GOTO 200
DO 80 J = 1,MX+1
I0 = I0 + 1
80 LIST(I0) = '0'
I0 = I0 + 1
LIST(I0) = '.'
200 IF (NFD.GE.I0) GOTO 240
IF (NSIG.LE.1) GOTO 220
NUSED = -1
IF (NTRY.GT.1) GOTO 203
201 NTRY = 2
IF (KFLAG.EQ.1) GOTO 22
203 NPL = NPL - 1
NTRY = NTRY + 1
GOTO 22
220 DO 225 J = I1+1,NFD
225 LIST(J) = '*'
NUSED = -2
GOTO 300
240 IF (I0.EQ.NFD) GOTO 300
IF (NFIELD) 245,280,260
245 DO 250 J = I0 + 1,NFD
250 LIST(J) = ' '
GOTO 300
260 I = NFD - I0
DO 265 J = NFD,I+1,-1
265 LIST(J) = LIST(J-I)
DO 270 J = 1,I
270 LIST(J) = ' '
GOTO 300
280 NFD = I0
NUSED = I0
300 DO 305 I = 1,NFD
ntxt = ntxt + 1
if (ntxt.le.ntxtx) text(ntxt) = list(i)
305 CONTINUE
if (ntxt.le.ntxtx) goto 320
write (6,308) ntxtx,(list(i),i=1,nfd)
308 format(/,1x,'[STRFT1] ntxtx =',i6,' too small; number is',
1 /,1x,1000a1)
ntxt = ntxtx
320 RETURN
END
C THIS SUBROUTINE GIVES A FLEXIBLE-FORMAT OUTPUT FOR INTEGERS.
C INPUT (INTEG,NSPACE,NFILL) WHERE: INTEG IS THE INTEGER TO BE
C TYPED; NSPACE IS THE NUMBER OF SPACES AVAILABLE (+VE FOR RIGHT-
C JUSTIFIED, -VE FOR LEFT-JUSTIFIED, 0 FOR EXACT NUMBER NEEDED) AND
C NFILL IS RETURNED AS THE NUMBER N OF SPACES USED OR -N IF
C THERE WASN'T ENOUGH ROOM. Output is added to the
c array text() , truncated if ntxtx is too small.
subroutine strit1(integ,nspace,nfill,text,ntxtx,ntxt)
dimension text(ntxtx)
character*1 t, text
NREC = 1
KQ = 1
I = INTEG
IF (I.LT.0) I = -I
10 IF (I.LT.10) GOTO 15
NREC = NREC + 1
KQ = KQ*10
I = I/10
GOTO 10
15 IF (INTEG.LT.0) NREC = NREC + 1
NFILL = NSPACE
IF (NSPACE.LT.0) NFILL = -NSPACE
IF (NSPACE.EQ.0) NFILL = NREC
IF (NREC.GT.NFILL) GOTO 90
IF (NSPACE.LE.0) GOTO 30
IF (NREC.EQ.NFILL) GOTO 30
DO 25 K = 1,NFILL-NREC
ntxt = ntxt + 1
25 if (ntxt.le.ntxtx) text(ntxt) = ' '
30 I = INTEG
IF (I.GE.0) GOTO 35
I = -I
ntxt = ntxt + 1
if (ntxt.le.ntxtx) text(ntxt) = '-'
32 FORMAT(1H+,'-',$)
35 J = I/KQ
ntxt = ntxt + 1
if (ntxt.le.ntxtx) t = char(j+48)
text(ntxt) = t
37 FORMAT(1H+,I1,$)
IF (KQ.EQ.1) GOTO 50
I = I - J*KQ
KQ = KQ/10
GOTO 35
50 IF ((NSPACE.GE.0).OR.(NREC.EQ.NFILL)) GOTO 100
DO 60 K = 1,NFILL-NREC
ntxt = ntxt + 1
60 if (ntxt.le.ntxtx) text(ntxt) = ' '
GOTO 100
90 DO 95 K = 1,NFILL
ntxt = ntxt + 1
95 if (ntxt.le.ntxtx) text(ntxt) = '*'
97 FORMAT(1H+,'*',$)
NFILL = -NFILL
100 CONTINUE
if (ntxt.le.ntxtx) goto 120
write (6,107) ntxtx,integ
107 format(/,1x,'[STRIT1] ntxtx =',i6,' too small; number is',
. /,1x,i25)
ntxt = ntxtx
120 RETURN
END
c subroutine to take a section from text(k0+1)-text(k1)
c and put it into taux(l0+1)-taux(l1) . If there's not
c room the section is truncated.
subroutine strshi(text,ntxtx,k0,k1,taux,mtxtx,l0,l1)
dimension text(ntxtx), taux(mtxtx)
character*1 text, taux
m0 = k1 - k0
m = m0
if ((m+l0).gt.mtxtx) m = mtxtx - l0
if (m.lt.0) m = 0
if (m.eq.0) goto 10
do 5 i = 1,m
5 taux(l0+i) = text(k0+i)
10 l1 = l0 + m
if (m.eq.m0) goto 15
write (6,12) mtxtx
12 format(/,1x,'[STRSHI] mtxtx =',i6,'; too small')
15 return
end
c this subroutine takes a string & tries to read it as
c a number (if possible). It returns mode = 1 if the string is
c blank; 2 if the string appears to start as an integer
c (returned as i ); 3 if the string appears to start as
c a decimal number (returned as x ); 4 if the string does not
c appear to represent a number. It will accept the notation
c 0.987E14.
subroutine strrd1(text,ntxtx,mode,x,i)
dimension text(ntxtx)
character*1 t, text
mode = 1
node = 1
ksign = 1
i = 0
x = 0
if (ntxtx.le.0) goto 200
do 50 n = 1,ntxtx
t = text(n)
if (t.eq.' ') goto (50,190,190) mode
if (t.eq.'-') goto (10,180,180) mode
if (t.eq.'.') goto (15,15,180) mode
if (t.eq.'E') goto (180,60) node
m = ichar(t)
if ((m.lt.48).or.(m.gt.57)) goto 180
node = 2
goto (20,25,30) mode
10 ksign = -1
mode = 2
goto 50
15 xset = 1.
mode = 3
goto 50
20 mode = 2
25 i = 10*i + m - 48
x = i
goto 50
30 xset = xset/10.
x = x + xset*(m-48)
goto 50
50 continue
goto 190
60 if (n.eq.ntxtx) goto 180
mode = 3
j = 0
jsign = 1
node = 1
do 80 n1 = n+1,ntxtx
t = text(n1)
if (t.eq.' ') goto (180,180,85) node
if (t.eq.'-') goto (65,180,180) node
m1 = ichar(t)
if ((m1.lt.48).or.(m1.gt.57)) goto 180
goto (70,70,75) node
65 jsign = -1
node = 2
goto 80
70 node = 3
75 j = 10*j + m1 - 48
goto 80
80 continue
85 j = j*jsign
x = x*(10.**j)
goto 190
180 mode = 4
190 x = x *ksign
i = i*ksign
200 return
end
c this subroutine takes a string text(nt0+1),...,text(ntz)
c & tries to read it
c as a number (if possible). It returns mode = 1 if the string is
c blank; 2 if the string appears to start as an integer
c (returned as i ); 3 if the string appears to start as
c a decimal number (returned as x ).
c It will accept the notation
c 0.987E14.
c The last character from the string used as part of the number is
c text(nt1).
c NOTE: the input ntz must be the last
c position in the string which is to be read (not the total length of
c the array!)
subroutine strtf1(text,ntz,nt0,nt1,mode,x,i)
dimension text(ntz)
character*1 t, text
mode = 1
node = 1
ksign = 1
i = 0
x = 0
nt1 = nt0
if (ntz.le.nt0) goto 200
do 50 n = nt0+1,ntz
t = text(n)
if (t.eq.' ') goto (50,190,190) mode
if (t.eq.'-') goto (10,180,180) mode
if (t.eq.'.') goto (15,15,180) mode
if (t.eq.'E') goto (180,60) node
m = ichar(t)
if ((m.lt.48).or.(m.gt.57)) goto 180
node = 2
goto (20,25,30) mode
10 ksign = -1
mode = 2
goto 50
15 xset = 1.
mode = 3
goto 50
20 mode = 2
25 i = 10*i + m - 48
x = i
goto 50
30 xset = xset/10.
x = x + xset*(m-48)
goto 50
50 nt1 = n
goto 190
60 if (n.eq.ntz) goto 180
mode = 3
j = 0
jsign = 1
node = 1
do 80 n1 = n+1,ntz
t = text(n1)
if (t.eq.' ') goto (180,180,85) node
if (t.eq.'-') goto (65,180,180) node
m1 = ichar(t)
if ((m1.lt.48).or.(m1.gt.57)) goto 180
goto (70,70,75) node
65 jsign = -1
node = 2
goto 80
70 node = 3
75 j = 10*j + m1 - 48
goto 80
80 nt1 = n1
85 j = j*jsign
x = x*(10.**j)
goto 190
180 continue
190 x = x *ksign
i = i*ksign
200 return
end
c this subroutine seeks a text string in store() , using hash-coding
c to accelerate the search.
subroutine strhsh(text,ntxt,ntxtx,store,nstorx,
. text1,ntx1, index, nindx, n, kh, iflag)
c The text to be looked for is in text(1)-text(ntxt) ; the hash-codes
c of the texts already stored are in khash() , and their listings are
c indexed in index() .
c The subroutine returns
c n as the position in index() corresponding to that
c string. If the string is recognised then iflag is returned as 1 ;
c if it is a new string, iflag is returned as 0 ; if an error is
c detected, iflag is returned as 2 .
c Before first use of this subroutine the array index() must be
c initialized with index(i)=-1 for all i .
dimension text(ntxtx), store(nstorx), text1(ntx1),
. index(nindx)
character*1 t,text, store, text1
kh = 0
12 if ((ntxt.gt.ntx1).or.(ntxt.gt.ntxtx)) goto 91
do 15 i = 1,ntxt
t = text(i)
j = ichar(t)
kh = 97*kh + j
15 continue
2 m = kh/nindx
n = kh - m*nindx
if (n.le.0) n = n + nindx
c write (6,2987) kh, n
c2987 format(/,i15,i10)
n0 = n
3 if (index(n).eq.-1) goto 6
nt1 = 0
call strget(store,nstorx,text1,ntx1,index(n),nt1)
if (nt1.ne.ntxt) goto 4
do 35 i = 1,ntxt
if (text(i).ne.text1(i)) goto 4
35 continue
goto 5
4 n = n + 1
if (n.gt.nindx) n = 1
if (n.eq.n0) goto 92
goto 3
5 iflag = 1
goto 99
6 iflag = 0
goto 99
91 write (6,911) ntx1, ntxt
911 format(/,'[STRHSH] error: ntx1 =',i6,' ntxtx =',i6,' but',i6,
. ' needed')
iflag = 2
goto 99
92 write (6,921) nindx
921 format(/,'[STRHSH] error: nindx =',i8,' which is too little')
iflag = 2
goto 99
99 return
end
c program to write a number in roman numerals (lower-case)
subroutine romno(text,ntx,nt,i)
dimension text(ntx)
character*1 text
if (i.le.0) goto 91
j = i
1 if (j.lt.1000) goto 2
call strch1(text,nt,ntx,'m')
j = j - 1000
goto 1
2 if (j.eq.0) goto 99
if (j.ge.4) goto 21
call strch1(text,nt,ntx,'i')
j = j - 1
goto 2
21 if (j.ge.5) goto 22
call strch1(text,nt,ntx,'i')
call strch1(text,nt,ntx,'v')
goto 99
22 if (j.ge.9) goto 23
call strch1(text,nt,ntx,'v')
j = j - 5
goto 2
23 if (j.ge.10) goto 24
call strch1(text,nt,ntx,'i')
call strch1(text,nt,ntx,'x')
goto 99
24 if (j.ge.40) goto 25
call strch1(text,nt,ntx,'x')
j = j - 10
goto 2
25 if (j.ge.49) goto 26
call strch1(text,nt,ntx,'x')
call strch1(text,nt,ntx,'l')
j = j - 40
goto 2
26 if (j.ge.50) goto 27
call strch1(text,nt,ntx,'i')
call strch1(text,nt,ntx,'l')
goto 99
27 if (j.ge.90) goto 28
call strch1(text,nt,ntx,'l')
j = j - 50
goto 2
28 if (j.ge.99) goto 29
call strch1(text,nt,ntx,'x')
call strch1(text,nt,ntx,'c')
j = j - 90
goto 2
29 if (j.ge.100) goto 30
call strch1(text,nt,ntx,'i')
call strch1(text,nt,ntx,'c')
goto 99
30 if (j.ge.400) goto 31
call strch1(text,nt,ntx,'c')
j = j - 100
goto 2
31 if (j.ge.490) goto 32
call strch1(text,nt,ntx,'c')
call strch1(text,nt,ntx,'d')
j = j - 400
goto 2
32 if (j.ge.499) goto 33
call strch1(text,nt,ntx,'x')
call strch1(text,nt,ntx,'d')
j = j - 490
goto 2
33 if (j.ge.500) goto 34
call strch1(text,nt,ntx,'i')
call strch1(text,nt,ntx,'d')
goto 99
34 if (j.ge.900) goto 35
call strch1(text,nt,ntx,'d')
j = j - 500
goto 2
35 if (j.ge.990) goto 36
call strch1(text,nt,ntx,'c')
call strch1(text,nt,ntx,'m')
j = j - 900
goto 2
36 if (j.ge.999) goto 37
call strch1(text,nt,ntx,'x')
call strch1(text,nt,ntx,'m')
j = j - 990
goto 2
37 call strch1(text,nt,ntx,'i')
call strch1(text,nt,ntx,'m')
goto 99
91 write (*,911) i
911 format(/,'[ROMNO] error: input zero or negative,',i10)
goto 99
99 continue
return
end
c to read a string of integers from keyboard. Characters not in the
c range 0..9 or - are treated as blanks.
subroutine readin(list,nlx,nl)
dimension text(255), list(nlx)
character*1 text, t
ntx = 255
read (*,1) (text(i),i=1,255)
1 format(255a1)
call strlen(text,ntx,nt)
nl = 0
nt0 = 0
2 if (nt0.eq.nt) goto 9
t = text(nt0+1)
m = ichar(t)
if ((m.eq.45).or.(m.eq.32)) goto 3
if ((m.ge.48).and.(m.le.57)) goto 3
nt0 = nt0 + 1
goto 2
3 call strtf1(text,ntx,nt0,nt1,mode,x,m)
if (mode.eq.1) goto 4
nl = nl + 1
list(nl) = m
if (nl.eq.nlx) goto 9
4 nt0 = nt1
goto 2
9 continue
return
end
c This program sorts the numbers in the array poj(1)-poj(nn)
c AND PUTS THE ORDER IN THE VECTOR jlist: IF THE LARGEST
C NUMBER IS POJ(10) THEN jlist(1) WILL BE 10.
C THE LENGTH OF THE VECTOR POJ IS LIMITED BY THE "DIMENSION"
C STATEMENT. THE DIMENSIONS GIVEN TO "jlist" AND "POJ" MUST BE AT
C LEAST NN AND THE DIMENSION GIVEN TO THE AUXILIARY VECTOR "KIST"
C MUST BE AT LEAST THE GREATER OF NN-M AND M/2 WHERE M IS
C THE LARGEST POWER OF 2 NOT GREATER THAN NN . ( (NN+1)/2 WILL
C ALWAYS DO.)
C THE PROGRAM PROCEEDS BY PUTTING SEGMENTS IN ORDER; AT STAGE
C K THE ORDERED SEGMENTS WILL BE OF LENGTH 2**K (STARTING FROM
C K = 0 ). IT TAKES ADJACENT SEGMENTS, COPIES THE LOWER ONE
C INTO KIST (LABEL 30), AND THEN PUTS THEM BOTH INTO jlist ,
C STARTING FROM THE BOTTOM.
C NOTE THAT THERE IS NO CONTROL AND IF NN IS TOO LARGE
C FOR THE VECTORS THERE WILL BE ERRORS WHICH MAY OR MAY NOT BE
C NOTIFIED.
SUBROUTINE SORTi(jlist,POJ,NN,KIST,IDIM1,IDIM2)
DIMENSION jlist(IDIM1),POJ(IDIM1),KIST(IDIM2)
integer poj
IF (NN.LE.0) GOTO 90
IF (NN.GT.IDIM1) GOTO 92
M = 1
DO 2 I = 1,34
M = 2**I
IF (M.GT.NN) GOTO 3
2 CONTINUE
GOTO 90
3 M = M/2
IF (((NN-M).GT.IDIM2).OR.((M/2).GT.IDIM2)) GOTO 94
NNBY2 = NN/2
IF (NNBY2.EQ.0) GOTO 17
DO 15 L = 1,NNBY2
M = 2*L
M0 = M
IF (POJ(M).GT.POJ(M-1)) M = M-1
M1 = 4*L - M - 1
jlist(M0-1) = M1
jlist(M0) = M
15 CONTINUE
17 IF (NN.NE.(2*NNBY2)) jlist(NN) = NN
KSTEP = 2
C KSTEP IS THE LENGTH OF THE SEGMENTS ORDERED SO FAR
20 IF (KSTEP.GE.NN) GOTO 100
KSTOP = 0
C KSTOP IS JUST ABOVE THE TOP OF THE UPPER SEGMENT TO BE WORKED ON
25 KSTART = KSTOP + KSTEP
C KSTART IS THE BOTTOM OF THE UPPER SEGMENT
C JSTART THE LENGTH OF THE LOWER SEGMENT
C KST1 THE BOTTOM OF THE LOWER SEGMENT
IF (KSTART.GE.NN) GOTO 60
JSTART = NN - KSTART
IF (JSTART.GT.KSTEP) JSTART = KSTEP
KST1 = JSTART + KSTART
C NOW COPY THE LOWER SEGMENT INTO KIST
DO 30 L = 1,JSTART
30 KIST(L) = jlist(L+KSTART)
C NOW ORDER THE TWO SEGMENTS TOGETHER; KNOW IS THE BOTTOM OF THE
C SO-FAR-UNTOUCHED PART OF THE UPPER SEGMENT; JNOW THE CORESPONDING
C ELEMENT OF THE LOWER SEGMENT (IN KIST ); NR THE LOCATION IN
C jlist WHERE THE INFERIOR OF THESE IS TO BE PUT
NR = KST1
KNOW = KSTART
JNOW = JSTART
35 LK = jlist(KNOW)
LJ = KIST(JNOW)
POJL = POJ(LK)
POJJ = POJ(LJ)
IF (POJL - POJJ) 40,50,50
40 jlist(NR) = LK
KNOW = KNOW - 1
NR = NR - 1
IF (KNOW.GT.KSTOP) GOTO 35
C IF THE UPPER SEGMENT IS EXHAUSTED, COPY THE LOWER SEGMENT BACK IN
DO 45 L = 1,JNOW
45 jlist(L+KSTOP) = KIST(L)
GOTO 55
50 jlist(NR) = LJ
JNOW = JNOW - 1
NR = NR - 1
IF (JNOW.GT.0) GOTO 35
55 KSTOP = KST1
GOTO 25
60 KSTEP = KSTEP + KSTEP
GOTO 20
90 write (6,91)NN
91 FORMAT(1X,'ERROR IN SORT - NN =',I6)
GOTO 98
92 write (6,93)NN,IDIM1
93 FORMAT(1X,'ERROR IN SORT - NN =',I6,' IDIM1 =',I6)
GOTO 98
94 write (6,95)NN,IDIM2
95 FORMAT(1X,'ERROR IN SORT - NN =',I6,' IDIM2 = ',I6)
GOTO 98
98 write (6,99)
99 FORMAT(1X,'SORT ABORTED')
NN = -1
100 CONTINUE
RETURN
END
|