1 #if defined(HAVE_CONFIG_H) && !defined(DISABLE_DAKOTA_CONFIG_H)
2 #include "ddace_config.h"
3 #endif
4 /**
5
6 These programs construct and manipulate orthogonal
7 arrays. They were prepared by
8
9 Art Owen
10 Department of Statistics
11 Sequoia Hall
12 Stanford CA 94305
13
14 They may be freely used and shared. This code comes
15 with no warranty of any kind. Use it at your own
16 risk.
17
18 I thank the Semiconductor Research Corporation and
19 the National Science Foundation for supporting this
20 work.
21
22 */
23
24
25 #include <stdio.h>
26 /* for exit: */
27 #include <stdlib.h>
28
29 #include "oa.h"
30
31 int **imatrix(), *ivector();
32
33 /* BMA, 20151007: Prototypes added for functions in this file to avoid
34 implicit declaration warnings, rather than reordering. */
35 void OA_fput(FILE *stream, int **A, int n, int k, int q);
36 int OA_fget(FILE *stream, int **A, int n, int k, int q, int eof_assert);
37 int OA_fread(FILE *stream, int ***A, int *n, int *k, int*q);
38 int OA_str0(int q, int nrow, int ncol, int **A, int verbose);
39 int OA_str1(int q, int nrow, int ncol, int **A, int verbose);
40 int OA_strt(int q, int nrow, int ncol, int **A, int t, int verbose);
41
42 /* BMA, 20151007: Prototypes added for functions in other files to
43 avoid implicit declaration warnings. Keeping out of headers to
44 avoid conflicts. */
45 void free_ivector(int *v, int nl, int nh);
46 int grow_imatrix_byrows(int ***imat, int oldrowsize, int newrowsize, int colsize);
47 int ipow(int a, int b);
48
49
50 /**
51
52 Glossary:
53
54 OA_put write OA to standard output
55 OA_fput write OA to stream
56 OA_get get OA from standard input
57 OA_fget get OA from stream
58
59 OA_parsein read arguments q,nrow,ncol to OA "filter programs"
60 OA_strworkcheck warn about large work loads in strength checking programs
61
62 */
63
64
65 /* OUTPUT OUTPUT OUTPUT OUTPUT OUTPUT OUTPUT OUTPUT */
66
OA_put(A,n,k,q)67 void OA_put( A,n,k,q )
68 int **A, n, k, q;
69 {
70 OA_fput( stdout,A,n,k,q );
71 }
72
OA_fput(stream,A,n,k,q)73 void OA_fput( stream,A,n,k,q )
74 FILE *stream;
75 int **A, n, k, q;
76 {
77 int i,j;
78 char* format;
79
80 format = "%d%s";
81 if( q < 1000 )format = "%3d%s";
82 if( q < 100 )format = "%2d%s";
83 if( q < 10 )format = "%1d%s";
84 if( q < 0 )format = "%d%s";
85
86 for( i=0; i<n; i++ )
87 for( j=0; j<k; j++ )
88 fprintf(stream,format,A[i][j],(j==k-1)?"\n":" ");
89 }
90
91
92 /* INPUT INPUT INPUT INPUT INPUT INPUT INPUT INPUT */
93
OA_get(A,n,k,q,eof_assert)94 int OA_get( A,n,k,q, eof_assert )
95 int **A, n, k, q, eof_assert;
96 {
97 return OA_fget( stdin,A,n,k,q,eof_assert );
98 }
99
OA_fget(stream,A,n,k,q,eof_assert)100 int OA_fget( stream,A,n,k,q,eof_assert )
101 FILE *stream;
102 int **A, n, k, q, eof_assert;
103 {
104 int i,j;
105 for( i=0; i<n; i++ )
106 for( j=0; j<k; j++ )
107 if( fscanf(stream,"%d",&A[i][j]) == EOF ){
108 fprintf(stderr,"Unexpected end of input encountered. Wanted to read\n");
109 fprintf(stderr,"%d rows of %d cols. Failed trying for row %d, col %d.\n",
110 n,k,i,j);
111 return 0;
112 }else if( A[i][j] >= q ){
113 fprintf(stderr,"Invalid array element %d. All elements should be\n",
114 A[i][j]);
115 fprintf(stderr,"strictly less than q = %d.\n",q);
116 return 0;
117 }else if( A[i][j] <0 ){
118 fprintf(stderr,"Invalid array element %d, should be >= 0.\n",A[i][j]);
119 return 0;
120 }
121
122 if( eof_assert && fscanf(stream,"%d",&eof_assert) != EOF ){
123 fprintf(stderr,"Input has more integers than expected.\n");
124 fprintf(stderr,"Perhaps the number of rows and/or columns is incorrect.\n");
125 return 0;
126 }
127 return 1;
128 }
129
130 /* READ READ READ READ READ READ READ READ READ */
131
132 /* Read array and determine dimensions. Inspired by approach
133 taken in Xgobi. */
134
135 #define MAXK 5000
136 #define ROWINC 1000
137 int line0[ MAXK ];
138
OA_read(A,n,k,q)139 int OA_read( A,n,k,q )
140 int ***A, *n, *k, *q;
141 {
142 return OA_fread( stdin,A,n,k,q );
143 }
144
145
OA_fread(stream,A,n,k,q)146 int OA_fread( stream,A,n,k,q )
147 FILE *stream;
148 int ***A, *n, *k, *q;
149 {
150 int i, j;
151 char c;
152
153 *k = 0;
154
155 while( (c = getc(stream)) != '\n' ){
156 if( c=='\t' || c==' ' )
157 ;
158 else{
159 ungetc(c,stream);
160 if( *k >= MAXK ){
161 fprintf(stderr,"Error: Input appears to have more than %d columns.\n",*k);
162 fprintf(stderr,"If such large input is desired, increase MAXK in oa.c\n");
163 fprintf(stderr,"and re-install the software.\n");
164 return 0;
165 }
166 if( fscanf(stream,"%d",&line0[(*k)++]) <= 0 ){
167 fprintf(stderr,"Error: no newline character found. Could be empty\n");
168 fprintf(stderr,"input or matrix all on one line.\n");
169 return 0;
170 }
171 }
172 }
173
174 *A = imatrix( 0,ROWINC-1, 0,*k-1 );
175 if( !(*A) ){
176 fprintf(stderr,"Unable to allocate memory to read the array.\n");
177 return 0;
178 }
179
180 *n = 0;
181 for( j=0; j<*k; j++ )
182 A[0][*n][j] = line0[j];
183
184 while( 1 ){
185 (*n)++;
186 if( (*n % ROWINC)==0 )
187 if( !grow_imatrix_byrows( A, (*n), (*n)+ROWINC, *k ) ){
188 fprintf(stderr,"Unable to allocate extra memory for row %d of the array.\n"
189 ,*n);
190 return 0;
191 }
192
193 if( fscanf(stream,"%d",&(A[0][*n][0])) == EOF )
194 break;
195
196 for( j=1; j<*k; j++ )
197 if( fscanf(stream,"%d",&(A[0][*n][j])) == EOF ){
198 fprintf(stderr,"Unexpected end of input encountered. Row %d only has\n",*n);
199 fprintf(stderr,"%d elements of %d expected.\n",j,*k);
200 return 0;
201 }
202 }
203 *q = 0; /* Assume that q = max(A)+1 */
204
205 for( i=0; i<*n; i++ )
206 for( j=0; j<*k; j++ )
207 if( A[0][i][j] > *q )
208 *q = A[0][i][j];
209 *q = *q +1;
210
211 return 1;
212 }
213
214 /* PARSE PARSE PARSE PARSE PARSE PARSE PARSE PARSE */
215
OA_parsein(argc,argv,q,nrow,ncol,A)216 void OA_parsein( argc,argv, q,nrow,ncol, A )
217 int argc;
218 char *argv[];
219 int *q, *nrow, *ncol, ***A;
220 {
221 if( argc<=1 ){
222 if( !OA_read( A,nrow,ncol,q ) ){
223 fprintf(stderr,"Fatal error while reading the array.\n");
224 exit(1);
225 }
226 }
227 else if( argc==2 ){
228 sscanf(argv[1],"%d",q);
229 scanf("%d %d",nrow,ncol);
230 }else if( argc==3 ){
231 sscanf(argv[1],"%d",q);
232 sscanf(argv[2],"%d",nrow);
233 scanf("%d",ncol);
234 }else{
235 sscanf(argv[1],"%d",q);
236 sscanf(argv[2],"%d",nrow);
237 sscanf(argv[3],"%d",ncol);
238 }
239
240 if( *q<1 ){
241 fprintf(stderr,"Array has only %d symbol(s). At least one\n",*q);
242 fprintf(stderr,"symbol is necessary in an orthogonal array.\n");
243 exit(1);
244 }
245
246 if( *ncol<1 ){
247 fprintf(stderr,"Array has only %d column(s). At least one\n",*ncol);
248 fprintf(stderr,"column is necessary in an orthogonal array.\n");
249 exit(1);
250 }
251
252 if( *nrow<1 ){
253 fprintf(stderr,"Array has only %d rows. At least one\n",*nrow);
254 fprintf(stderr,"row is necessary in an orthogonal array.\n");
255 exit(1);
256 }
257
258
259 if( argc >1 ){
260 *A = imatrix( 0,*nrow-1,0,*ncol-1 );
261 if( !(*A) ){
262 fprintf(stderr,"The array is too large (%d by %d) to fit in memory.\n",nrow,ncol);
263 exit(1);
264 }
265
266 if( !OA_get( *A,*nrow,*ncol,*q, 1 ) ){/* Read 'em all */
267 fprintf(stderr,"Read error getting the orthogonal array.\n");
268 exit(1);
269 }
270 }
271 }
272
273
274 /* WORK WORK WORK WORK WORK WORK WORK WORK WORK */
275
OA_strworkcheck(work,str)276 void OA_strworkcheck( work,str )
277 double work;
278 int str;
279 {
280 if( work > BIGWORK ){
281 fprintf(stderr,"If the array has strength %d, %g comparisons will\n",
282 str, work);
283 fprintf(stderr,"be required to prove it. This might take a long time.\n");
284 fprintf(stderr,"This warning is triggered when more than %d comparisons\n",
285 BIGWORK);
286 fprintf(stderr,"are required. To avoid this warning increase BIGWORK in\n");
287 fprintf(stderr,"oa.h. Intermediate results will be printed.\n\n");
288 fflush(stderr);
289 }else if( work > MEDWORK ){
290 fprintf(stderr,"Since more than %d comparisons may be required to\n",MEDWORK);
291 fprintf(stderr,"to check whether the array has strength %d, intermediate\n",
292 str);
293 fprintf(stderr,"results will be printed. To avoid this warning increase\n");
294 fprintf(stderr,"MEDWORK in oa.h\n\n");
295 fflush(stderr);
296 }
297 }
298
299
300 /* STRENGTH STRENGTH STRENGTH STRENGTH STRENGTH */
301
302
OA_strength(q,nrow,ncol,A,str,verbose)303 void OA_strength( q,nrow,ncol,A,str,verbose )
304 int q,nrow,ncol,**A,*str, verbose;
305 /*
306 Calculate and return the strength of the array A.
307
308 verbose = 0 => No printed output
309 verbose = 1 => Only stderr output
310 verbose = 2 => Output to both stdout and stderr
311
312 */
313 {
314 *str = -1;
315
316 if( OA_str0( q,nrow,ncol,A,verbose) )
317 *str = 0;
318 else
319 return;
320 if( OA_str1( q,nrow,ncol,A,verbose) )
321 *str = 1;
322 else
323 return;
324 while( OA_strt( q,nrow,ncol,A,*str+1,verbose ) )
325 (*str)++;
326 return;
327 }
328
329
330 /* Check strength 0 */
OA_str0(q,nrow,ncol,A,verbose)331 int OA_str0( q,nrow,ncol,A,verbose )
332 int q,nrow,ncol,**A, verbose;
333 {
334 int i, j1;
335
336 for( j1=0; j1<ncol; j1++ )
337 for( i=0; i<nrow; i++ )
338 if( A[i][j1] < 0 || A[i][j1] >= q ){
339 if( verbose >= 2 ){
340 printf("Array is not even of strength 0, that is there are elements\n");
341 // BMA 20141204; added q-1 as argument to silence compiler warning
342 //printf("other than integers 0 through %d inclusive in it.\n");
343 printf("other than integers 0 through %d inclusive in it.\n", (q-1));
344 printf("The first exception is A[%d][%d] = %d.\n",i,j1,A[i][j1]);
345 }
346 return 0;
347 }
348 if( verbose >=2 )
349 printf("The array has strength (at least) 0.\n");
350 return 1;
351 }
352
353
354 /* Check strength 1 */
OA_str1(q,nrow,ncol,A,verbose)355 int OA_str1( q,nrow,ncol,A,verbose )
356 int q,nrow,ncol,**A, verbose;
357 {
358 int i, j1, q1;
359 int lambda, count;
360 double work;
361
362 if( nrow%q ){
363 if( verbose >= 2 ){
364 printf("The array cannot have strength 1, because the number\n");
365 printf("of rows %d is not a multiple of q = %d.\n",nrow,q);
366 }
367 return 0;
368 }
369
370 lambda = nrow/q;
371 work = nrow * ncol * q * 1.0;
372 OA_strworkcheck( work,1 );
373 for( j1=0; j1<ncol; j1++ ){
374 for( q1=0; q1<q; q1++ ){
375 count = 0;
376 for( i=0; i<nrow; i++ )
377 count += (A[i][j1]==q1);
378 if( count != lambda ){
379 if( verbose >= 2 ){
380 printf("Array is not of strength 1. The first violation arises for\n");
381 printf("the number of times A[,%d] = %d.\n",
382 j1, q1);
383 printf("This happened in %d rows, it should have happened in %d rows.\n",
384 count, lambda);
385 }
386 return 0;
387 }
388 }
389 if( work > MEDWORK && verbose > 0 )
390 fprintf(stderr,"No violation of strength 1 involves column %d.\n",j1);
391 }
392 if( verbose >=2 )
393 printf("The array has strength (at least) 1.\n");
394 return 1;
395 }
396
397 /* Check strength 2 */
OA_str2(q,nrow,ncol,A,verbose)398 int OA_str2( q,nrow,ncol,A,verbose )
399 int q,nrow,ncol,**A, verbose;
400 {
401 int i, j1,j2, q1,q2;
402 int lambda, count;
403 double work;
404
405 if( ncol<2 ){
406 if( verbose > 0 ){
407 fprintf(stderr,"Array has only %d column(s). At least two\n",ncol);
408 fprintf(stderr,"columns are necessary for strength 2 to make sense.\n");
409 }
410 return 0;
411 }
412 if( nrow % (q*q) ){
413 if( verbose > 0 ){
414 fprintf(stderr,"The array cannot have strength 2, because the number\n");
415 fprintf(stderr,"of rows %d is not a multiple of q^2 = %d^2 = %d.\n",nrow,q,q*q);
416 }
417 return 0;
418 }
419
420 lambda = nrow/(q*q);
421 work = nrow*ncol*(ncol-1.0)*q*q/2.0;
422 OA_strworkcheck( work,2 );
423
424 for( j1=0; j1<ncol; j1++ ){
425 for( j2=j1+1; j2<ncol; j2++ ){
426 for( q1=0; q1<q; q1++ )
427 for( q2=0; q2<q; q2++ ){
428 count = 0;
429 for( i=0; i<nrow; i++ )
430 count += (A[i][j1]==q1)&&(A[i][j2]==q2);
431 if( count != lambda ){
432 if( verbose >= 2 ){
433 printf("Array is not of strength 2. The first violation arises for\n");
434 printf("the number of times (A[,%d],A[,%d]) = (%d,%d).\n",
435 j1,j2, q1,q2 );
436 printf("This happened in %d rows, it should have happened in %d rows.\n",
437 count, lambda);
438 }
439 return 0;
440 }
441 }
442 }
443 if( work > MEDWORK && verbose > 0 )
444 fprintf(stderr,"No violation of strength 2 involves column %d.\n",j1);
445 }
446
447 if( verbose >=2 )
448 printf("The array has strength (at least) 2.\n");
449 return 1;
450 }
451
452
453
454 /* Check strength 3 */
OA_str3(q,nrow,ncol,A,verbose)455 int OA_str3( q,nrow,ncol,A,verbose )
456 int q,nrow,ncol,**A, verbose;
457 {
458 int i, j1,j2,j3, q1,q2,q3;
459 int lambda, count;
460 double work;
461
462 if( ncol<3 ){
463 if( verbose > 0 ){
464 fprintf(stderr,"Array has only %d column(s). At least three\n",ncol);
465 fprintf(stderr,"columns are necessary for strength 3 to make sense.\n");
466 }
467 return 0;
468 }
469 if( nrow % (q*q*q) ){
470 if( verbose > 0 ){
471 fprintf(stderr,"The array cannot have strength 3, because the number\n");
472 fprintf(stderr,"of rows %d is not a multiple of q^3 = %d^3 = %d.\n",nrow,q,q*q*q);
473 }
474 return 0;
475 }
476
477 lambda = nrow/(q*q*q);
478 work = nrow*ncol*(ncol-1.0)*(ncol-2.0)*q*q*q/6.0;
479 OA_strworkcheck( work,3 );
480
481 for( j1=0; j1<ncol; j1++ ){
482 for( j2=j1+1; j2<ncol; j2++ )
483 for( j3=j2+1; j3<ncol; j3++ ){
484 for( q1=0; q1<q; q1++ )
485 for( q2=0; q2<q; q2++ )
486 for( q3=0; q3<q; q3++ ){
487 count = 0;
488 for( i=0; i<nrow; i++ )
489 count += (A[i][j1]==q1)&&(A[i][j2]==q2)&&(A[i][j3]==q3);
490 if( count != lambda ){
491 if( verbose >= 2 ){
492 printf("Array is not of strength 3. The first violation arises for\n");
493 printf("the number of times (A[,%d],A[,%d],A[,%d]) = (%d,%d,%d).\n",
494 j1,j2,j3, q1,q2,q3 );
495 printf("This happened in %d rows, it should have happened in %d rows.\n",
496 count, lambda);
497 }
498 return 0;
499 }
500 }
501 }
502 if( work > MEDWORK && verbose > 0 )
503 fprintf(stderr,"No violation of strength 3 involves column %d.\n",j1);
504 }
505 if( verbose >=2 )
506 printf("The array has strength (at least) 3.\n");
507 return 1;
508 }
509
510
511 /* Check strength 4 */
OA_str4(q,nrow,ncol,A,verbose)512 int OA_str4( q,nrow,ncol,A,verbose )
513 int q,nrow,ncol,**A, verbose;
514 {
515 int i, j1,j2,j3,j4, q1,q2,q3,q4;
516 int lambda, count;
517 double work;
518
519 if( ncol<4 ){
520 if( verbose > 0 ){
521 fprintf(stderr,"Array has only %d column(s). At least four\n",ncol);
522 fprintf(stderr,"columns are necessary for strength 4 to make sense.\n");
523 }
524 return 0;
525 }
526 if( nrow % (q*q*q*q) ){
527 if( verbose > 0 ){
528 fprintf(stderr,"The array cannot have strength 4, because the number\n");
529 fprintf(stderr,"of rows %d is not a multiple of q^4 = %d^4 = %d.\n",nrow,q,q*q*q*q);
530 }
531 return 0;
532 }
533
534 lambda = nrow/(q*q*q*q);
535 work = nrow*ncol*(ncol-1.0)*(ncol-2.0)*(ncol-3.0)*q*q*q*q/24.0;
536 OA_strworkcheck( work,4 );
537
538 for( j1=0; j1<ncol; j1++ ){
539 for( j2=j1+1; j2<ncol; j2++ )
540 for( j3=j2+1; j3<ncol; j3++ )
541 for( j4=j3+1; j4<ncol; j4++ ){
542 for( q1=0; q1<q; q1++ )
543 for( q2=0; q2<q; q2++ )
544 for( q3=0; q3<q; q3++ )
545 for( q4=0; q4<q; q4++ ){
546 count = 0;
547 for( i=0; i<nrow; i++ )
548 count += (A[i][j1]==q1)&&(A[i][j2]==q2)&&(A[i][j3]==q3)&&(A[i][j4]==q4);
549 if( count != lambda ){
550 if( verbose >= 2 ){
551 printf("Array is not of strength 4. The first violation arises for\n");
552 printf("the number of times (A[,%d],A[,%d],A[,%d],A[,%d]) = (%d,%d,%d,%d).\n",
553 j1,j2,j3,j4, q1,q2,q3,q4 );
554 printf("This happened in %d rows, it should have happened in %d rows.\n",
555 count, lambda);
556 }
557 return 0;
558 }
559 }
560 }
561 }
562 if( work > MEDWORK && verbose > 0 )
563 fprintf(stderr,"No violation of strength 4 involves column %d.\n",j1);
564
565 if( verbose >=2 )
566 printf("The array has strength (at least) 4.\n");
567 return 1;
568 }
569
570
571 /* Check strength t */
OA_strt(q,nrow,ncol,A,t,verbose)572 int OA_strt( q,nrow,ncol,A,t,verbose )
573 int q,nrow,ncol,**A,t,verbose;
574 {
575 int row, i, ic, iq, *clist, *qlist, ctuples, qtuples;
576 int lambda, count, match;
577 double work;
578
579 if( t<0 ){
580 if( verbose > 0 ){
581 fprintf(stderr,"Don't know how to verify strength %d. It doesn't\n",t);
582 fprintf(stderr,"make sense.\n");
583 }
584 return 0;
585 }
586 if( ncol<t ){
587 if( verbose > 0 ){
588 fprintf(stderr,"Array has only %d column(s). At least %d\n",ncol,t);
589 fprintf(stderr,"columns are necessary for strength %d to make sense.\n",t);
590 }
591 return 0;
592 }
593 if( t==0 )
594 return OA_str0( q,nrow,ncol,A,verbose );
595 if( nrow % ipow(q,t) ){
596 if( verbose > 0 ){
597 fprintf(stderr,"The array cannot have strength %d, because the number\n",t);
598 fprintf(stderr,"of rows %d is not a multiple of q^%d = %d^%d = %d.\n",
599 nrow,t,q,t,ipow(q,t));
600 }
601 return 0;
602 }
603
604 lambda = nrow/ipow(q,t);
605 work = nrow*ipow(q,t);
606 ctuples = 1;
607
608 clist = ivector( 0,t-1 );
609 qlist = ivector( 0,t-1 );
610
611 for( i=0; i<t; i++ ){
612 work *= (ncol-i)/(i+1.0);
613 ctuples *= ncol-i;
614 qlist[i] = 0;
615 clist[i] = i;
616 }
617 for( i=0; i<t; i++ )
618 ctuples /= (i+1);
619 qtuples = ipow(q,t);
620
621 OA_strworkcheck( work,t );
622
623 for( ic=0; ic<ctuples; ic++ ){ /* Loop over ordered tuples of columns */
624 /* for( i=0; i<t; i++ )
625 printf("%s %d%s",(i==0)?"Col":"",clist[i],(i==t-1)?"\n":" ");*/
626
627 for( iq=0; iq<qtuples; iq++ ){ /* Loop over unordered tuples of symbols */
628 /* for( i=0; i<t; i++ )
629 printf(" %s %d%s",(i==0)?"Sym":"",qlist[i],(i==t-1)?"\n":" ");*/
630 count = 0;
631 for( row=0; row<nrow; row++ ){
632 match = 1;
633 for( i=0; i<t && match; i++ )
634 match *= A[row][clist[i]] == qlist[i];
635 count += match;
636 }
637 if( count != lambda ){
638 if( verbose >= 2 ){
639 printf("Array is not of strength %d. The first violation arises for\n",t);
640 printf("the number of times (");
641 for( i=0; i<t; i++ )
642 printf("A[,%d]%s",clist[i],(i==t-1)?")":",");
643 printf(" = (");
644 for( i=0; i<t; i++ )
645 printf("%d%s",qlist[i],(i==t-1)?").\n":",");
646 printf("This happened in %d rows, it should have happened in %d rows.\n",
647 count, lambda);
648 }
649 return 0;
650 }
651 for( i=t-1; i>=0; i-- ){
652 qlist[i] = (qlist[i]+1) % q;
653 if( qlist[i] )break;
654 }
655 }
656
657 for( i= t-1; i>=0; i-- ){
658 clist[i] = (clist[i]+1) % (ncol+i-t+1);
659 if( clist[i] )break;
660 }
661
662 if( work > MEDWORK && verbose > 0 && (t==1||t>1 && clist[1]==0) )
663 fprintf(stderr,"No violation of strength %d involves column %d.\n",
664 t,(clist[0]+ncol-1)%ncol);
665
666 for( i=1; i< t; i++ )
667 if( clist[i] <= clist[i-1] )
668 clist[i] = clist[i-1]+1;
669 }
670
671 if( verbose >=2 )
672 printf("The array has strength (at least) %d.\n",t);
673
674 /* Added by J Cramp on 07.02.04 to
675 * fix memory leak problems. */
676 free_ivector(clist,0,t-1);
677 free_ivector(qlist,0,t-1);
678
679 return 1;
680 }
681