1 /* twohamg.c */
2
3 #define USAGE "twohamg [-sgvq] [-L#] [infile [outfile]]"
4
5 #define HELPTEXT \
6 " Partition quartic graphs into two hamiltonian cycles.\n\
7 Output those which cannot be partitioned.\n\
8 \n\
9 -s force output to sparse6 format\n\
10 -g force output to graph6 format\n\
11 If neither -s or -g are given, the output format is\n\
12 determined by the header or, if there is none, by the\n\
13 format of the first input graph. Also see -S. \n\
14 \n\
15 The output file will have a header if and only if the input file does.\n\
16 \n\
17 -p Read a cubic graph and use its prism. Vertex i of the input becomes\n\
18 vertices 2*i,2*i+1 in the prism.\n\
19 -x Test for decompositions using each 2-path\n\
20 -X As -x but only output if two 2-paths are missed at some vertex\n\
21 -y Test for decompositions using each non-triangular 3-path\n\
22 -t# With -x and -X, consider only paths with center #\n\
23 With -y, consider only paths starting at #\n\
24 -Y With -p, only consider paths whose central edge is vertical\n\
25 -v Give a partition for those graphs who have one and a message\n\
26 for those which don't. With -x, list exceptional 2-paths.\n\
27 -L# Limit to 1000*# iterations; write with message if timeout.\n\
28 Graphs that time out are written to the output.\n\
29 \n\
30 -q suppress auxiliary information\n"
31
32 #define DEBUG 0
33 #define FISHTAIL 1 /* Whether to use fish-tail rule */
34
35 /*************************************************************************/
36
37 #include "gtools.h"
38
39 #define WHITE 0 /* unknown colour */
40 #define BLUE 1
41 #define RED 2
42
43 #define NO 0
44 #define YES 1
45 #define TIMEOUT 2
46 #define NOTHING 3
47
48 typedef struct addrval_struct
49 {
50 int *addr;
51 int value;
52 } addrval;
53
54 static long seed = 314159265;
55
56 /* valstack is a stack to save and restore degrees, farends, and colours */
57 DYNALLSTAT(addrval,valstack,valstack_sz);
58 static addrval *valstacktop;
59 #define SETVAL(ad,val) { valstacktop->addr = &(ad); \
60 valstacktop->value = ad; ++valstacktop; ad = val; }
61 #define UNSETVAL { --valstacktop; *(valstacktop->addr) = valstacktop->value; }
62
63 typedef struct p4
64 {
65 int e1,e2,e3;
66 int v1,v2,v3,v4;
67 boolean ok;
68 } p4;
69
70 DYNALLSTAT(int,bluefarend,bluefarend_sz); /* Parallel to sg.v */
71 DYNALLSTAT(int,redfarend,redfarend_sz); /* Parallel to sg.v */
72 DYNALLSTAT(int,eno,eno_sz); /* Edge numbers, parallel to sg.e */
73 DYNALLSTAT(int,cross,cross_sz); /* Parallel to sg.eno */
74 DYNALLSTAT(int,colour,colour_sz); /* indexed by edge number */
75 DYNALLSTAT(int,v1,v1_sz);
76 DYNALLSTAT(int,v2,v2_sz);
77 DYNALLSTAT(int,bluedeg,bluedeg_sz); /* Parallel to sg.v */
78 DYNALLSTAT(int,reddeg,reddeg_sz); /* Parallel to sg.v */
79 DYNALLSTAT(int,beste,beste_sz); /* List of possible best edges */
80 DYNALLSTAT(p4,p4list,p4list_sz); /* List of non-triangular p4s */
81
82 /* vstack is a stack of interesting vertices; onstack says whether a
83 vertex is on vstack so that we don't put it there twice */
84 DYNALLSTAT(int,vstack,vstack_sz);
85 DYNALLSTAT(boolean,onstack,onstack_sz);
86 static int *top;
87 #define PUSH(v) if (!onstack[v]) { *(top++) = (v); onstack[v] = TRUE; }
88 #define POP(v) { v = *(--top); onstack[v] = FALSE; }
89 #define STACKISEMPTY (top == vstack)
90
91 static nauty_counter nodes,limit,totallimit;
92 static nauty_counter locallimit[]
93 = {5,7,10,20,30,40,50,60,70,100,200,400,600,1000,2000,4000,8000,
94 20000,50000,100000,200000,5000000,0};
95 #define NUMLIMITS (sizeof(locallimit)/sizeof(locallimit[0]))
96 static nauty_counter nphases[100]; /* Must be greater than NUMLIMITS */
97 static nauty_counter ntimeouts;
98
99 /**************************************************************************/
100
101 static void
makeprism_sg(sparsegraph * sg,sparsegraph * sh)102 makeprism_sg(sparsegraph *sg, sparsegraph *sh)
103 /* Make the cartesian product sh = sg x K2.
104 sh must be initialized. */
105 {
106 size_t *v,*vv,vi,vvi;
107 int *d,*dd,*e,*ee;
108 int i,j,n;
109
110 SG_VDE(sg,v,d,e);
111 n = sg->nv;
112
113 SG_ALLOC(*sh,2*sg->nv,2*(n+sg->nde),"makeprism_sg");
114 sh->nv = 2*n;
115 sh->nde = 2*(n+sg->nde);
116 SG_VDE(sh,vv,dd,ee);
117
118 vvi = 0;
119 for (i = 0; i < n; ++i)
120 {
121 vi = v[i];
122 dd[2*i] = dd[2*i+1] = d[i] + 1;
123 vv[2*i] = vvi;
124 for (j = 0; j < d[i]; ++j) ee[vvi++] = 2*e[vi++];
125 ee[vvi++] = 2*i+1;
126 vi = v[i];
127 vv[2*i+1] = vvi;
128 for (j = 0; j < d[i]; ++j) ee[vvi++] = 2*e[vi++]+1;
129 ee[vvi++] = 2*i;
130 }
131 }
132
133 /**************************************************************************/
134
135 static void
dumpdata(int id,int nblue,int nred,int n)136 dumpdata(int id, int nblue, int nred, int n)
137 {
138 int i;
139
140 printf("%d: nblue=%d nred=%d -------------------------\n",id,nblue,nred);
141
142 #if DEBUG > 1
143 for (i = 0; i < n; ++i)
144 {
145 printf("%2d: ",i);
146 for (j = 0; j < 4; ++j)
147 printf(" %2d%c",eno[4*i+j],"wbr"[colour[eno[4*i+j]]]);
148 printf(" b=%d r=%d bfe=%d rfe=%d\n",
149 bluedeg[i],reddeg[i],bluefarend[i],redfarend[i]);
150 }
151 #endif
152 }
153
154 /**************************************************************************/
155
156 static void
initialise_g(int n,int * e)157 initialise_g(int n, int *e)
158 /* Allocate all and initialise eno[],v1[],v2[]
159 e is a vector like sg.e */
160 {
161 int ne,i,j,k,l;
162
163 DYNALLOC1(addrval,valstack,valstack_sz,10*n,"malloc");
164 DYNALLOC1(int,bluefarend,bluefarend_sz,n,"malloc");
165 DYNALLOC1(int,redfarend,redfarend_sz,n,"malloc");
166 DYNALLOC1(int,eno,eno_sz,4*n,"malloc");
167 DYNALLOC1(int,v1,v1_sz,2*n,"malloc");
168 DYNALLOC1(int,v2,v2_sz,2*n,"malloc");
169 DYNALLOC1(int,colour,colour_sz,2*n,"malloc");
170 DYNALLOC1(int,bluedeg,bluedeg_sz,n,"malloc");
171 DYNALLOC1(int,reddeg,reddeg_sz,n,"malloc");
172 DYNALLOC1(int,vstack,vstack_sz,n,"malloc");
173 DYNALLOC1(boolean,onstack,onstack_sz,n,"malloc");
174 DYNALLOC1(int,beste,beste_sz,2*n,"malloc");
175
176 /* Randomize e; seems to be no purpose for this any more. */
177
178 for (i = 0; i < n; ++i)
179 {
180 j = KRAN(4);
181 k = e[4*i+j]; e[4*i+j] = e[4*i+3]; e[4*i+3] = k;
182 j = KRAN(3);
183 k = e[4*i+j]; e[4*i+j] = e[4*i+2]; e[4*i+2] = k;
184 j = KRAN(2);
185 k = e[4*i+j]; e[4*i+j] = e[4*i+1]; e[4*i+1] = k;
186 }
187
188 ne = 0;
189 for (i = 0; i < n; ++i)
190 {
191 for (j = 0; j < 4; ++j)
192 {
193 k = e[4*i+j];
194 if (k > i)
195 {
196 v1[ne] = i;
197 v2[ne] = k;
198 eno[4*i+j] = ne++;
199 }
200 else /* Note: this code assumes a simple graph */
201 {
202 for (l = 0; l < 4; ++l)
203 if (e[4*k+l] == i) break;
204 eno[4*i+j] = eno[4*k+l];
205 }
206 }
207 }
208 if (ne != 2*n) gt_abort(">E ne is incorrect");
209
210 #if DEBUG
211 { int ii;
212 printf("===== n=%d === ne=%d ===================\n",n,ne);
213
214 for (ii = 0; ii < n; ++ii)
215 printf("%2d: %2d %2d %2d %2d %2d %2d %2d %2d\n",
216 ii,e[4*ii],e[4*ii+1],e[4*ii+2],e[4*ii+3],
217 eno[4*ii],eno[4*ii+1],eno[4*ii+2],eno[4*ii+3]);
218 }
219 #endif
220 }
221
222 /**************************************************************************/
223
224 static void
initialise_colouring(int n)225 initialise_colouring(int n)
226 /* Make all edges white */
227 {
228 int i;
229
230 for (i = 0; i < 2*n; ++i) colour[i] = WHITE;
231
232 for (i = 0; i < n; ++i)
233 {
234 bluefarend[i] = i;
235 redfarend[i] = i;
236 bluedeg[i] = 0;
237 reddeg[i] = 0;
238 onstack[i] = FALSE;
239 }
240
241 top = vstack;
242 valstacktop = valstack;
243 }
244
245 /**************************************************************************/
246
247 static boolean
makeblue(int edge,boolean lastok)248 makeblue(int edge, boolean lastok)
249 /* Colour WHITE edge BLUE, return success.
250 lastok indicates if it is ok to add the final blue edge */
251 {
252 int w1,w2,f1,f2;
253
254 w1 = v1[edge];
255 w2 = v2[edge];
256 if (bluedeg[w1] == 2 || bluedeg[w2] == 2) return FALSE;
257
258 f1 = bluefarend[w1];
259 f2 = bluefarend[w2];
260 if (f1 == w2 && !lastok) return FALSE;
261
262 SETVAL(colour[edge],BLUE);
263 SETVAL(bluedeg[w1],bluedeg[w1]+1);
264 SETVAL(bluedeg[w2],bluedeg[w2]+1);
265 SETVAL(bluefarend[f1],f2);
266 SETVAL(bluefarend[f2],f1);
267 PUSH(w1);
268 PUSH(w2);
269 PUSH(f1);
270 // PUSH(f2); not sure this one is needed
271
272 return TRUE;
273 }
274
275 /**************************************************************************/
276
277 static boolean
makered(int edge,boolean lastok)278 makered(int edge, boolean lastok)
279 /* Colour WHITE edge RED, return success.
280 lastok indicates if it is ok to add the final red edge */
281 {
282 int w1,w2,f1,f2;
283
284 w1 = v1[edge];
285 w2 = v2[edge];
286 if (reddeg[w1] == 2 || reddeg[w2] == 2) return FALSE;
287
288 f1 = redfarend[w1];
289 f2 = redfarend[w2];
290 if (f1 == w2 && !lastok) return FALSE;
291
292 SETVAL(colour[edge],RED);
293 SETVAL(reddeg[w1],reddeg[w1]+1);
294 SETVAL(reddeg[w2],reddeg[w2]+1);
295 SETVAL(redfarend[f1],f2);
296 SETVAL(redfarend[f2],f1);
297 PUSH(w1);
298 PUSH(w2);
299 PUSH(f1);
300 // PUSH(f2); not sure this one is needed
301
302 return TRUE;
303 }
304
305 /**************************************************************************/
306
307 static boolean
propagate(int n,int * e,int * nblue,int * nred)308 propagate(int n, int *e, int *nblue, int *nred)
309 /* Look at active vertices and propagate colourings */
310 {
311 int i,j,v,w;
312
313 while (!STACKISEMPTY)
314 {
315 POP(v);
316
317 if (reddeg[v] == 2 && bluedeg[v] < 2)
318 {
319 for (i = 0; i < 4; ++i)
320 {
321 j = eno[4*v+i];
322 if (colour[j] == WHITE)
323 {
324 if (!makeblue(j,*nblue==n-1)) return FALSE;
325 ++*nblue;
326 }
327 }
328 }
329 else if (bluedeg[v] == 2 && reddeg[v] < 2)
330 {
331 for (i = 0; i < 4; ++i)
332 {
333 j = eno[4*v+i];
334 if (colour[j] == WHITE)
335 {
336 if (!makered(j,*nred==n-1)) return FALSE;
337 ++*nred;
338 }
339 }
340 }
341
342 if (bluedeg[v] == 1)
343 {
344 w = bluefarend[v];
345 for (i = 0; i < 4; ++i)
346 if (e[4*v+i] == w) break;
347 if (i < 4)
348 {
349 j = eno[4*v+i];
350 if (colour[j] == WHITE)
351 {
352 if (*nblue == n-1)
353 {
354 if (!makeblue(j,TRUE)) return FALSE;
355 ++*nblue;
356 }
357 else
358 {
359 if (!makered(j,*nred==n-1)) return FALSE;
360 ++*nred;
361 }
362 }
363 }
364 }
365
366 if (reddeg[v] == 1)
367 {
368 w = redfarend[v];
369 for (i = 0; i < 4; ++i)
370 if (e[4*v+i] == w) break;
371 if (i < 4)
372 {
373 j = eno[4*v+i];
374 if (colour[j] == WHITE)
375 {
376 if (*nred == n-1)
377 {
378 if (!makered(j,TRUE)) return FALSE;
379 ++*nred;
380 }
381 else
382 {
383 if (!makeblue(j,*nblue==n-1)) return FALSE;
384 ++*nblue;
385 }
386 }
387 }
388 }
389 }
390
391 return TRUE;
392 }
393
394 /**************************************************************************/
395
396 static int
fishtail(int n,int * nblue,int * nred)397 fishtail(int n, int *nblue, int *nred)
398 /* Try to colour edges by the fishtail method
399 Return NO (contradiction), YES (did some good), or NOTHING. */
400 {
401 int i,j;
402 int e1,w1,e2,w2,e3,w3;
403 int edge,col;
404 int ans;
405
406 ans = NOTHING;
407
408 for (i = 0; i < n; ++i)
409 {
410 edge = -1;
411
412 if (reddeg[i] == 1 && bluedeg[i] == 0 && *nblue < n-2)
413 {
414 j = 4*i;
415 if (colour[eno[j]] != WHITE) ++j;
416 e1 = eno[j];
417 w1 = (v1[e1] == i ? v2[e1] : v1[e1]);
418 ++j;
419 if (colour[eno[j]] != WHITE) ++j;
420 e2 = eno[j];
421 w2 = (v1[e2] == i ? v2[e2] : v1[e2]);
422 ++j;
423 if (colour[eno[j]] != WHITE) ++j;
424 e3 = eno[j];
425 w3 = (v1[e3] == i ? v2[e3] : v1[e3]);
426
427 if (bluedeg[w1] == 1 && bluefarend[w1] == w2)
428 {
429 edge = e3;
430 col = BLUE;
431 }
432 else if (bluedeg[w1] == 1 && bluefarend[w1] == w3)
433 {
434 edge = e2;
435 col = BLUE;
436 }
437 else if (bluedeg[w2] == 1 && bluefarend[w2] == w3)
438 {
439 edge = e1;
440 col = BLUE;
441 }
442 }
443 else if (reddeg[i] == 0 && bluedeg[i] == 1 && *nred < n-2)
444 {
445 j = 4*i;
446 if (colour[eno[j]] != WHITE) ++j;
447 e1 = eno[j];
448 w1 = (v1[e1] == i ? v2[e1] : v1[e1]);
449 ++j;
450 if (colour[eno[j]] != WHITE) ++j;
451 e2 = eno[j];
452 w2 = (v1[e2] == i ? v2[e2] : v1[e2]);
453 ++j;
454 if (colour[eno[j]] != WHITE) ++j;
455 e3 = eno[j];
456 w3 = (v1[e3] == i ? v2[e3] : v1[e3]);
457
458 if (reddeg[w1] == 1 && redfarend[w1] == w2)
459 {
460 edge = e3;
461 col = RED;
462 }
463 else if (reddeg[w1] == 1 && redfarend[w1] == w3)
464 {
465 edge = e2;
466 col = RED;
467 }
468 else if (reddeg[w2] == 1 && redfarend[w2] == w3)
469 {
470 edge = e1;
471 col = RED;
472 }
473 }
474 if (edge >= 0)
475 {
476 if (col == BLUE)
477 {
478 if (!makeblue(edge,*nblue == n-1)) return NO;
479 ++*nblue;
480 }
481 else
482 {
483 if (!makered(edge,*nred == n-1)) return NO;
484 ++*nred;
485 }
486 ans = YES;
487 #if DEBUG
488 printf("Fishtail: %d=(%d-%d) -> %c\n",
489 edge,v1[edge],v2[edge],"wbr"[col]);
490 #endif
491 }
492 }
493
494 return ans;
495 }
496
497 /**************************************************************************/
498
499 static int
searchnode(int level,int n,int * e,int nblue,int nred)500 searchnode(int level, int n, int *e, int nblue, int nred)
501 {
502 boolean ok;
503 int i,status,nbest;
504 addrval *valptr;
505 int best,score,bestscore;
506 long ran;
507
508 ok = propagate(n,e,&nblue,&nred);
509
510 #if DEBUG
511 if (ok) dumpdata(level,nblue,nred,n);
512 else { printf("FAIL\n"); return NO; }
513 #else
514 if (!ok) return NO;
515 #endif
516
517 if (nblue == n && nred == n) return YES;
518
519 #if FISHTAIL
520 status = fishtail(n,&nblue,&nred);
521 if (status == NO) return NO;
522 if (status == YES) return searchnode(level+1,n,e,nblue,nred);
523 #endif
524
525 valptr = valstacktop;
526
527 bestscore = -1;
528 nbest = 0;
529 for (i = 0; i < 2*n; ++i)
530 if (colour[i] == WHITE)
531 {
532 score = (bluedeg[v1[i]] == 1) + (bluedeg[v2[i]] == 1) +
533 (reddeg[v1[i]] == 1) + (reddeg[v2[i]] == 1);
534 if (score > bestscore)
535 {
536 bestscore = score;
537 beste[0] = i;
538 nbest = 1;
539 }
540 else if (score == bestscore)
541 beste[nbest++] = i;
542 }
543
544 if (bestscore == 0 && nred + nblue > 0)
545 return FALSE; /* Disconnected */
546
547 ran = KRAN(2*nbest);
548 best = beste[ran/2];
549
550 if ((ran&1))
551 {
552 if (makeblue(best,nblue==n-1))
553 {
554 #if DEBUG
555 printf(" setting %d(%d-%d) blue\n",best,v1[best],v2[best]);
556 #endif
557 status = searchnode(level+1,n,e,nblue+1,nred);
558 if (status != NO) return status;
559 while (valstacktop > valptr) UNSETVAL;
560 }
561
562 if (++nodes == limit) return TIMEOUT;
563
564 if (makered(best,nred==n-1))
565 {
566 #if DEBUG
567 printf(" setting %d(%d-%d) red\n",best,v1[best],v2[best]);
568 #endif
569 status = searchnode(level+1,n,e,nblue,nred+1);
570 if (status != NO) return status;
571 while (valstacktop > valptr) UNSETVAL;
572 }
573 }
574 else
575 {
576 if (makered(best,nred==n-1))
577 {
578 #if DEBUG
579 printf(" setting %d(%d-%d) red\n",best,v1[best],v2[best]);
580 #endif
581 status = searchnode(level+1,n,e,nblue,nred+1);
582 if (status != NO) return status;
583 while (valstacktop > valptr) UNSETVAL;
584 }
585
586 if (++nodes == limit) return TIMEOUT;
587
588 if (makeblue(best,nblue==n-1))
589 {
590 #if DEBUG
591 printf(" setting %d(%d-%d) blue\n",best,v1[best],v2[best]);
592 #endif
593 status = searchnode(level+1,n,e,nblue+1,nred);
594 if (status != NO) return status;
595 while (valstacktop > valptr) UNSETVAL;
596 }
597 }
598
599 return NO;
600 }
601
602 /**************************************************************************/
603
604 static int
dispatchsearch(int n,int * e,int nblue,int nred)605 dispatchsearch(int n, int *e, int nblue, int nred)
606 /* Multiple tries at solution */
607 {
608 int i,status;
609 addrval *valptr;
610 boolean ok;
611 nauty_counter remaininglimit;
612
613 ok = propagate(n,e,&nblue,&nred);
614
615 #if DEBUG
616 if (ok) dumpdata(0,nblue,nred,n);
617 else { printf("FAIL\n"); ++nphases[0]; return NO; }
618 #else
619 if (!ok) { ++nphases[0]; return NO; }
620 #endif
621
622 valptr = valstacktop;
623
624 remaininglimit = totallimit;
625
626 for (i = 0; i < NUMLIMITS; ++i)
627 {
628 if (totallimit > 0)
629 {
630 limit = locallimit[i];
631 if (limit == 0 || limit > remaininglimit) limit = remaininglimit;
632 remaininglimit -= limit;
633 if (limit == 0) limit = 1;
634 }
635 else
636 limit = locallimit[i];
637
638 nodes = 0;
639
640 status = searchnode(1,n,e,nblue,nred);
641 if (status != TIMEOUT)
642 {
643 ++nphases[i+1];
644 return status;
645 }
646 while (valstacktop > valptr) UNSETVAL;
647 }
648
649 ++ntimeouts;
650 return TIMEOUT;
651 }
652
653 /**************************************************************************/
654
655 static int
isdecomposable(sparsegraph sg,int * ham1,int * ham2)656 isdecomposable(sparsegraph sg, int *ham1, int *ham2)
657 /* Test if sg is decomposable as two hamiltonian cycles */
658 {
659 int nblue,nred,n,e0,status;
660 int i0,i,j,lastv,laste;
661
662 n = sg.nv;
663 initialise_g(n,sg.e);
664 initialise_colouring(n);
665
666 e0 = KRAN(2*n);
667
668 nblue = nred = 0;
669 if (!makeblue(e0,nblue==n-1)) /* make edge e0 blue */
670 return NO;
671 ++nblue;
672
673 #if DEBUG
674 printf(" setting %d(%d-%d) blue\n",e0,v1[0],v2[0]);
675 #endif
676
677 // status = searchnode(1,n,sg.e,nblue,nred);
678 status = dispatchsearch(n,sg.e,nblue,nred);
679
680 if (status != YES) return status;
681
682 if (!ham1 || !ham2) return YES;
683
684 for (i0 = 0; colour[i0] != BLUE; ++i0) {}
685 ham1[0] = v1[i0];
686 laste = i0;
687 lastv = ham1[1] = v2[i0];
688 for (i = 2; i < n; ++i)
689 {
690 for (j = 0; j < 4; ++j)
691 if (eno[4*lastv+j] != laste && colour[eno[4*lastv+j]] == BLUE)
692 break;
693 laste = eno[4*lastv+j];
694 lastv = ham1[i] = (v1[laste] == lastv ? v2[laste] : v1[laste]);
695 }
696
697 for (i0 = 0; colour[i0] != RED; ++i0) {}
698 ham2[0] = v1[i0];
699 laste = i0;
700 lastv = ham2[1] = v2[i0];
701 for (i = 2; i < n; ++i)
702 {
703 for (j = 0; j < 4; ++j)
704 if (eno[4*lastv+j] != laste && colour[eno[4*lastv+j]] == RED)
705 break;
706 laste = eno[4*lastv+j];
707 lastv = ham2[i] = (v1[laste] == lastv ? v2[laste] : v1[laste]);
708 }
709
710 return YES;
711 }
712
713 /**************************************************************************/
714
715 static int
iscrossdecomposable(sparsegraph sg,int vertex)716 iscrossdecomposable(sparsegraph sg, int vertex)
717 /* Test if sg is decomposable as two hamiltonian cycles
718 in all ways through each vertex (or the given vertex if it
719 is >= 0. Details left in cross[].
720 Return -1: not decomposable at all
721 0: misses some 2-paths including two at some vertices
722 1: misses some 2-paths but at most one per vertex
723 2: fully decomposable
724 3: timeout
725 */
726 {
727 int i,j,k,l,n,c;
728 int ans,status;
729 int imin,imax;
730
731 n = sg.nv;
732 DYNALLOC1(int,cross,cross_sz,4*n,"malloc");
733
734 if (vertex >= n) gt_abort(">E vertex given to -t is too large");
735
736 for (i = 0; i < 4*n; ++i) cross[i] = 0;
737
738 status = isdecomposable(sg,FALSE,FALSE);
739 if (status == NO) return -1;
740 else if (status == TIMEOUT) return 3;
741
742 for (i = 0; i < n; ++i)
743 {
744 c = colour[eno[4*i]];
745 for (j = 1; j < 4; ++j)
746 if (colour[eno[4*i+j]] == c) cross[4*i+j] = 1;
747 }
748
749 ans = 2;
750
751 if (vertex >= 0) { imin = imax = vertex; }
752 else { imin = 0; imax = n-1; }
753
754 for (i = imin; i <= imax; ++i)
755 for (j = 1; j < 4; ++j)
756 if (!cross[4*i+j])
757 {
758 initialise_colouring(n);
759 if (makeblue(eno[4*i],FALSE) && makeblue(eno[4*i+j],n==2))
760 {
761 status = dispatchsearch(n,sg.e,2,0);
762 if (status == TIMEOUT) return 3;
763 if (status == YES)
764 {
765 for (k = 0; k < n; ++k)
766 {
767 c = colour[eno[4*k]];
768 for (l = 1; l < 4; ++l)
769 if (colour[eno[4*k+l]] == c) cross[4*k+l] = 1;
770 }
771 }
772 else
773 ans = 1;
774 }
775 else
776 ans = 1;
777 }
778
779 if (ans == 1)
780 for (i = 0; i < n; ++i)
781 {
782 if (cross[4*i+1] + cross[4*i+2] + cross[4*i+3] <= 1)
783 ans = 0;
784 }
785
786 return ans;
787 }
788
789 /**************************************************************************/
790
791 static int
p4decomposition(sparsegraph sg,int vertex,boolean vertical)792 p4decomposition(sparsegraph sg, int vertex, boolean vertical)
793 /* Test which non-triangular P4s extend to a decomposition.
794 Return -2: timeout
795 -1: not decomposable at all
796 >=0: number of P4s missed
797 If 0 <= vertex < n, only P4s starting at that vertex are considered.
798 If vertical (only for prisms), the central edge must be vertical.
799 The paths missed can be found in p4list[].
800 */
801 {
802 int i,k,n;
803 int ans,status;
804 int imin,imax,nump4;
805 int v1,v2,v3,v4,e1,e2,e3,j1,j2,j3;
806
807 n = sg.nv;
808 if (vertex >= n) gt_abort(">E vertex given to -t is too large");
809
810 status = isdecomposable(sg,FALSE,FALSE);
811 if (status == NO) return -1;
812 else if (status == TIMEOUT) return -2;
813
814 if (vertex >= 0)
815 {
816 imin = imax = vertex;
817 DYNALLOC1(p4,p4list,p4list_sz,36,"malloc");
818 }
819 else
820 {
821 imin = 0; imax = n-1;
822 DYNALLOC1(p4,p4list,p4list_sz,18*n,"malloc");
823 }
824
825 nump4 = 0;
826 for (v1 = imin; v1 <= imax; ++v1)
827 for (j1 = 0; j1 < 4; ++j1)
828 {
829 v2 = sg.e[4*v1+j1];
830 e1 = eno[4*v1+j1];
831 for (j2 = 0; j2 < 4; ++j2)
832 {
833 v3 = sg.e[4*v2+j2];
834 if (v3 == v1) continue;
835 if (vertical && (v2|1) != (v3|1)) continue;
836 e2 = eno[4*v2+j2];
837 for (j3 = 0; j3 < 4; ++j3)
838 {
839 v4 = sg.e[4*v3+j3];
840 if (v4 == v1 || v4 == v2) continue;
841 e3 = eno[4*v3+j3];
842 if (vertex >= 0 || v1 < v4)
843 {
844 p4list[nump4].v1 = v1;
845 p4list[nump4].v2 = v2;
846 p4list[nump4].v3 = v3;
847 p4list[nump4].v4 = v4;
848 p4list[nump4].e1 = e1;
849 p4list[nump4].e2 = e2;
850 p4list[nump4].e3 = e3;
851 p4list[nump4].ok = FALSE;
852 ++nump4;
853 }
854 }
855 }
856 }
857
858 for (i = 0; i < nump4; ++i)
859 if (colour[p4list[i].e1] == colour[p4list[i].e2] &&
860 colour[p4list[i].e1] == colour[p4list[i].e3])
861 p4list[i].ok = TRUE;
862
863 for (i = 0; i < nump4; ++i)
864 if (!p4list[i].ok)
865 {
866 initialise_colouring(n);
867 if (makeblue(p4list[i].e1,FALSE) &&
868 makeblue(p4list[i].e2,n==2) &&
869 makeblue(p4list[i].e3,n==3))
870 {
871 status = dispatchsearch(n,sg.e,3,0);
872 if (status == TIMEOUT) return -2;
873 if (status == YES)
874 {
875 for (k = 0; k < nump4; ++k)
876 {
877 if (p4list[k].ok) continue;
878 if (colour[p4list[k].e1] == colour[p4list[k].e2]
879 && colour[p4list[k].e1] == colour[p4list[k].e3])
880 p4list[k].ok = TRUE;
881 }
882 }
883 }
884 }
885
886 ans = 0;
887 for (i = 0; i < nump4; ++i) if (!p4list[i].ok) ++ans;
888
889 return ans;
890 }
891
892 /**************************************************************************/
893
894 int
main(int argc,char * argv[])895 main(int argc, char *argv[])
896 {
897 sparsegraph sg,sh,*this;
898 int n,codetype;
899 int argnum,i,j,outcode;
900 char *arg,sw;
901 boolean badargs;
902 boolean sswitch,gswitch,qswitch,vswitch,xswitch,Xswitch;
903 boolean pswitch,Lswitch,tswitch,yswitch,Yswitch;
904 long Lvalue;
905 double t;
906 char *infilename,*outfilename;
907 FILE *infile,*outfile;
908 nauty_counter nin,nout;
909 int status,xstatus,ystatus,tvalue,imin,imax;
910 DYNALLSTAT(int,ham1,ham1_sz);
911 DYNALLSTAT(int,ham2,ham2_sz);
912
913 HELP; PUTVERSION;
914
915 INITSEED;
916 ran_init(seed);
917
918 sswitch = gswitch = yswitch = qswitch = vswitch = FALSE;
919 tswitch = xswitch = Xswitch = Lswitch = pswitch = FALSE;
920 Yswitch = FALSE;
921 infilename = outfilename = NULL;
922
923 argnum = 0;
924 badargs = FALSE;
925 for (j = 1; !badargs && j < argc; ++j)
926 {
927 arg = argv[j];
928 if (arg[0] == '-' && arg[1] != '\0')
929 {
930 ++arg;
931 while (*arg != '\0')
932 {
933 sw = *arg++;
934 SWBOOLEAN('s',sswitch)
935 else SWBOOLEAN('g',gswitch)
936 else SWBOOLEAN('q',qswitch)
937 else SWBOOLEAN('x',xswitch)
938 else SWBOOLEAN('y',yswitch)
939 else SWBOOLEAN('Y',Yswitch)
940 else SWBOOLEAN('X',Xswitch)
941 else SWBOOLEAN('v',vswitch)
942 else SWBOOLEAN('p',pswitch)
943 else SWLONG('L',Lswitch,Lvalue,"-L")
944 else SWINT('t',tswitch,tvalue,"-t")
945 else badargs = TRUE;
946 }
947 }
948 else
949 {
950 ++argnum;
951 if (argnum == 1) infilename = arg;
952 else if (argnum == 2) outfilename = arg;
953 else badargs = TRUE;
954 }
955 }
956
957 if (sswitch && gswitch)
958 gt_abort(">E twohamg: -s and -g are incompatible");
959 if ((xswitch != 0) + (Xswitch != 0)
960 + (yswitch != 0) + (Yswitch != 0) > 1)
961 gt_abort(">E twohamg: -x, -X, -y, -Y are incompatible");
962 if (Yswitch && !pswitch)
963 gt_abort(">E twohang: -Y is not allowed without -p");
964
965 if (!Lswitch) totallimit = 0;
966 else totallimit = Lvalue * (nauty_counter)1000;
967
968 if (!tswitch) tvalue = -1;
969
970 if (badargs || argnum > 2)
971 {
972 fprintf(stderr,">E Usage: %s\n",USAGE);
973 GETHELP;
974 exit(1);
975 }
976
977 if (!qswitch)
978 {
979 fprintf(stderr,">A twohamg");
980 if (sswitch || gswitch || vswitch || xswitch || Xswitch
981 || tswitch || pswitch || Lswitch || yswitch || Yswitch)
982 fprintf(stderr," -");
983 if (sswitch) fprintf(stderr,"s");
984 if (gswitch) fprintf(stderr,"g");
985 if (xswitch) fprintf(stderr,"x");
986 if (Xswitch) fprintf(stderr,"X");
987 if (yswitch) fprintf(stderr,"y");
988 if (Yswitch) fprintf(stderr,"Y");
989 if (pswitch) fprintf(stderr,"p");
990 if (vswitch) fprintf(stderr,"v");
991 if (Lswitch) fprintf(stderr,"L%ld",Lvalue);
992 if (tswitch) fprintf(stderr,"t%d",tvalue);
993 if (argnum > 0) fprintf(stderr," %s",infilename);
994 if (argnum > 1) fprintf(stderr," %s",outfilename);
995 fprintf(stderr,"\n");
996 fflush(stderr);
997 }
998
999 if (infilename && infilename[0] == '-') infilename = NULL;
1000 infile = opengraphfile(infilename,&codetype,FALSE,1);
1001 if (!infile) exit(1);
1002 if (!infilename) infilename = "stdin";
1003
1004 if (!outfilename || outfilename[0] == '-')
1005 {
1006 outfilename = "stdout";
1007 outfile = stdout;
1008 }
1009 else if ((outfile = fopen(outfilename,"w")) == NULL)
1010 {
1011 fprintf(stderr,"Can't open output file %s\n",outfilename);
1012 gt_abort(NULL);
1013 }
1014
1015 NODIGRAPHSYET(codetype);
1016
1017 if (sswitch || (!gswitch && (codetype&SPARSE6)))
1018 outcode = SPARSE6;
1019 else
1020 outcode = GRAPH6;
1021
1022 if (codetype&HAS_HEADER)
1023 {
1024 if (outcode == SPARSE6) writeline(outfile,SPARSE6_HEADER);
1025 else writeline(outfile,GRAPH6_HEADER);
1026 }
1027
1028 nin = nout = 0;
1029 t = CPUTIME;
1030 SG_INIT(sg);
1031 SG_INIT(sh);
1032
1033 while (read_sg(infile,&sg) != NULL)
1034 {
1035 ++nin;
1036 if (pswitch)
1037 {
1038 n = sg.nv;
1039 for (i = 0; i < n; ++i)
1040 if (sg.d[i] != 3) break;
1041 if (i != n)
1042 {
1043 fprintf(stderr,
1044 ">W input " COUNTER_FMT " is not cubic\n",nin);
1045 continue;
1046 }
1047 makeprism_sg(&sg,&sh);
1048 n = sh.nv;
1049 this = &sh;
1050 }
1051 else
1052 {
1053 n = sg.nv;
1054 for (i = 0; i < n; ++i)
1055 if (sg.d[i] != 4) break;
1056 if (i != n)
1057 {
1058 fprintf(stderr,
1059 ">W input " COUNTER_FMT " is not quartic\n",nin);
1060 continue;
1061 }
1062 for (i = 0; i < n; ++i)
1063 if (sg.v[i] != 4*i) break;
1064 if (i != n)
1065 gt_abort("readg6_sg() result is not in standard form");
1066 this = &sg;
1067 }
1068
1069 if (xswitch || Xswitch)
1070 {
1071 xstatus = iscrossdecomposable(*this,tvalue);
1072 if ((xswitch && xstatus < 2) || (Xswitch && xstatus < 1))
1073 {
1074 if (outcode == SPARSE6) writes6_sg(outfile,&sg);
1075 else if (outcode == GRAPH6) writeg6_sg(outfile,&sg);
1076 ++nout;
1077
1078 if (vswitch)
1079 {
1080 if (xstatus < 0)
1081 {
1082 fprintf(stderr,">H Graph " COUNTER_FMT
1083 " is indecomposable\n",nin);
1084 }
1085 else
1086 {
1087 fprintf(stderr,">X" COUNTER_FMT ": ",nin);
1088 if (tswitch) { imin = imax = tvalue; }
1089 else { imin = 0; imax = n-1; }
1090 for (i = imin; i <= imax; ++i)
1091 for (j = 1; j < 4; ++j)
1092 if (!cross[4*i+j])
1093 {
1094 fprintf(stderr," %d-%d-%d",
1095 (v1[eno[4*i]] == i
1096 ? v2[eno[4*i]] : v1[eno[4*i]]),
1097 i,
1098 (v1[eno[4*i+j]] == i
1099 ? v2[eno[4*i+j]] : v1[eno[4*i+j]]));
1100 }
1101 fprintf(stderr,"\n");
1102 }
1103 }
1104 }
1105 else if (xstatus == 3)
1106 {
1107 fprintf(stderr,">H Graph " COUNTER_FMT
1108 " timed out\n",nin);
1109 if (outcode == SPARSE6) writes6_sg(outfile,&sg);
1110 else if (outcode == GRAPH6) writeg6_sg(outfile,&sg);
1111 ++nout;
1112 }
1113 }
1114 else if (yswitch || Yswitch)
1115 {
1116 ystatus = p4decomposition(*this,tvalue,Yswitch);
1117 if (ystatus != 0)
1118 {
1119 if (outcode == SPARSE6) writes6_sg(outfile,&sg);
1120 else if (outcode == GRAPH6) writeg6_sg(outfile,&sg);
1121 ++nout;
1122 }
1123
1124 if (ystatus == -2)
1125 {
1126 fprintf(stderr,">H Graph " COUNTER_FMT
1127 " timed out\n",nin);
1128 }
1129 else if (vswitch)
1130 {
1131 if (ystatus == -1)
1132 {
1133 fprintf(stderr,">H Graph " COUNTER_FMT
1134 " is indecomposable\n",nin);
1135 }
1136 else if (ystatus > 0)
1137 {
1138 fprintf(stderr,">X" COUNTER_FMT ": ",nin);
1139 for (i = j = 0; j < ystatus; ++i)
1140 if (!p4list[i].ok)
1141 {
1142 ++j;
1143 fprintf(stderr," %d-%d-%d-%d",
1144 p4list[i].v1,p4list[i].v2,
1145 p4list[i].v3,p4list[i].v4);
1146 }
1147 fprintf(stderr,"\n");
1148 }
1149 }
1150 }
1151 else
1152 {
1153 if (vswitch)
1154 {
1155 DYNALLOC1(int,ham1,ham1_sz,n,"malloc");
1156 DYNALLOC1(int,ham2,ham2_sz,n,"malloc");
1157 status = isdecomposable(*this,ham1,ham2);
1158 }
1159 else
1160 status = isdecomposable(*this,NULL,NULL);
1161
1162 if (status != YES)
1163 {
1164 if (outcode == SPARSE6) writes6_sg(outfile,&sg);
1165 else if (outcode == GRAPH6) writeg6_sg(outfile,&sg);
1166
1167 if (status == TIMEOUT)
1168 fprintf(stderr,">H Graph " COUNTER_FMT
1169 " timed out\n",nin);
1170 else if (vswitch)
1171 fprintf(stderr,">H Graph " COUNTER_FMT
1172 " is indecomposable\n",nin);
1173 ++nout;
1174 }
1175 else if (vswitch)
1176 {
1177 #if DEBUG
1178 dumpdata(0,0,0,n);
1179 #endif
1180 fprintf(stderr,">H" COUNTER_FMT ": ",nin);
1181 for (i = 0; i < n; ++i) fprintf(stderr," %d",ham1[i]);
1182 fprintf(stderr,"\n ");
1183 for (i = 0; i < n; ++i) fprintf(stderr," %d",ham2[i]);
1184 fprintf(stderr,"\n");
1185 }
1186 }
1187 }
1188 t = CPUTIME - t;
1189
1190 if (!qswitch)
1191 {
1192 fprintf(stderr,">T to=" COUNTER_FMT " phases=",ntimeouts);
1193 for (i = 0; i <= NUMLIMITS; ++i)
1194 fprintf(stderr," " COUNTER_FMT,nphases[i]);
1195 fprintf(stderr,"\n");
1196
1197 fprintf(stderr,
1198 ">Z " COUNTER_FMT " graphs read from %s, "
1199 COUNTER_FMT " written to %s; %3.2f sec.\n",
1200 nin,infilename,nout,outfilename,t);
1201 }
1202
1203 exit(0);
1204 }
1205