      program icosaminx
c
c     Mike Collins -- cinclodes@yahoo.com
c
c     This FORTRAN code simulates on the screen a generalization of 
c     Rubik's Cube to the icosahedron. I call this puzzle the Icosaminx. 
c     It apparently does not exist outside the computer. It has about 
c     8.33x10^56 configurations. 
c
c     The output file icosaminx.ps contains the display. After compiling 
c     the code, run it and enter zeros when prompted. The code will 
c     quit after you enter the second zero. This initializes the 
c     display file, which you can now open with your Postscript 
c     displayer. Keep the display open and run the code again. This 
c     time, you can enter a random integer to mix up the code or zero 
c     again if you want to start with the puzzle solved. Then you can 
c     make a series of moves by entering one integer per line. Enter
c     negatives to twist in the opposite directions. You can reopen 
c     the display after each move.  
c
      real vx(20),vy(20),vz(20),fx(12,5),fy(12,5),fz(12,5),
     >   nx(20),ny(20),nz(20),tilex(140,4),tiley(140,4),
     >   gx(20,3),gy(20,3),gz(20,3)
      integer iclr(140),jclr(20)
      character*15 cname(20)
      pi=4.0*atan(1.0)
      theta=10.0
    1 cname(1)=' orange'
      cname(2)=' purple'
      cname(3)=' blue'
      cname(4)=' green'
      cname(5)=' red' 
      cname(6)=' yellow'
      cname(7)=' pine'
      cname(8)=' brown'
      cname(9)=' aqua'
      cname(10)=' white'
      cname(11)=' gray'
      cname(12)=' burgandy'
      cname(13)=' lavender'
      cname(14)=' tan'
      cname(15)=' blue-green'
      cname(16)=' lime'
      cname(17)=' cream' 
      cname(18)=' pink'
      cname(19)=' slate'
      cname(20)=' turquoise'
c
      jclr(1)=1
      jclr(2)=2
      jclr(3)=3
      jclr(4)=4
      jclr(5)=5
      jclr(6)=6
      jclr(7)=7
      jclr(8)=8
      jclr(9)=9
      jclr(10)=10
      jclr(11)=11
      jclr(12)=12
      jclr(13)=13
      jclr(14)=14
      jclr(15)=15
      jclr(16)=16
      jclr(17)=17
      jclr(18)=18
      jclr(19)=19
      jclr(20)=20
c
      n=0
      do 3 i=1,20
      do 2 j=1,7
      n=n+1
      iclr(n)=jclr(i)
    2 continue
    3 continue
c
      call vrtx(vx,vy,vz,pi)
      call faces(vx,vy,vz,fx,fy,fz,nx,ny,nz)
      call gclr(fx,fy,fz,gx,gy,gz,nx,ny,nz)
      call fclr(gx,gy,gz,nx,ny,nz,theta,pi,tilex,tiley,1.0)
      call fclr(gx,gy,gz,nx,ny,nz,theta,pi,tilex,tiley,-1.0)
c
c     Mix up the cube.
c
      write(*,*)' '
      write(*,*)'   Randomly mix up the puzzle. Enter an integer seed'
      write(*,*)'   for the random number generator. Enter 0 if you'
      write(*,*)'   dont want to mix up the cube.'
      write(*,*)' '
      read(*,*)iseed
c
      if(iseed.ne.0)then
      do 4 i=1,1000
      xran=1.0+19.999*ran0(iseed)
      n=ifix(xran)
      call opera(iclr,n,1,jclr,cname)
    4 continue
      end if
c
      open(unit=1,status='unknown',file='icosaminx.ps')
      call start(jclr,cname)
      call tclr(tilex,tiley,iclr)
      call frame(tilex,tiley)
      write(1,*)' showpage'
      close(1)
c
c     Solve the cube.
c
      write(*,*)' '
      write(*,*)' Enter an operator number one at a time.'
      write(*,*)' Definition of the operators is given in the display.'
      write(*,*)' Use negative numbers for inverse operators.'
      write(*,*)' Enter 99 to reset the cube. Enter 0 to quit.'
      write(*,*)' '
c
    5 read(*,*)i
      if(i.eq.0)go to 6
c
      k=1
      if(i.lt.0)then
      i=-i
      if(i.ne.21)k=2
      if(i.eq.21)k=4
      end if
c
      call opera(iclr,i,k,jclr,cname)
      open(unit=1,status='unknown',file='icosaminx.ps')
      call start(jclr,cname)
      call tclr(tilex,tiley,iclr)
      call frame(tilex,tiley)
      write(1,*)' showpage'
      close(1)
      if(i.eq.99)go to 1
      go to 5
c
    6 stop
      end
c
      subroutine opera(iclr,i,k,jclr,cname)
      integer iclr(140),jclr(20)
      character*15 cname(20),ctemp
c
      do 1 j=1,k
c
      if(i.eq.1)then
      call cycle(iclr,29,31,33)
      call cycle(iclr,30,32,34)
      call cycle(iclr,23,83,4)
      call cycle(iclr,3,22,82)
      call cycle(iclr,10,75,87)
      call cycle(iclr,17,101,36)
      call cycle(iclr,24,78,5)
      end if
c
      if(i.eq.2)then
      call cycle(iclr,22,24,26)
      call cycle(iclr,23,25,27)
      call cycle(iclr,32,16,76)
      call cycle(iclr,31,15,75)
      call cycle(iclr,3,54,101)
      call cycle(iclr,17,71,33)
      call cycle(iclr,10,94,78)
      end if
c
      if(i.eq.3)then
      call cycle(iclr,36,38,40)
      call cycle(iclr,37,39,41)
      call cycle(iclr,6,86,60)
      call cycle(iclr,59,5,85)
      call cycle(iclr,43,29,136)
      call cycle(iclr,12,82,131)
      call cycle(iclr,1,87,61)
      end if
c
      if(i.eq.4)then
      call cycle(iclr,2,4,6)
      call cycle(iclr,1,3,5)
      call cycle(iclr,11,30,41)
      call cycle(iclr,12,31,36)
      call cycle(iclr,43,24,87)
      call cycle(iclr,59,17,82)
      call cycle(iclr,40,10,29)
      end if
c
      if(i.eq.5)then
      call cycle(iclr,10,12,8)
      call cycle(iclr,11,13,9)
      call cycle(iclr,3,43,19)
      call cycle(iclr,31,59,50)
      call cycle(iclr,24,40,66)
      call cycle(iclr,17,1,47)
      call cycle(iclr,2,48,18)
      end if
c
      if(i.eq.6)then
      call cycle(iclr,15,17,19)
      call cycle(iclr,16,18,20)
      call cycle(iclr,25,9,55)
      call cycle(iclr,24,8,54)
      call cycle(iclr,31,47,94)
      call cycle(iclr,3,66,71)
      call cycle(iclr,10,50,26)
      end if
c
      if(i.eq.7)then
      call cycle(iclr,50,52,54)
      call cycle(iclr,51,53,55)
      call cycle(iclr,20,65,95)
      call cycle(iclr,19,64,94)
      call cycle(iclr,8,122,71)
      call cycle(iclr,47,117,26)
      call cycle(iclr,66,96,15)
      end if
c
      if(i.eq.8)then
      call cycle(iclr,57,59,61)
      call cycle(iclr,58,60,62)
      call cycle(iclr,44,39,130)
      call cycle(iclr,129,43,38)
      call cycle(iclr,124,12,85)
      call cycle(iclr,68,1,136)
      call cycle(iclr,45,40,131)
      end if
c
      if(i.eq.9)then
      call cycle(iclr,43,45,47)
      call cycle(iclr,44,46,48)
      call cycle(iclr,13,58,67)
      call cycle(iclr,12,57,66)
      call cycle(iclr,1,129,50)
      call cycle(iclr,40,124,19)
      call cycle(iclr,59,68,8)
      end if
c
      if(i.eq.10)then
      call cycle(iclr,64,66,68)
      call cycle(iclr,65,67,69)
      call cycle(iclr,51,46,123)
      call cycle(iclr,122,50,45)
      call cycle(iclr,117,19,57)
      call cycle(iclr,96,8,129)
      call cycle(iclr,52,47,124)
      end if
c
      if(i.eq.11)then
      call cycle(iclr,78,80,82)
      call cycle(iclr,79,81,83)
      call cycle(iclr,100,88,34)
      call cycle(iclr,33,99,87)
      call cycle(iclr,22,108,36)
      call cycle(iclr,75,138,5)
      call cycle(iclr,101,89,29)
      end if
c
      if(i.eq.12)then
      call cycle(iclr,99,101,103)
      call cycle(iclr,100,102,104)
      call cycle(iclr,109,79,74)
      call cycle(iclr,108,78,73)
      call cycle(iclr,138,33,92)
      call cycle(iclr,89,22,115)
      call cycle(iclr,80,75,110)
      end if
c
      if(i.eq.13)then
      call cycle(iclr,71,73,75)
      call cycle(iclr,72,74,76)
      call cycle(iclr,93,102,27)
      call cycle(iclr,26,92,101)
      call cycle(iclr,15,115,78)
      call cycle(iclr,54,110,33)
      call cycle(iclr,94,103,22)
      end if
c
      if(i.eq.14)then
      call cycle(iclr,85,87,89)
      call cycle(iclr,86,88,90)
      call cycle(iclr,37,81,137)
      call cycle(iclr,136,36,80)
      call cycle(iclr,131,5,99)
      call cycle(iclr,61,29,108)
      call cycle(iclr,38,82,138)
      end if
c
      if(i.eq.15)then
      call cycle(iclr,134,136,138)
      call cycle(iclr,135,137,139)
      call cycle(iclr,90,107,132)
      call cycle(iclr,106,131,89)
      call cycle(iclr,113,61,80)
      call cycle(iclr,120,38,99)
      call cycle(iclr,127,85,108)
      end if
c
      if(i.eq.16)then
      call cycle(iclr,106,108,110)
      call cycle(iclr,107,109,111)
      call cycle(iclr,139,104,114)
      call cycle(iclr,113,138,103)
      call cycle(iclr,120,89,73)
      call cycle(iclr,127,80,92)
      call cycle(iclr,134,99,115)
      end if
c
      if(i.eq.17)then
      call cycle(iclr,113,115,117)
      call cycle(iclr,114,116,118)
      call cycle(iclr,111,97,121)
      call cycle(iclr,120,110,96)
      call cycle(iclr,127,103,52)
      call cycle(iclr,134,73,64)
      call cycle(iclr,106,92,122)
      end if
c
      if(i.eq.18)then
      call cycle(iclr,92,94,96)
      call cycle(iclr,93,95,97)
      call cycle(iclr,116,72,53)
      call cycle(iclr,115,71,52)
      call cycle(iclr,110,26,64)
      call cycle(iclr,103,15,122)
      call cycle(iclr,73,54,117)
      end if
c
      if(i.eq.19)then
      call cycle(iclr,127,129,131)
      call cycle(iclr,128,130,132)
      call cycle(iclr,135,125,62)
      call cycle(iclr,134,124,61)
      call cycle(iclr,106,68,38)
      call cycle(iclr,113,45,85)
      call cycle(iclr,120,57,136)
      end if
c
      if(i.eq.20)then
      call cycle(iclr,120,122,124)
      call cycle(iclr,121,123,125)
      call cycle(iclr,128,118,69)
      call cycle(iclr,127,117,68)
      call cycle(iclr,134,96,45)
      call cycle(iclr,106,52,57)
      call cycle(iclr,113,64,129)
      end if
c
      if(i.eq.21)then
      call cycle5(iclr,35,28,77,105,84)
      call cycle5(iclr,33,22,75,101,78)
      call cycle5(iclr,34,23,76,102,79)
      call cycle5(iclr,29,24,71,103,80)
      call cycle5(iclr,30,25,72,104,81)
      call cycle5(iclr,31,26,73,99,82)
      call cycle5(iclr,32,27,74,100,83)
      call cycle5(iclr,7,21,98,112,91)
      call cycle5(iclr,1,19,96,106,85)
      call cycle5(iclr,2,20,97,107,86)
      call cycle5(iclr,3,15,92,108,87)
      call cycle5(iclr,4,16,93,109,88)
      call cycle5(iclr,5,17,94,110,89)
      call cycle5(iclr,6,18,95,111,90)
      call cycle5(iclr,42,14,56,119,140)
      call cycle5(iclr,36,10,54,115,138)
      call cycle5(iclr,37,11,55,116,139)
      call cycle5(iclr,38,12,50,117,134)
      call cycle5(iclr,39,13,51,118,135)
      call cycle5(iclr,40,8,52,113,136)
      call cycle5(iclr,41,9,53,114,137)
      call cycle5(iclr,63,49,70,126,133)
      call cycle5(iclr,57,45,68,124,129)
      call cycle5(iclr,58,46,69,125,130)
      call cycle5(iclr,59,47,64,120,131)
      call cycle5(iclr,60,48,65,121,132)
      call cycle5(iclr,61,43,66,122,127)
      call cycle5(iclr,62,44,67,123,128)
      ctemp=cname(2)
      cname(2)=cname(1)
      cname(1)=cname(11)
      cname(11)=cname(12)
      cname(12)=cname(13)
      cname(13)=ctemp
      ctemp=cname(4)
      cname(4)=cname(14)
      cname(14)=cname(16)
      cname(16)=cname(18)
      cname(18)=cname(6)
      cname(6)=ctemp
      ctemp=cname(7)
      cname(7)=cname(5)
      cname(5)=cname(3)
      cname(3)=cname(15)
      cname(15)=cname(17)
      cname(17)=ctemp
      ctemp=cname(10)
      cname(10)=cname(9)
      cname(9)=cname(8)
      cname(8)=cname(19)
      cname(19)=cname(20)
      cname(20)=ctemp
      end if
c
      if(i.eq.22)then
      call cycle(iclr,29,31,33)
      call cycle(iclr,30,32,34)
      call cycle(iclr,120,124,122)
      call cycle(iclr,121,125,123)
      call cycle(iclr,7,28,84)
      call cycle(iclr,1,26,80)
      call cycle(iclr,2,27,81)
      call cycle(iclr,3,22,82)
      call cycle(iclr,4,23,83)
      call cycle(iclr,5,24,78)
      call cycle(iclr,6,25,79)
      call cycle(iclr,14,77,91)
      call cycle(iclr,10,75,87)
      call cycle(iclr,11,76,88)
      call cycle(iclr,12,71,89)
      call cycle(iclr,13,72,90)
      call cycle(iclr,8,73,85)
      call cycle(iclr,9,74,86)
      call cycle(iclr,21,105,42)
      call cycle(iclr,17,101,36)
      call cycle(iclr,18,102,37)
      call cycle(iclr,19,103,38)
      call cycle(iclr,20,104,39)
      call cycle(iclr,15,99,40)
      call cycle(iclr,16,100,41)
      call cycle(iclr,70,119,133)
      call cycle(iclr,64,113,129)
      call cycle(iclr,65,114,130)
      call cycle(iclr,66,115,131)
      call cycle(iclr,67,116,132)
      call cycle(iclr,68,117,127)
      call cycle(iclr,69,118,128)
      call cycle(iclr,49,98,140)
      call cycle(iclr,45,96,134)
      call cycle(iclr,46,97,135)
      call cycle(iclr,47,92,136)
      call cycle(iclr,48,93,137)
      call cycle(iclr,43,94,138)
      call cycle(iclr,44,95,139)
      call cycle(iclr,56,112,63)
      call cycle(iclr,52,106,57)
      call cycle(iclr,53,107,58)
      call cycle(iclr,54,108,59)
      call cycle(iclr,55,109,60)
      call cycle(iclr,50,110,61)
      call cycle(iclr,51,111,62)
      ctemp=cname(2)
      cname(2)=cname(4)
      cname(4)=cname(11)
      cname(11)=ctemp
      ctemp=cname(6)
      cname(6)=cname(3)
      cname(3)=cname(12)
      cname(12)=ctemp
      ctemp=cname(5)
      cname(5)=cname(14)
      cname(14)=cname(13)
      cname(13)=ctemp
      ctemp=cname(7)
      cname(7)=cname(8)
      cname(8)=cname(16)
      cname(16)=ctemp
      ctemp=cname(9)
      cname(9)=cname(15)
      cname(15)=cname(18)
      cname(18)=ctemp
      ctemp=cname(17)
      cname(17)=cname(10)
      cname(10)=cname(19)
      cname(19)=ctemp
      end if
c
    1 continue
c
      return
      end
c
      subroutine cycle(iclr,i1,i2,i3)
      integer iclr(140)
c
      temp=iclr(i3)
      iclr(i3)=iclr(i2)
      iclr(i2)=iclr(i1)
      iclr(i1)=temp
c
      return
      end
c
      subroutine cycle5(iclr,i1,i2,i3,i4,i5)
      integer iclr(140)
c
      temp=iclr(i5)
      iclr(i5)=iclr(i4)
      iclr(i4)=iclr(i3)
      iclr(i3)=iclr(i2)
      iclr(i2)=iclr(i1)
      iclr(i1)=temp
c
      return
      end
c
      subroutine gclr(fx,fy,fz,gx,gy,gz,nx,ny,nz)
      real fx(12,5),fy(12,5),fz(12,5),gx(20,3),gy(20,3),
     >   gz(20,3),ux(12),uy(12),uz(12),nx(20),ny(20),nz(20)
c
c     Icosahedron vertices from the dodecahedron vertices.
c
      do 2 i=1,12
      ux(i)=0.0
      uy(i)=0.0
      uz(i)=0.0
      do 1 j=1,5
      ux(i)=ux(i)+0.2*fx(i,j)
      uy(i)=uy(i)+0.2*fy(i,j)
      uz(i)=uz(i)+0.2*fz(i,j)
    1 continue
    2 continue
c
c     The faces.
c
      gx(1,1)=ux(7)
      gx(1,2)=ux(2)
      gx(1,3)=ux(11)
      gy(1,1)=uy(7)
      gy(1,2)=uy(2)
      gy(1,3)=uy(11)
      gz(1,1)=uz(7)
      gz(1,2)=uz(2)
      gz(1,3)=uz(11)
c
      gx(2,1)=ux(2)
      gx(2,2)=ux(6)
      gx(2,3)=ux(11)
      gy(2,1)=uy(2)
      gy(2,2)=uy(6)
      gy(2,3)=uy(11)
      gz(2,1)=uz(2)
      gz(2,2)=uz(6)
      gz(2,3)=uz(11)
c
      gx(3,1)=ux(6)
      gx(3,2)=ux(10)
      gx(3,3)=ux(11)
      gy(3,1)=uy(6)
      gy(3,2)=uy(10)
      gy(3,3)=uy(11)
      gz(3,1)=uz(6)
      gz(3,2)=uz(10)
      gz(3,3)=uz(11)
c
      gx(4,1)=ux(10)
      gx(4,2)=ux(12)
      gx(4,3)=ux(11)
      gy(4,1)=uy(10)
      gy(4,2)=uy(12)
      gy(4,3)=uy(11)
      gz(4,1)=uz(10)
      gz(4,2)=uz(12)
      gz(4,3)=uz(11)
c
      gx(5,1)=ux(12)
      gx(5,2)=ux(7)
      gx(5,3)=ux(11)
      gy(5,1)=uy(12)
      gy(5,2)=uy(7)
      gy(5,3)=uy(11)
      gz(5,1)=uz(12)
      gz(5,2)=uz(7)
      gz(5,3)=uz(11)
c
      gx(6,1)=ux(2)
      gx(6,2)=ux(7)
      gx(6,3)=ux(3)
      gy(6,1)=uy(2)
      gy(6,2)=uy(7)
      gy(6,3)=uy(3)
      gz(6,1)=uz(2)
      gz(6,2)=uz(7)
      gz(6,3)=uz(3)
c
      gx(7,1)=ux(6)
      gx(7,2)=ux(2)
      gx(7,3)=ux(1)
      gy(7,1)=uy(6)
      gy(7,2)=uy(2)
      gy(7,3)=uy(1)
      gz(7,1)=uz(6)
      gz(7,2)=uz(2)
      gz(7,3)=uz(1)
c
      gx(8,1)=ux(10)
      gx(8,2)=ux(6)
      gx(8,3)=ux(5)
      gy(8,1)=uy(10)
      gy(8,2)=uy(6)
      gy(8,3)=uy(5)
      gz(8,1)=uz(10)
      gz(8,2)=uz(6)
      gz(8,3)=uz(5)
c
      gx(9,1)=ux(12)
      gx(9,2)=ux(10)
      gx(9,3)=ux(9)
      gy(9,1)=uy(12)
      gy(9,2)=uy(10)
      gy(9,3)=uy(9)
      gz(9,1)=uz(12)
      gz(9,2)=uz(10)
      gz(9,3)=uz(9)
c
      gx(10,1)=ux(7)
      gx(10,2)=ux(12)
      gx(10,3)=ux(8)
      gy(10,1)=uy(7)
      gy(10,2)=uy(12)
      gy(10,3)=uy(8)
      gz(10,1)=uz(7)
      gz(10,2)=uz(12)
      gz(10,3)=uz(8)
c
      gx(11,1)=ux(8)
      gx(11,2)=ux(3)
      gx(11,3)=ux(7)
      gy(11,1)=uy(8)
      gy(11,2)=uy(3)
      gy(11,3)=uy(7)
      gz(11,1)=uz(8)
      gz(11,2)=uz(3)
      gz(11,3)=uz(7)
c
      gx(12,1)=ux(3)
      gx(12,2)=ux(1)
      gx(12,3)=ux(2)
      gy(12,1)=uy(3)
      gy(12,2)=uy(1)
      gy(12,3)=uy(2)
      gz(12,1)=uz(3)
      gz(12,2)=uz(1)
      gz(12,3)=uz(2)
c
      gx(13,1)=ux(1)
      gx(13,2)=ux(5)
      gx(13,3)=ux(6)
      gy(13,1)=uy(1)
      gy(13,2)=uy(5)
      gy(13,3)=uy(6)
      gz(13,1)=uz(1)
      gz(13,2)=uz(5)
      gz(13,3)=uz(6)
c
      gx(14,1)=ux(5)
      gx(14,2)=ux(9)
      gx(14,3)=ux(10)
      gy(14,1)=uy(5)
      gy(14,2)=uy(9)
      gy(14,3)=uy(10)
      gz(14,1)=uz(5)
      gz(14,2)=uz(9)
      gz(14,3)=uz(10)
c
      gx(15,1)=ux(9)
      gx(15,2)=ux(8)
      gx(15,3)=ux(12)
      gy(15,1)=uy(9)
      gy(15,2)=uy(8)
      gy(15,3)=uy(12)
      gz(15,1)=uz(9)
      gz(15,2)=uz(8)
      gz(15,3)=uz(12)
c
      gx(16,1)=ux(9)
      gx(16,2)=ux(4)
      gx(16,3)=ux(8)
      gy(16,1)=uy(9)
      gy(16,2)=uy(4)
      gy(16,3)=uy(8)
      gz(16,1)=uz(9)
      gz(16,2)=uz(4)
      gz(16,3)=uz(8)
c
      gx(17,1)=ux(5)
      gx(17,2)=ux(4)
      gx(17,3)=ux(9)
      gy(17,1)=uy(5)
      gy(17,2)=uy(4)
      gy(17,3)=uy(9)
      gz(17,1)=uz(5)
      gz(17,2)=uz(4)
      gz(17,3)=uz(9)
c
      gx(18,1)=ux(1)
      gx(18,2)=ux(4)
      gx(18,3)=ux(5)
      gy(18,1)=uy(1)
      gy(18,2)=uy(4)
      gy(18,3)=uy(5)
      gz(18,1)=uz(1)
      gz(18,2)=uz(4)
      gz(18,3)=uz(5)
c
      gx(19,1)=ux(3)
      gx(19,2)=ux(4)
      gx(19,3)=ux(1)
      gy(19,1)=uy(3)
      gy(19,2)=uy(4)
      gy(19,3)=uy(1)
      gz(19,1)=uz(3)
      gz(19,2)=uz(4)
      gz(19,3)=uz(1)
c
      gx(20,1)=ux(8)
      gx(20,2)=ux(4)
      gx(20,3)=ux(3)
      gy(20,1)=uy(8)
      gy(20,2)=uy(4)
      gy(20,3)=uy(3)
      gz(20,1)=uz(8)
      gz(20,2)=uz(4)
      gz(20,3)=uz(3)
c
      do 3 i=1,20
      x1=gx(i,2)-gx(i,1)
      y1=gy(i,2)-gy(i,1)
      z1=gz(i,2)-gz(i,1)
      x2=gx(i,3)-gx(i,2)
      y2=gy(i,3)-gy(i,2)
      z2=gz(i,3)-gz(i,2)
      nx(i)=y1*z2-y2*z1
      ny(i)=-x1*z2+x2*z1
      nz(i)=x1*y2-x2*y1
    3 continue
c
      return
      end
c
      subroutine fclr(fx,fy,fz,nx,ny,nz,theta,pi,tilex,tiley,sign)
      real fx(20,3),fy(20,3),fz(20,3),nx(20),ny(20),nz(20),
     >   tilex(140,4),tiley(140,4),x(20),y(20)
c
      wx=cos(theta*pi/180.0)
      wy=0.0
      wz=sin(theta*pi/180.0)
c
      if(sign.gt.0.0)itile=0
      if(sign.lt.0.0)itile=70
      do 1 i=1,20
      dot=wx*nx(i)+wy*ny(i)+wz*nz(i)
      if(sign*dot.gt.0.0)then
      x(1)=fy(i,1)
      y(1)=fz(i,1)*cos(theta*pi/180.0)-fx(i,1)*sin(theta*pi/180.0)
      x(2)=fy(i,2)
      y(2)=fz(i,2)*cos(theta*pi/180.0)-fx(i,2)*sin(theta*pi/180.0)
      x(3)=fy(i,3)
      y(3)=fz(i,3)*cos(theta*pi/180.0)-fx(i,3)*sin(theta*pi/180.0)
c
      alp=0.2
      bet=1.0-alp
      x(4)=bet*x(1)+alp*x(2)
      y(4)=bet*y(1)+alp*y(2)
      x(5)=alp*x(1)+bet*x(2)
      y(5)=alp*y(1)+bet*y(2)
      x(6)=bet*x(2)+alp*x(3)
      y(6)=bet*y(2)+alp*y(3)
      x(7)=alp*x(2)+bet*x(3)
      y(7)=alp*y(2)+bet*y(3)
      x(8)=bet*x(3)+alp*x(1)
      y(8)=bet*y(3)+alp*y(1)
      x(9)=alp*x(3)+bet*x(1)
      y(9)=alp*y(3)+bet*y(1)
c
      xa=x(7)
      ya=y(7)
      xb=x(4)
      yb=y(4)
      xc=x(9)
      yc=y(9)
      xd=x(6)
      yd=y(6)
      call lines(xa,ya,xb,yb,xc,yc,xd,yd,xint,yint)
      x(10)=xint
      y(10)=yint
c
      xa=x(9)
      ya=y(9)
      xb=x(6)
      yb=y(6)
      xc=x(5)
      yc=y(5)
      xd=x(8)
      yd=y(8)
      call lines(xa,ya,xb,yb,xc,yc,xd,yd,xint,yint)
      x(11)=xint
      y(11)=yint
c
      xa=x(5)
      ya=y(5)
      xb=x(8)
      yb=y(8)
      xc=x(7)
      yc=y(7)
      xd=x(4)
      yd=y(4)
      call lines(xa,ya,xb,yb,xc,yc,xd,yd,xint,yint)
      x(12)=xint
      y(12)=yint
c
      call scale(x,y,sign)
c
      itile=itile+1
      tilex(itile,1)=x(2)
      tiley(itile,1)=y(2)
      tilex(itile,2)=x(6)
      tiley(itile,2)=y(6)
      tilex(itile,3)=x(11)
      tiley(itile,3)=y(11)
      tilex(itile,4)=x(5)
      tiley(itile,4)=y(5)
c
      itile=itile+1
      tilex(itile,1)=x(6)
      tiley(itile,1)=y(6)
      tilex(itile,2)=x(7)
      tiley(itile,2)=y(7)
      tilex(itile,3)=x(12)
      tiley(itile,3)=y(12)
      tilex(itile,4)=x(11)
      tiley(itile,4)=y(11)
c
      itile=itile+1
      tilex(itile,1)=x(3)
      tiley(itile,1)=y(3)
      tilex(itile,2)=x(8)
      tiley(itile,2)=y(8)
      tilex(itile,3)=x(12)
      tiley(itile,3)=y(12)
      tilex(itile,4)=x(7)
      tiley(itile,4)=y(7)
c
      itile=itile+1
      tilex(itile,1)=x(8)
      tiley(itile,1)=y(8)
      tilex(itile,2)=x(9)
      tiley(itile,2)=y(9)
      tilex(itile,3)=x(10)
      tiley(itile,3)=y(10)
      tilex(itile,4)=x(12)
      tiley(itile,4)=y(12)
c
      itile=itile+1
      tilex(itile,1)=x(9)
      tiley(itile,1)=y(9)
      tilex(itile,2)=x(1)
      tiley(itile,2)=y(1)
      tilex(itile,3)=x(4)
      tiley(itile,3)=y(4)
      tilex(itile,4)=x(10)
      tiley(itile,4)=y(10)
c
      itile=itile+1
      tilex(itile,1)=x(4)
      tiley(itile,1)=y(4)
      tilex(itile,2)=x(5)
      tiley(itile,2)=y(5)
      tilex(itile,3)=x(11)
      tiley(itile,3)=y(11)
      tilex(itile,4)=x(10)
      tiley(itile,4)=y(10)
c
      itile=itile+1
      tilex(itile,1)=x(11)
      tiley(itile,1)=y(11)
      tilex(itile,2)=x(12)
      tiley(itile,2)=y(12)
      tilex(itile,3)=0.5*(x(12)+x(10))
      tiley(itile,3)=0.5*(y(12)+y(10))
      tilex(itile,4)=x(10)
      tiley(itile,4)=y(10)
c
      end if
    1 continue
c
      return
      end
c
      subroutine lines(xa,ya,xb,yb,xc,yc,xd,yd,xint,yint)
      a=xd-xc
      b=xa-xb
      c=yd-yc
      d=ya-yb
      r1=xa-xc
      r2=ya-yc
      s=(d*r1-b*r2)/(a*d-b*c)
      xint=xc+s*(xd-xc)
      yint=yc+s*(yd-yc)
      return
      end
c
      subroutine scale(x,y,sign)
      real x(20),y(20)
c
      do 1 i=1,12
      if(sign.gt.0.0)then
      x(i)=410.0+112.0*x(i)
      y(i)=525.0+112.0*y(i)
      else
      x(i)=410.0+112.0*x(i)
      y(i)=220.0+112.0*y(i)
      end if
    1 continue
      return
      end
c
      subroutine tclr(tilex,tiley,iclr)
      real tilex(140,4),tiley(140,4)
      integer iclr(140)
c
c
      do 1 i=1,140
      write(1,*)tilex(i,1),tiley(i,1),' moveto'
      write(1,*)tilex(i,2),tiley(i,2),' lineto'
      write(1,*)tilex(i,3),tiley(i,3),' lineto'
      write(1,*)tilex(i,4),tiley(i,4),' lineto'
      write(1,*)' closepath'
      n=iclr(i)
      call color(n)
      write(1,*)' fill'
    1 continue
    2 continue
c
      return
      end
c
      subroutine frame(tilex,tiley)
      real tilex(140,4),tiley(140,4)
c
      do 1 i=1,140
      write(1,*)tilex(i,1),tiley(i,1),' moveto'
      write(1,*)tilex(i,2),tiley(i,2),' lineto'
      write(1,*)tilex(i,3),tiley(i,3),' lineto'
      write(1,*)tilex(i,4),tiley(i,4),' lineto'
      write(1,*)' closepath'
      call color(99)
      write(1,*)' stroke'
    1 continue
    2 continue
c
      return
      end
c
      subroutine color(n)
c
      if(n.eq.1)write(1,*)' 0.333 1 1 sethsbcolor'
      if(n.eq.2)write(1,*)' 0 1 1 sethsbcolor'
      if(n.eq.3)write(1,*)' 0.167 1 1 sethsbcolor'
      if(n.eq.4)write(1,*)' 0.75 1 1 sethsbcolor'
      if(n.eq.5)write(1,*)' 0.1 1 1 sethsbcolor'
      if(n.eq.6)write(1,*)' 0.667 1 1  sethsbcolor'
      if(n.eq.7)write(1,*)' 0.55 1 1  sethsbcolor'
      if(n.eq.8)write(1,*)' 0.333 1 0.5 sethsbcolor'
      if(n.eq.9)write(1,*)' 0.1 1 0.6 sethsbcolor'
      if(n.eq.10)write(1,*)' 1 1 1  setrgbcolor'
      if(n.eq.16)write(1,*)' 0.85 1.0 0.4 setrgbcolor'
      if(n.eq.12)write(1,*)' 0.7 0.7 0.7  setrgbcolor'
      if(n.eq.11)write(1,*)' 0.75 0.4 1 sethsbcolor'
      if(n.eq.20)write(1,*)' 0.45 0.4 1 sethsbcolor'
      if(n.eq.15)write(1,*)' 0 1 0.8 sethsbcolor'
      if(n.eq.13)write(1,*)' 1 0.7 0.5 setrgbcolor'
      if(n.eq.17)write(1,*)' 0.166 0.2 1 sethsbcolor'
      if(n.eq.18)write(1,*)' 0 1 1 setrgbcolor'
      if(n.eq.19)write(1,*)' 0.35 0.35 0.35 setrgbcolor'
      if(n.eq.14)write(1,*)' 1 0.15 0.7 setrgbcolor'
      if(n.eq.44)write(1,*)' 1 1 1 setrgbcolor'
      if(n.eq.99)write(1,*)' 0 0 0 setrgbcolor'
c
      return
      end
c
      subroutine color2(n)
c
      if(n.eq.1)write(1,*)' 0.333 1 1 sethsbcolor'
      if(n.eq.2)write(1,*)' 0 1 1 sethsbcolor'
      if(n.eq.3)write(1,*)' 0.167 1 1 sethsbcolor'
      if(n.eq.4)write(1,*)' 0.75 1 1 sethsbcolor'
      if(n.eq.5)write(1,*)' 0.1 1 1 sethsbcolor'
      if(n.eq.6)write(1,*)' 0.667 1 1  sethsbcolor'
      if(n.eq.7)write(1,*)' 0.55 1 1  sethsbcolor'
      if(n.eq.8)write(1,*)' 0.333 1 0.5 sethsbcolor'
      if(n.eq.9)write(1,*)' 0.1 1 0.6 sethsbcolor'
      if(n.eq.10)write(1,*)' 1 1 1  setrgbcolor'
      if(n.eq.11)write(1,*)' 0.1 1 0.6 sethsbcolor'
      if(n.eq.12)write(1,*)' 1 1 1  setrgbcolor'
      if(n.eq.13)write(1,*)' 0.333 1 0.5 sethsbcolor'
      if(n.eq.14)write(1,*)' 0.667 1 1  sethsbcolor'
      if(n.eq.15)write(1,*)' 0.55 1 1  sethsbcolor'
      if(n.eq.16)write(1,*)' 0 1 1 sethsbcolor'
      if(n.eq.17)write(1,*)' 0.333 1 1 sethsbcolor'
      if(n.eq.18)write(1,*)' 0.1 1 1 sethsbcolor'
      if(n.eq.19)write(1,*)' 0.75 1 1 sethsbcolor'
      if(n.eq.20)write(1,*)' 0.167 1 1 sethsbcolor'
      if(n.eq.44)write(1,*)' 1 1 1 setrgbcolor'
      if(n.eq.99)write(1,*)' 0 0 0 setrgbcolor'
c
      return
      end
c
      subroutine vrtx(vx,vy,vz,pi)
      real vx(20),vy(20),vz(20)
      cost=cos(72.0*pi/180.0)
      sint=sin(72.0*pi/180.0)
      cost2=cos(36.0*pi/180.0)
      sint2=sin(36.0*pi/180.0)
c
      vx(1)=-0.5/sin(36.0*pi/180.0)
      vy(1)=0.0
      vz(1)=0.0
      do 1 i=2,5
      vx(i)=cost*vx(i-1)-sint*vy(i-1)
      vy(i)=sint*vx(i-1)+cost*vy(i-1)
      vz(i)=0.0
    1 continue
c
      alpha=cos(108.0*pi/180.0)/cos(54.0*pi/180.0)
      dz1=sqrt(1.0-alpha**2)
      vx(6)=vx(1)+alpha
      vy(6)=0.0
      vz(6)=dz1
      do 2 i=7,10
      vx(i)=cost*vx(i-1)-sint*vy(i-1)
      vy(i)=sint*vx(i-1)+cost*vy(i-1)
      vz(i)=vz(6)
    2 continue
c
      dz2=dz1*cos(54.0*pi/180.0)/sin(72.0*pi/180.0)
      vx(11)=cost2*vx(6)-sint2*vy(6)
      vy(11)=sint2*vx(6)+cost2*vy(6)
      vz(11)=dz1+dz2
      do 3 i=12,15
      vx(i)=cost*vx(i-1)-sint*vy(i-1)
      vy(i)=sint*vx(i-1)+cost*vy(i-1)
      vz(i)=vz(11)
    3 continue
c
      vx(16)=cost2*vx(1)-sint2*vy(1)
      vy(16)=sint2*vx(1)+cost2*vy(1)
      vz(16)=2.0*dz1+dz2
      do 4 i=17,20
      vx(i)=cost*vx(i-1)-sint*vy(i-1)
      vy(i)=sint*vx(i-1)+cost*vy(i-1)
      vz(i)=vz(16)
    4 continue
c
      return
      end
c
      subroutine faces(vx,vy,vz,fx,fy,fz,nx,ny,nz)
      real vx(20),vy(20),vz(20),fx(12,5),fy(12,5),fz(12,5),
     >   nx(20),ny(20),nz(20)
c
      call setfc(vx,vy,vz,fx,fy,fz,1,1,1)
      call setfc(vx,vy,vz,fx,fy,fz,1,2,2)
      call setfc(vx,vy,vz,fx,fy,fz,1,3,3)
      call setfc(vx,vy,vz,fx,fy,fz,1,4,4)
      call setfc(vx,vy,vz,fx,fy,fz,1,5,5)
c
      call setfc(vx,vy,vz,fx,fy,fz,2,1,1)
      call setfc(vx,vy,vz,fx,fy,fz,2,2,6)
      call setfc(vx,vy,vz,fx,fy,fz,2,3,11)
      call setfc(vx,vy,vz,fx,fy,fz,2,4,7)
      call setfc(vx,vy,vz,fx,fy,fz,2,5,2)
c
      call setfc(vx,vy,vz,fx,fy,fz,3,1,2)
      call setfc(vx,vy,vz,fx,fy,fz,3,2,7)
      call setfc(vx,vy,vz,fx,fy,fz,3,3,12)
      call setfc(vx,vy,vz,fx,fy,fz,3,4,8)
      call setfc(vx,vy,vz,fx,fy,fz,3,5,3)
c
      call setfc(vx,vy,vz,fx,fy,fz,4,1,3)
      call setfc(vx,vy,vz,fx,fy,fz,4,2,8)
      call setfc(vx,vy,vz,fx,fy,fz,4,3,13)
      call setfc(vx,vy,vz,fx,fy,fz,4,4,9)
      call setfc(vx,vy,vz,fx,fy,fz,4,5,4)
c
      call setfc(vx,vy,vz,fx,fy,fz,5,1,4)
      call setfc(vx,vy,vz,fx,fy,fz,5,2,9)
      call setfc(vx,vy,vz,fx,fy,fz,5,3,14)
      call setfc(vx,vy,vz,fx,fy,fz,5,4,10)
      call setfc(vx,vy,vz,fx,fy,fz,5,5,5)
c
      call setfc(vx,vy,vz,fx,fy,fz,6,1,5)
      call setfc(vx,vy,vz,fx,fy,fz,6,2,10)
      call setfc(vx,vy,vz,fx,fy,fz,6,3,15)
      call setfc(vx,vy,vz,fx,fy,fz,6,4,6)
      call setfc(vx,vy,vz,fx,fy,fz,6,5,1)
c
      call setfc(vx,vy,vz,fx,fy,fz,7,1,16)
      call setfc(vx,vy,vz,fx,fy,fz,7,2,17)
      call setfc(vx,vy,vz,fx,fy,fz,7,3,12)
      call setfc(vx,vy,vz,fx,fy,fz,7,4,7)
      call setfc(vx,vy,vz,fx,fy,fz,7,5,11)
c
      call setfc(vx,vy,vz,fx,fy,fz,8,1,17)
      call setfc(vx,vy,vz,fx,fy,fz,8,2,18)
      call setfc(vx,vy,vz,fx,fy,fz,8,3,13)
      call setfc(vx,vy,vz,fx,fy,fz,8,4,8)
      call setfc(vx,vy,vz,fx,fy,fz,8,5,12)
c
      call setfc(vx,vy,vz,fx,fy,fz,9,1,18)
      call setfc(vx,vy,vz,fx,fy,fz,9,2,19)
      call setfc(vx,vy,vz,fx,fy,fz,9,3,14)
      call setfc(vx,vy,vz,fx,fy,fz,9,4,9)
      call setfc(vx,vy,vz,fx,fy,fz,9,5,13)
c
      call setfc(vx,vy,vz,fx,fy,fz,10,1,19)
      call setfc(vx,vy,vz,fx,fy,fz,10,2,20)
      call setfc(vx,vy,vz,fx,fy,fz,10,3,15)
      call setfc(vx,vy,vz,fx,fy,fz,10,4,10)
      call setfc(vx,vy,vz,fx,fy,fz,10,5,14)
c
      call setfc(vx,vy,vz,fx,fy,fz,11,1,20)
      call setfc(vx,vy,vz,fx,fy,fz,11,2,16)
      call setfc(vx,vy,vz,fx,fy,fz,11,3,11)
      call setfc(vx,vy,vz,fx,fy,fz,11,4,6)
      call setfc(vx,vy,vz,fx,fy,fz,11,5,15)
c
      call setfc(vx,vy,vz,fx,fy,fz,12,1,20)
      call setfc(vx,vy,vz,fx,fy,fz,12,2,19)
      call setfc(vx,vy,vz,fx,fy,fz,12,3,18)
      call setfc(vx,vy,vz,fx,fy,fz,12,4,17)
      call setfc(vx,vy,vz,fx,fy,fz,12,5,16)
c
      return
      end
c
      subroutine setfc(vx,vy,vz,fx,fy,fz,i,j,k)
      real vx(20),vy(20),vz(20),fx(12,5),fy(12,5),fz(12,5)
      fx(i,j)=vx(k)
      fy(i,j)=vy(k)
      fz(i,j)=vz(k)
      return
      end
c
      subroutine start(jclr,cname)
      integer jclr(20)
      character*15 cname(20)
c
      x=40
      z=800
      write(1,*)'%!'
      write(1,*)' /Times-Roman findfont'
      write(1,*)' 14 scalefont setfont'
      z=z-40
      write(1,*)x,z,' moveto'
      write(1,*)' (Rotate a front face:) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  1 =',cname(jclr(1)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  2 =',cname(jclr(2)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  3 =',cname(jclr(3)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  4 =',cname(jclr(4)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  5 =',cname(jclr(5)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  6 =',cname(jclr(6)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  7 =',cname(jclr(7)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  8 =',cname(jclr(8)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  9 =',cname(jclr(9)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  10 =',cname(jclr(10)),' ) show'
      z=z-40
      write(1,*)x,z,' moveto'
      write(1,*)' (Rotate a back face:) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  11 =',cname(jclr(11)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  12 =',cname(jclr(12)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  13 =',cname(jclr(13)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  14 =',cname(jclr(14)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  15 =',cname(jclr(15)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  16 =',cname(jclr(16)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  17 =',cname(jclr(17)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  18 =',cname(jclr(18)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  19 =',cname(jclr(19)),' ) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  20 =',cname(jclr(20)),' ) show'
c
      z=z-40
      write(1,*)x,z,' moveto'
      write(1,*)' (Rotate the entire icosahedron:) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  21 = about the vertical axis) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (  22 = about the 1-20 axis) show'
c
      z=z-40
      write(1,*)x,z,' moveto'
      write(1,*)' (All rotations are counterclockwise.) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (Rotations of back faces appear clockwise.) show'
c
      x=300
      z=500
      write(1,*)x,z,' moveto'
      write(1,*)' (Above: Direct view of the front faces.) show'
      z=z-20
      write(1,*)x,z,' moveto'
      write(1,*)' (Below: Mirror view of the back faces.) show'
c
      return
      end
c
      function ran0(iseed)
      if(iseed.eq.0)iseed=178544878
      m=2147483647
      n=16807
      k=iseed/127773
      iseed=n*(iseed-k*127773)-k*2836
      if(iseed.lt.0)iseed=iseed+m
      ran0=float(iseed)/float(m)
      return
      end
