1 /*
2 * pstep-recursive.c
3 *
4 *
5 * Part of TREE-PUZZLE 5.2 (July 2004)
6 *
7 * (c) 2003-2004 by Heiko A. Schmidt, Korbinian Strimmer, and Arndt von Haeseler
8 * (c) 1999-2003 by Heiko A. Schmidt, Korbinian Strimmer,
9 * M. Vingron, and Arndt von Haeseler
10 * (c) 1995-1999 by Korbinian Strimmer and Arndt von Haeseler
11 *
12 * All parts of the source except where indicated are distributed under
13 * the GNU public licence. See http://www.opensource.org for details.
14 *
15 * ($Id$)
16 *
17 */
18
19 #ifdef HAVE_CONFIG_H
20 # include <config.h>
21 #endif
22
23 #define EXTERN extern
24
25 #include<stdio.h>
26 #include<stdlib.h>
27 #include"puzzle.h"
28 /* #include"util.h" */
29 /* #include"pstep.h" */
30
31
32 /******************************************************************************/
33 /* functions for representing and building puzzling step trees */
34 /******************************************************************************/
35
36 /* initialize tree with the following starting configuration
37 * A(=0)
38 * [up =NULL]
39 * [downleft =1]
40 * [downright=2]
41 *
42 * o
43 * |
44 * |
45 * +-------------+--------------+
46 * | |
47 * | |
48 * o o
49 *
50 * C(=1) B(=2)
51 * [up =0] [up =0]
52 * [downleft =NULL] [downleft =NULL]
53 * [downright=NULL] [downright=NULL]
54 *
55 * nextedge = 3
56 * nextleaf = 3
57 * and set according edge maps
58 *
59 */
inittree_orig(ONEEDGE_ORIG ** edge,int ** edgeofleaf,int Maxspc,int Maxbrnch,int * nextedge,int * nextleaf)60 void inittree_orig(ONEEDGE_ORIG **edge, /* out: new array of edges */
61 int **edgeofleaf,/* out: array of external edge ptrs */
62 int Maxspc, /* in: Number of species (n) */
63 int Maxbrnch, /* in: Number of branches (2n-3) */
64 int *nextedge, /* out: next free edge index (=3) */
65 int *nextleaf) /* out: next free leaf index (=3) */
66 {
67 int i;
68 ONEEDGE_ORIG *tmpedge;
69 int *tmpedgeofleaf;
70
71 /* allocate the memory for the whole tree */
72
73 /* allocate memory for vector with all the edges of the tree */
74 tmpedge = (ONEEDGE_ORIG *) calloc((size_t) Maxbrnch, sizeof(ONEEDGE_ORIG) );
75 if (tmpedge == NULL) maerror("edge in inittree");
76 *edge = tmpedge;
77
78 /* allocate memory for vetmpctor with edge numbers of leaves */
79 tmpedgeofleaf = (int *) calloc((size_t) Maxspc, sizeof(int) );
80 if (tmpedgeofleaf == NULL) maerror("edgeofleaf in inittree");
81 *edgeofleaf = tmpedgeofleaf;
82
83 /* allocate memory for all the edges the edge map */
84 for (i = 0; i < Maxbrnch; i++) {
85 (tmpedge)[i].edgemap = (int *) calloc((size_t) Maxbrnch, sizeof(int) );
86 if (tmpedge[i].edgemap == NULL) maerror("edgemap in inittree");
87 }
88
89 /* number all edges */
90 for (i = 0; i < Maxbrnch; i++) tmpedge[i].numedge = i;
91
92 /* initialize tree */
93
94 *nextedge = 3;
95 *nextleaf = 3;
96
97 /* edge maps */
98 (tmpedge[0].edgemap)[0] = 0; /* you are on the right edge */
99 (tmpedge[0].edgemap)[1] = 4; /* go down left for leaf 1 */
100 (tmpedge[0].edgemap)[2] = 5; /* go down right for leaf 2 */
101 (tmpedge[1].edgemap)[0] = 1; /* go up for leaf 0 */
102 (tmpedge[1].edgemap)[1] = 0; /* you are on the right edge */
103 (tmpedge[1].edgemap)[2] = 3; /* go up/down right for leaf 2 */
104 (tmpedge[2].edgemap)[0] = 1; /* go up for leaf 0 */
105 (tmpedge[2].edgemap)[1] = 2; /* go up/down left for leaf 1 */
106 (tmpedge[2].edgemap)[2] = 0; /* you are on the right edge */
107
108 /* interconnection */
109 tmpedge[0].up = NULL;
110 tmpedge[0].downleft = &tmpedge[1];
111 tmpedge[0].downright = &tmpedge[2];
112 tmpedge[1].up = &tmpedge[0];
113 tmpedge[1].downleft = NULL;
114 tmpedge[1].downright = NULL;
115 tmpedge[2].up = &tmpedge[0];
116 tmpedge[2].downleft = NULL;
117 tmpedge[2].downright = NULL;
118
119 /* edges of leaves */
120 tmpedgeofleaf[0] = 0;
121 tmpedgeofleaf[1] = 1;
122 tmpedgeofleaf[2] = 2;
123 } /* inittree */
124
125 /******************/
126 /* add next leaf on the specified edge */
addnextleaf_orig(int dockedge,ONEEDGE_ORIG * edge,int * edgeofleaf,int Maxspc,int Maxbrnch,int * in_nextedge,int * in_nextleaf)127 void addnextleaf_orig(int dockedge, /* insert here */
128 ONEEDGE_ORIG *edge, /* edge array */
129 int *edgeofleaf, /* ext. edge idx array */
130 int Maxspc, /* No. of species */
131 int Maxbrnch, /* No. of branches */
132 int *in_nextedge, /* next free edge idx */
133 int *in_nextleaf) /* next free leaf idx */
134 {
135 int i;
136 int nextedge;
137 int nextleaf;
138
139 nextedge=*in_nextedge;
140 nextleaf=*in_nextleaf;
141
142 if (dockedge >= nextedge) {
143 /* Trying to add leaf nextleaf to nonexisting edge dockedge */
144 fprintf(STDOUT, "\n\n\nHALT: PLEASE REPORT ERROR F TO DEVELOPERS\n\n\n");
145 # if PARALLEL
146 PP_Finalize();
147 # endif
148 exit(1);
149 }
150
151 if (nextleaf >= Maxspc) {
152 /* Trying to add leaf nextleaf to a tree with Maxspc leaves */
153 fprintf(STDOUT, "\n\n\nHALT: PLEASE REPORT ERROR G TO DEVELOPERS\n\n\n");
154 # if PARALLEL
155 PP_Finalize();
156 # endif
157 exit(1);
158 }
159
160 /* necessary change in edgeofleaf if dockedge == edgeofleaf[0] */
161 if (edgeofleaf[0] == dockedge) edgeofleaf[0] = nextedge;
162
163 /* adding nextedge to the tree */
164 edge[nextedge].up = edge[dockedge].up;
165 edge[nextedge].downleft = &edge[dockedge];
166 edge[nextedge].downright = &edge[nextedge+1];
167 edge[dockedge].up = &edge[nextedge];
168
169 if (edge[nextedge].up != NULL) {
170 if ( ((edge[nextedge].up)->downleft) == &edge[dockedge] )
171 (edge[nextedge].up)->downleft = &edge[nextedge];
172 else
173 (edge[nextedge].up)->downright = &edge[nextedge];
174 }
175
176 /* adding nextedge + 1 to the tree */
177 edge[nextedge+1].up = &edge[nextedge];
178 edge[nextedge+1].downleft = NULL;
179 edge[nextedge+1].downright = NULL;
180 edgeofleaf[nextleaf] = nextedge+1;
181
182 /* the two new edges get info about the old edges */
183 /* nextedge */
184 for (i = 0; i < nextedge; i++) {
185 switch ( (edge[dockedge].edgemap)[i] ) {
186
187 /* down right changes to down left */
188 case 5: (edge[nextedge].edgemap)[i] = 4;
189 break;
190
191 /* null changes to down left */
192 case 0: (edge[nextedge].edgemap)[i] = 4;
193 break;
194
195 default:(edge[nextedge].edgemap)[i] = (edge[dockedge].edgemap)[i];
196 break;
197 }
198 }
199
200 /* nextedge + 1 */
201 for (i = 0; i < nextedge; i++) {
202 switch ( (edge[dockedge].edgemap)[i] ) {
203
204 /* up/down left changes to up */
205 case 2: (edge[nextedge+1].edgemap)[i] = 1;
206 break;
207
208 /* up/down right changes to up */
209 case 3: (edge[nextedge+1].edgemap)[i] = 1;
210 break;
211
212 /* down left changes to up/down left */
213 case 4: (edge[nextedge+1].edgemap)[i] = 2;
214 break;
215
216 /* down right changes to up/down left */
217 case 5: (edge[nextedge+1].edgemap)[i] = 2;
218 break;
219
220 /* null changes to up/down left */
221 case 0: (edge[nextedge+1].edgemap)[i] = 2;
222 break;
223
224 /* up stays up */
225 default:(edge[nextedge+1].edgemap)[i] =
226 (edge[dockedge].edgemap)[i];
227 break;
228 }
229 }
230
231 /* dockedge */
232 for (i = 0; i < nextedge; i++) {
233 switch ( (edge[dockedge].edgemap)[i] ) {
234
235 /* up/down right changes to up */
236 case 3: (edge[dockedge].edgemap)[i] = 1;
237 break;
238
239 /* up/down left changes to up */
240 case 2: (edge[dockedge].edgemap)[i] = 1;
241 break;
242
243 default: break;
244 }
245 }
246
247 /* all edgemaps are updated for the two new edges */
248 /* nextedge */
249 (edge[nextedge].edgemap)[nextedge] = 0;
250 (edge[nextedge].edgemap)[nextedge+1] = 5; /* down right */
251
252 /* nextedge + 1 */
253 (edge[nextedge+1].edgemap)[nextedge] = 1; /* up */
254 (edge[nextedge+1].edgemap)[nextedge+1] = 0;
255
256 /* all other edges */
257 for (i = 0; i < nextedge; i++) {
258 (edge[i].edgemap)[nextedge] = (edge[i].edgemap)[dockedge];
259 (edge[i].edgemap)[nextedge+1] = (edge[i].edgemap)[dockedge];
260 }
261
262 /* an extra for dockedge */
263 (edge[dockedge].edgemap)[nextedge] = 1; /* up */
264 (edge[dockedge].edgemap)[nextedge+1] = 3; /* up/down right */
265
266 nextleaf++;
267 nextedge = nextedge + 2;
268
269 *in_nextedge=nextedge;
270 *in_nextleaf=nextleaf;
271 } /* addnextleaf */
272
273
274 /******************/
275 /* free memory (to be called after inittree) */
freetree_orig(ONEEDGE_ORIG * edge,int * edgeofleaf,int Maxspc)276 void freetree_orig(ONEEDGE_ORIG *edge, /* edge array */
277 int *edgeofleaf, /* ext. edge idx array */
278 int Maxspc) /* No. of species */
279 {
280 int i;
281
282 for (i = 0; i < 2 * Maxspc - 3; i++) free(edge[i].edgemap);
283 free(edge);
284 free(edgeofleaf);
285 } /* freetree */
286
287
288 /******************/
289 /* writes OTU sitting on edge ed */
writeOTU_orig(FILE * outfp,int ed,ONEEDGE_ORIG * edge,int * edgeofleaf,int nextedge,int nextleaf,int * column,int * trueID)290 void writeOTU_orig(FILE *outfp, /* output file */
291 int ed, /* edge to subtree */
292 ONEEDGE_ORIG *edge, /* edge array */
293 int *edgeofleaf, /* ext. edge idx array */
294 int nextedge, /* next free edge idx */
295 int nextleaf, /* next free leaf idx */
296 int *column, /* current screen depth */
297 int *trueID) /* species permutation */
298 {
299 int i;
300
301 /* test whether we are on a leaf */
302 if (edge[ed].downright == NULL && edge[ed].downleft == NULL) {
303 for (i = 1; i < nextleaf; i++) {
304 if (edgeofleaf[i] == ed) { /* i is the leaf of ed */
305 *column += fputid(outfp, trueID[i]);
306 return;
307 }
308 }
309 }
310
311 /* we are NOT on a leaf */
312 fprintf(outfp, "(");
313 (*column)++;
314 writeOTU_orig(outfp, edge[ed].downleft->numedge, edge, edgeofleaf, nextedge, nextleaf, column, trueID);
315 fprintf(outfp, ",");
316 (*column)++;
317 (*column)++;
318 if (*column > 55) {
319 *column = 2;
320 fprintf(outfp, "\n ");
321 }
322 writeOTU_orig(outfp, edge[ed].downright->numedge, edge, edgeofleaf, nextedge, nextleaf, column, trueID);
323 fprintf(outfp, ")");
324 (*column)++;
325 } /* writeOTU */
326
327 /******************/
328 /* write tree */
writetree_orig(FILE * outfp,ONEEDGE_ORIG * edge,int * edgeofleaf,int nextedge,int nextleaf,int * trueID)329 void writetree_orig(FILE *outfp, /* output file */
330 ONEEDGE_ORIG *edge, /* edge array */
331 int *edgeofleaf, /* ext. edge idx array */
332 int nextedge, /* next free edge idx */
333 int nextleaf, /* next free leaf idx */
334 int *trueID) /* species permutation */
335 {
336 int column;
337
338 column = 1;
339 fprintf(outfp, "(");
340 column += fputid(outfp, trueID[0]) + 3;
341 fprintf(outfp, ",");
342 writeOTU_orig(outfp, edge[edgeofleaf[0]].downleft->numedge, edge, edgeofleaf, nextedge, nextleaf, &column, trueID);
343 column++;
344 column++;
345 fprintf(outfp, ",");
346 writeOTU_orig(outfp, edge[edgeofleaf[0]].downright->numedge, edge, edgeofleaf, nextedge, nextleaf, &column, trueID);
347 fprintf(outfp, ");\n");
348 } /* writetree */
349
350
351 /******************/
352 /* clear all edgeinfos */
resetedgeinfo_orig(ONEEDGE_ORIG * edge,int nextedge)353 void resetedgeinfo_orig(ONEEDGE_ORIG *edge, /* edge array */
354 int nextedge) /* next free edge idx */
355 {
356 int i;
357
358 for (i = 0; i < nextedge; i++)
359 edge[i].edgeinfo = 0;
360 } /* resetedgeinfo */
361
362 /******************/
363 /* increment all edgeinfo between leaf A and B */
incrementedgeinfo_orig(int A,int B,ONEEDGE_ORIG * edge,int * edgeofleaf)364 void incrementedgeinfo_orig(int A, /* start leaf of penalty path */
365 int B, /* start leaf of penalty path */
366 ONEEDGE_ORIG *edge, /* edge array */
367 int *edgeofleaf) /* ext. edge idx array */
368 {
369 int curredge, finaledge, nextstep;
370
371 if (A == B) return;
372
373 finaledge = edgeofleaf[B];
374
375 curredge = edgeofleaf[A];
376 edge[curredge].edgeinfo = edge[curredge].edgeinfo + 1;
377
378 while (curredge != finaledge) {
379 nextstep = (edge[curredge].edgemap)[finaledge];
380 switch (nextstep) {
381
382 /* up */
383 case 1: curredge = (edge[curredge].up)->numedge;
384 break;
385
386 /* up/down left */
387 case 2: curredge = ((edge[curredge].up)->downleft)->numedge;
388 break;
389
390 /* up/down right */
391 case 3: curredge = ((edge[curredge].up)->downright)->numedge;
392 break;
393
394 /* down left */
395 case 4: curredge = (edge[curredge].downleft)->numedge;
396 break;
397
398 /* down right */
399 case 5: curredge = (edge[curredge].downright)->numedge;
400 break;
401
402 }
403 edge[curredge].edgeinfo = edge[curredge].edgeinfo + 1;
404 }
405 } /* incrementedgeinfo */
406
407
408 /******************/
409 /* checks which edge has the lowest edgeinfo
410 if there are several edges with the same lowest edgeinfo,
411 one of them will be selected randomly */
minimumedgeinfo_orig(ONEEDGE_ORIG * edge,int * edgeofleaf,int nextedge,int nextleaf,int * out_minedge,uli * out_mininfo)412 void minimumedgeinfo_orig(ONEEDGE_ORIG *edge, /* edge array */
413 int *edgeofleaf, /* ext. edge idx array */
414 int nextedge, /* next free edge idx */
415 int nextleaf, /* next free leaf idx */
416 int *out_minedge, /* minimum edge set */
417 uli *out_mininfo) /* minumum penalty */
418 {
419 int i, k, howmany, randomnum;
420 int minedge; /* minimum edge set */
421 uli mininfo; /* minumum penalty */
422
423 howmany = 1;
424 minedge = 0;
425 mininfo = edge[0].edgeinfo;
426 for (i = 1; i < nextedge; i++)
427 if (edge[i].edgeinfo <= mininfo) {
428 if (edge[i].edgeinfo == mininfo) {
429 howmany++;
430 } else {
431 minedge = i;
432 mininfo = edge[i].edgeinfo;
433 howmany = 1;
434 }
435 }
436
437 if (howmany > 1) { /* draw random edge */
438 randomnum = randominteger(howmany) + 1; /* 1 to howmany */
439 i = -1;
440 for (k = 0; k < randomnum; k++) {
441 do {
442 i++;
443 } while (edge[i].edgeinfo != mininfo);
444 minedge = i;
445 }
446 }
447
448 *out_minedge=minedge;
449 *out_mininfo=mininfo;
450 } /* minimumedgeinfo */
451
452
453
454
455
456 /*************************************************************************/
457 /* global functions of the puzzling step */
458 /*************************************************************************/
459
460 /* perform one single puzzling step to produce one intermediate tree */
onepstep(ONEEDGE_ORIG ** edge,int ** edgeofleaf,unsigned char * quartettopols,int Maxspc,ivector permutation)461 void onepstep( /* PStep (intermediate) tree topol: */
462 ONEEDGE_ORIG **edge, /* out: new array of edges */
463 int **edgeofleaf, /* out: array of extern edge ptrs */
464 unsigned char *quartettopols, /* in: quartetblock with all topols */
465 int Maxspc, /* in: Number of species (n) */
466 ivector permutation) /* in: species permutation (trueID) */
467 {
468 int Maxbrnch = (2*Maxspc)-3; /* Number of branches (2n-3) */
469 int nextedge; /* next free edge index (=3) */
470 int nextleaf; /* next free leaf index (=3) */
471 int minedge; /* insertion edge of minimum penalty */
472 uli mininfo; /* minimum penalty */
473 int a, b, c, i; /* quartet leaves, i to be inserted */
474 int chooseX, chooseY; /* end leaves of penalty path */
475
476 /* allocate and initialize new tree topology */
477 inittree_orig(edge, edgeofleaf, Maxspc, Maxbrnch, &nextedge, &nextleaf);
478
479 /* adding all other leafs */
480 for (i = 3; i < Maxspc; i++) {
481
482 /* clear all edgeinfos */
483 resetedgeinfo_orig(*edge, nextedge);
484
485 /*
486 * core of quartet puzzling algorithm
487 */
488
489 for (a = 0; a < nextleaf - 2; a++)
490 for (b = a + 1; b < nextleaf - 1; b++)
491 for (c = b + 1; c < nextleaf; c++) {
492
493 /* check which two leaves out of
494 a,b,c are closer related to each
495 other than to leaf i according to a
496 least squares fit of the continous
497 Bayesian weights to the seven
498 trivial "attractive regions". We
499 assign a score of 1 to all edges
500 on the penalty path between these
501 two leaves chooseX and chooseY */
502
503 checkquartet(a, b, c, i, permutation, &chooseX, &chooseY);
504 incrementedgeinfo_orig(chooseX, chooseY, *edge, *edgeofleaf);
505
506 } /* for all quartets q=(a,b,c,i) */
507
508 /* find out which edge has lowest edgeinfo */
509 minimumedgeinfo_orig(*edge, *edgeofleaf, nextedge, nextleaf, &minedge, &mininfo);
510
511 /* add the next leaf on minedge */
512 addnextleaf_orig(minedge, *edge, *edgeofleaf, Maxspc, Maxbrnch, &nextedge, &nextleaf);
513 } /* adding all other leafs */
514 } /* onepstep */
515
516
517
518 /*************************************************************************/
519
520 /* perform Numtrial single puzzling steps constructing Numtrial intermediate */
521 /* trees, sort each of them and extract splits for consensus step */
allpstep(uli Numtrial,unsigned char * quartetinfo,int Maxspc,int fixedorder_optn)522 void allpstep(uli Numtrial, /* in: no. psteps to perform */
523 unsigned char *quartetinfo, /* in: quartetblock with all topols*/
524 int Maxspc, /* in: Number of species (n) */
525 int fixedorder_optn) /* in: 'fixed' anchored RNG (debug)*/
526 {
527 /* misc variables */
528 ivector trueIDtmp; /* species permutation (trueID) */
529 uli Currtrial; /* step counter */
530
531 /* PStep (intermediate) tree topol: */
532 ONEEDGE_ORIG *tmpedgeset; /* array of edges */
533 int *tmpedgeofleaf; /* array of extern edge pointers */
534
535 /* for unique sorting of tree topologies */
536 int *ctree; /* sort tree */
537 int startnode; /* root in sort tree */
538 char *trstr; /* phylip tree string */
539 treelistitemtype *treeitem; /* list entry of tree */
540
541 # if PARALLEL
542 cmatrix *bp; /* temporary buffer for splits in slaves */
543 int n; /* step count for init/free split buffer */
544 # endif /* PARALLEL */
545
546
547 # if PARALLEL
548 /* add times to arry */
549 addtimes(GENERAL, &tarr);
550
551 /* alloc temporary split memory */
552 bp = (cmatrix *) calloc((size_t) Numtrial, sizeof(void *));
553 for(n=0; n<Numtrial; n++) {
554 bp[n] = new_cmatrix(Maxspc - 3, Maxspc);
555 }
556 # endif /* PARALLEL */
557
558
559 /* do single puzzling steps Numtrial times */
560 trueIDtmp = new_ivector(Maxspc);
561 for (Currtrial = 0; Currtrial < Numtrial; Currtrial++) {
562 /* randomize input order */
563 /* todo: chooser for missing data (HAS) */
564 if (fixedorder_optn) {
565 fixedchooser(Maxspc, Maxspc, trueIDtmp);
566 } else {
567 chooser(Maxspc, Maxspc, trueIDtmp);
568 }
569
570 /* one single puzzling step */
571 onepstep(&tmpedgeset, &tmpedgeofleaf, quartetinfo, Maxspc, trueIDtmp);
572
573 /* compute bipartitions of current tree */
574 /* missing parameters (HAS) */
575 # if ! PARALLEL
576 computebiparts(tmpedgeset, tmpedgeofleaf, trueIDtmp,
577 biparts, outgroup, Maxspc);
578
579 makenewsplitentries(biparts, Maxspc-3,
580 &splitpatterns, &splitsizes,
581 &splitfreqs, &numbiparts,
582 &maxbiparts, Maxspc);
583 # else /* PARALLEL */
584 computebiparts(tmpedgeset, tmpedgeofleaf, trueIDtmp,
585 bp[Currtrial], outgroup, Maxspc);
586
587 /* file splits not (yet) needed in slaves */
588 /* filing is done centrally at the master */
589 /* maybe optimizable by presorting (HAS) */
590 # if 0
591 makenewsplitentries(bp[Currtrial], Maxspc-3,
592 &splitpatterns, &splitsizes,
593 &splitfreqs, &numbiparts,
594 &maxbiparts, Maxspc);
595 # endif
596 # endif /* else PARALLEL */
597
598
599 /* Sorting to find unique topologies: */
600 /****************************************************/
601 /* sort tree topology branches, so that within */
602 /* each pair of sister subtrees the right-most */
603 /* leaf's ID of the right subtree is smaller than */
604 /* the right-most leaf's ID of the lest subtree. */
605 /* Exception: The outgroup (root) is the tree's */
606 /* right-most leaf! */
607
608 /* init/alloc copy tree structure */
609 ctree = initctree();
610 /* copy tree structure */
611 copytree(ctree, trueIDtmp, tmpedgeset, tmpedgeofleaf, Maxspc);
612 /* sort tree, startnode=outgroup */
613 startnode = sortctree(ctree);
614 /* sorted tree -> phylip tree string */
615 trstr=sprintfctree(ctree, psteptreestrlen);
616
617 /* add sorted tree to unique tree list */
618 treeitem = addtree2list(&trstr, 1, &psteptreelist, &psteptreenum, &psteptreesum);
619
620 # if ! PARALLEL
621 /* output unique tree to trace list, if option set */
622 /* not done on slave processes */
623 if((listqptrees == PSTOUT_LIST) || (listqptrees == PSTOUT_LISTORDER)) {
624 /* print: order no/# topol per this id/tree id/sum of unique topologies/sum of trees so far */
625 fprintf(qptlist, "%ld.\t1\t%d\t%d\t%d\t%d\n",
626 Currtrial + 1, (*treeitem).count, (*treeitem).id, psteptreenum, psteptreesum);
627 }
628 # endif /* ! PARALLEL */
629
630 # ifdef VERBOSE1
631 printf("%s\n", trstr);
632 printfsortedpstrees(psteptreelist);
633 # endif
634
635 /* free sorting tree structure */
636 freectree(&ctree);
637
638 /* free tree before building the next tree */
639 freetree_orig(tmpedgeset, tmpedgeofleaf, Maxspc);
640
641 # if ! PARALLEL
642 /* generate message roughly every 15 minutes */
643 /* check timer */
644 time(&time2);
645 if ( (time2 - time1) > TIMECHECK_INTERVAL) {
646 double tc2, mintogo, minutes, hours;
647 /* every TIMECHECK_INTERVAL seconds */
648 /* percentage of completed trees */
649 if (mflag == 0) {
650 fprintf(STDOUT, "\n");
651 mflag = 1;
652 }
653 tc2 = 100.0*Currtrial/Numtrial;
654 /* excluded for less runtime, but less accuracy of the 15min */
655 /* + 100.0*nq/Numquartets/Numtrial; */
656 mintogo = (100.0-tc2) *
657 (double) (time2-time0)/60.0/tc2;
658 hours = floor(mintogo/60.0);
659 minutes = mintogo - 60.0*hours;
660 fprintf(STDOUT, "%2.2f%%", tc2);
661 fprintf(STDOUT, " completed (remaining");
662 fprintf(STDOUT, " time: %.0f", hours);
663 fprintf(STDOUT, " hours %.0f", minutes);
664 fprintf(STDOUT, " minutes)\n");
665 fflush(STDOUT);
666 time1 = time2;
667 } /* check timer */
668 # endif /* ! PARALLEL */
669
670 addtimes(PUZZLING, &tarr);
671 } /* for Numtrials */
672
673 # if PARALLEL
674 /* in slaves: send results: splits, trees */
675 PP_SendSplitsBlock(Maxspc, Numtrial, bp, psteptreenum, psteptreelist);
676
677 /* free temporary split memory */
678 for (Currtrial = 0; Currtrial < Numtrial; Currtrial++) {
679 free_cmatrix(bp[Currtrial]);
680 }
681 free(bp);
682
683 /* sent! Thus, in slave process not needed any more */
684 freetreelist(&psteptreelist, &psteptreenum, &psteptreesum);
685 # endif /* PARALLEL */
686
687 } /* allpstep */
688
689
690