/*************************************************************************** GenBracelet.cpp - description ------------------- begin : 9 March 2008 copyright : (C) 2008 by Ernesto Estevez Rams email : estevez@imre.oc.uh.cu last modified : 9 March 2008 ***************************************************************************/ /*************************************************************************** * * * This program is free software; you can redistribute it and/or modify * * it under the terms of the GNU General Public License as published by * * the Free Software Foundation; either version 2 of the License, or * * (at your option) any later version. * * * ***************************************************************************/ /***************************************************************************** * Based on the algorithm of Sawada: * Generating Bracelets in Constant Amortized Time. J. Sawada, SIAM J. Comput * * 31 (2001), pp. 259-268 * ******************************************************************************/ #include #include #include int a[101]; // buffer int n; // length of objects unsigned long nobject; //----------------------------------------------------- unsigned long LexicographicOrder(void) { int j=1; unsigned long lo=(2< 0 && a[pal+p/2-1]==1) return pal; } return -1; } //----------------------------------------------------------- int Pm6m2(int p) { int pal=0; if(n%2==0){ pal=palindrome(p, p-1); if(pal > 0 && a[pal+p/2-1]==0) return pal; } return -1; } //----------------------------------------------------------- void Print(int p, int ch, int ch2p, int ch3p) { int j=0; int pal=0; if(ch==5){ if (p == n) { nobject++; printf("%d ", LexicographicOrder()); if((pal=Pm3m1o(p)) < 0) if((pal=Pm3m1so(p)) < 0) if((pal=Pm3m1s(p)) < 0) if((pal=Pm6m2(p)) < 0){ pal=1; printf("P3m1 "); } else printf("P-6m2 "); else printf("P-3m1(s) "); else printf("P-3m1(so) "); else printf("P-3m1(o) "); for(j=1; j<=n; j++) printf("%d ",a[j]); printf("\n"); } else if (n == 3*p && ch3p != 5) { //Rombohedral nobject++; printf("%d ", LexicographicOrder()); if((pal=Rm3mo(p)) < 0) if((pal=Rm3ms(p)) < 0) if((pal=Rm3mso(p)) < 0){ pal=1; printf("R3m "); } else printf("R-3m(so) "); else printf("R-3 m(s) "); else printf("R-3m(o) "); for(j=1; j<=n; j++) printf("%d ",a[j]); printf("\n"); } else if (n == 2*p && ch2p != 5) { //Hexagonal with 6_3 nobject++; printf("%d ", LexicographicOrder()); if((pal=P63mmcs(p)) < 0) if((pal=P63mmco(p)) < 0){ pal=1; printf("P6_3mc "); } else printf("P6_3/mmc(o) "); else printf("P6_3/mmc(s) "); for(j=1; j<=n; j++) printf("%d ",a[j]); printf("\n"); } } } //----------------------------------------------------------------- int charge(int c, int op) { #define H 0 #define K 1 #define KK 2 #define KH 3 #define HK 4 #define E 5 // H K KK KH HK E static char oper[6][6]={{E, HK, KH, KK, K, H}, // H {KH, KK, E, HK, H, K}, // K {HK, E, K, H, KH, KK}, // KK {K, H, HK, E, KK, KH}, // KH {KK, KH, H, K, E, HK}, // HK {H, K, KK, KH, HK, E}}; // E int ch; ch=oper[c][op]; return ch; } //------------------------------------------------------------ int CheckRev(int t, int i) { int j=0; for(j=i+1; j<=(t+1)/2; ++j){ if(a[j] < a[t-j+1]) return 0; if(a[j] > a[t-j+1]) return -1; } return 1; } //--------------------------------------------------------------------------------------------- // The real vainille //--------------------------------------------------------------------------------------------- void GenNeutralBracelets(int t, int p, int r, int u, int v, int rs, int ch, int ch2p, int ch3p) { int j=0; int rev=0; int i=0; int ch1=0; int n2=n/2; int n3=n/3; if(t-1 > (n-r)/2+r){ if(a[t-1] > a[n-t+2+r]) rs=0; else if(a[t-1] < a[n-t+2+r]) rs=1; } if(t > n) { if(rs == 0 && n%p == 0) Print(p, ch, ch2p, ch3p); } else { a[t] = a[t-p]; ch1=ch; ch=charge(ch, a[t-p]); // keeping track of charge if(n%2==0 && t<=n2) // keeping track of charge for periodic bracelet ch2p=ch; if(n%3==0 && t<=n3) // keeping track of charge for periodic bracelet ch3p=ch; if(a[t] == a[1]) v++; else v=0; if(u == t-1 && a[t-1] == a[1]) u++; if(t == n && u != n && a[n] == a[1]){;} else if(u == v){ rev=CheckRev(t, u); if(rev == 0) GenNeutralBracelets(t+1, p, r, u, v, rs, ch, ch2p, ch3p); if(rev == 1) GenNeutralBracelets(t+1, p, t, u, v, 0, ch, ch2p, ch3p); } else GenNeutralBracelets(t+1, p, r, u, v, rs, ch, ch2p, ch3p); if(u == t) u--; if(a[t-p] == 0) { a[t] = 1; ch=charge(ch1, 1); if(n%2==0 && t<=n2) // keeping track of charge for periodic bracelet ch2p=ch; if(n%3==0 && t<=n3) // keeping track of charge for periodic bracelet ch3p=ch; if(t == 1) GenNeutralBracelets(t+1,t,r,1,1,rs, ch, ch2p, ch3p); else GenNeutralBracelets(t+1, t, r, u, 0, rs, ch, ch2p, ch3p); } } } //------------------------------------------------------------------------------- int main(int argc, char *argv[]) { int i=2; int l=0; if (argc<2){ printf("Polytype generation in the HK notation. Copyright GPL. IMRE (2008).\n\n"); printf("To few arguments:\n polytypes l\nwhere l stands for the maximum polytype length.\n\n"); return 0; } else l=atoi(argv[1]); printf("Polytype generation in the HK notation. Copyright GPL. IMRE (2008).\n\nOutput: lexicographic order Space group code\n0 stands for h, and 1 stands for k.\n\n"); for(i=2; i<=l; ++i){ n=i; a[0] = 0; nobject=0; GenNeutralBracelets(1,1,0,0,0,0,5,5,5); } }