C****************************************************************** C This program calculate the vertex positions of 2D quasiperiodic tilings C using the generalized dual method. It is just a simple program C from were you can start other applications like the 3D version (it is just C very easy to do,from the present code). C First, you give the star vectors e_{i}, C i.e. a set of unitary vectors that have the required symmetry. For example, C in the Penrose it is required 5 fold symmetry, thus, the star is a set C of 5 vectors that point to the vertex of a pentagon. Then, all intersections C of a family of lines perpendicular to each vector is calculated. Following C the formulas of G.G. Naumis, J.L. Aragón,"Analytic expressions for the C vertex coordinates of quasiperiodic lattices",Z. Kristallogr. 218, 397 (2003). C This program has been made by Gerardo G. Naumis and José Luis Aragon, C in 2002, at the Instituto de Fisica, National University of Mexico. C To draw the tiling, just make a program that joins all points separated C by a distance "1". See our mathematica version of the program in the C same webpage: www.fisica.unam.mx/naumis C If you have comments, contact: G. Naumis naumis@fisica.unam.mx C********************************************************************* C Variable declaration REAL EX(5),EY(5),UX(5),UY(5),GAMMA(5),rx,ry,arx0,ary0,aux,AL,PI REAL Aij,arx,ary COMMON H,NTP,U,V,PI,N0,NTG,N,NFR integer iseed,irand,i,ID(200),KD(200),KA double precision rand C++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ C Parameters PI=4.0*ATAN(1.0) AL=2.0*PI/5.0 C******************************************************************* C Following the notation of Steinhardt and Levine, on needs to define C phases that shifts the family of lines on each vector. C Different isomorphisms (different types of tilings) are obtained C if these phases are changed. If the sum of these phases is 0, C one obtain the Penrose isomorphism class, which is the most famous. C see: J Socolar, Steinhardt, Phys. Rev. B 34, 617 (1986). C However, many people forget to tell you that you must be careful C with singular tilings,i.e., if these shifts are so that C one has points where three of more lines cross, the tiles are no longer C rhombuses, one can have hexagons (duals of three ´lines that meet in a point) C octagons (four lines), decagon (5 lines). This last one ocuurs at the origin C if all phases are set to zero. C A FAQ is "why I don´t get the right Penrose?. First, you will always C get different portions of the Penrose tiling unless you know which are the C exact parameters of the tile that you want to compare. Usually,instead one C needs to check that both tilings are in the same local isomorphism class, C that means that each patch of a given radius, is observed in both tilings. C To check the Penrose, as first check, verify the proportions of each C vertex configuration and second neighbour. A table was published in: C G.G. Naumis, "Density of states and first spectral moments of a C quasiperiodic lattice", J. Phys: Condens. Matter 11 (1999) 7 143-7 154. C I an exotic symmetry is required, change the number of star vectors C************************************************************************ GAMMA(1)=0.0 GAMMA(2)=-0.50 GAMMA(3)=0.25 GAMMA(4)=0.25 GAMMA(5)=0.44 C write(*,*) AINT(0.1135) C write(*,*) AINT(2.5760) C write(*,*) AINT(-1.2563) C************************************************************************ C Definition of star (basis vectors), and generation of auxiliary C vectors that allow to get the intersection of each family of lines C************************************************************************ DO 10 i=1,5 ri=REAL(i)-1.0 EX(i)=cos(AL*ri) EY(i)=sin(AL*ri) UX(i)=-sin(AL*ri) UY(i)=cos(AL*ri) write(*,*) EX(i),EY(i) 10 CONTINUE C************************************************************************** C Just for the record, two files are generated as output, one with the C intersection points of each family ("pensorseint.dat"). This is known C as the skeleton. The other file, "penrsoe.dat" contains the vertex positions. C The Average lattice is also calculated and store in "penroseav.dat" C************************************************************************** OPEN(2,FILE='penroseav.dat',STATUS='UNKNOWN') OPEN(4,FILE='penrose.dat',STATUS='UNKNOWN') OPEN(6,FILE='penroseint.dat',STATUS='UNKNOWN') C************************************************************************** C A loop is carried over I and J which are couples of vectors that generate C all the intersections. For each copule I and J, we have a family of lines C indexed by m and n, which are the "numeration" of each line over the directions C I and J. C************************************************************************** DO 20 i=1,5 DO 20 j=i+1,5 if (i.NE.j) then Aij=((EX(i)*EY(j))-(EY(i)*EX(j))) DO 30 n=-10,10 DO 30 m=-10,10 rx=(REAL(n)*EX(i))+(REAL(m)*EX(j)) ry=(REAL(n)*EY(i))+(REAL(m)*EY(j)) arx0=(-(REAL(m)+GAMMA(j))*EY(i)+(REAL(n)+GAMMA(i))*EY(j))/Aij ary0=(((REAL(m)+GAMMA(j))*EX(i)-(REAL(n)+GAMMA(i))*EX(j)))/Aij co1=(arx0*EX(i))+(ary0*EY(i))-REAL(n)-GAMMA(i) co2=(arx0*EX(j))+(ary0*EY(j))-REAL(m)-GAMMA(j) C*************************************************************************** C Write the family and interseection coordinates of two lines, C which in fact correspond to the center of a rhombus. C*************************************************************************** C write(*,*) 'N',N,' M',M,' INT=',arx0,ary0 if ((ABS(co1).GT.0.0001).OR.(ABS(co2).GT.0.0001)) then write(*,*) 'error',I,J,n,m,co1,co2 endif C write(*,*) 'I',I,' J',J,' N',N,' M',M,co1 C**************************************************************************** C Dual transformation C Calculate ordinal coordinates with respect to each family of lines C of the generated intersection. Of course, coordinates I and J are known C already, since they were used to generate the point. Thus, such step C is not performed for lines of family I and J. C**************************************************************************** arxp=rx aryp=ry DO 40 L=1,5 if ((L.NE.i).AND.(L.NE.j)) then aux=(arx0*EX(L))+(ary0*EY(L))-GAMMA(L) if (aux.GE.0.0) then n1=AINT(aux) else n1=AINT(aux)-1 endif C****************************************************************** C Perform in a recursive way the dual transformation by making a linear C combination with the right coefficients of the star vectors. C****************************************************************** arx=rx+(REAL(n1)*EX(L)) ary=ry+(REAL(n1)*EY(L)) arxp=arxp+((aux-0.5)*EX(L)) aryp=aryp+((aux-0.5)*EY(L)) rx=arx ry=ary endif 40 CONTINUE C Once a point is found in the tiling, the corresponding other 3 points of C a rhombus are generated by adding the proper vectors if ((ABS(rx).LT.9.0).AND.(ABS(ry).LT.9.0)) then write(6,*) arx0,ary0 write(4,*) rx,ry,1.0 write(4,*) rx-EX(i),ry-EY(i),1.0 write(4,*) rx-EX(j),ry-EY(j),1.0 write(4,*) rx-EX(i)-EX(j),ry-EY(i)-EY(j),1.0 write(2,*) arxp,aryp,1.0 write(2,*) arxp-EX(i),aryp-EY(i),1.0 write(2,*) arxp-EX(j),aryp-EY(j),1.0 write(2,*) arxp-EX(i)-EX(j),aryp-EY(i)-EY(j),1.0 endif 30 CONTINUE endif 20 CONTINUE CLOSE(2) CLOSE(4) CLOSE(6) C************************************************************** C This feature makes a nice drawing of the family of lines C************************************************************** OPEN(8,FILE='lineas.dat',STATUS='UNKNOWN') AMUX=8.0 DO 120 i=1,5 C write(8,*) EX(i),EY(i) Do 119 n=1,AMUX DO 119 naux=1,1000 if (i.NE.1) then x=REAL(naux)*(0.01)-(AMUX/2.0) y=(REAL(n)-(AMUX/2.0)+GAMMA(i)-x*EX(i))/EY(i) else y=(REAL(naux)*(0.01)-(AMUX/2.0)+GAMMA(i)) x=REAL(n)-(AMUX/2.0)+GAMMA(i) endif write(8,*) x,y,1.0 119 CONTINUE 120 CONTINUE CLOSE(8) STOP END