/[escript]/trunk/dudley/src/Util.cpp
ViewVC logotype

Contents of /trunk/dudley/src/Util.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4789 - (show annotations)
Fri Mar 21 00:30:32 2014 UTC (5 years, 3 months ago) by sshaw
File size: 18773 byte(s)
changing more awkward logic into more obvious (and correct, though i doubt a negative length ever happens) logic
1
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2014 by University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Open Software License version 3.0
9 * http://www.opensource.org/licenses/osl-3.0.php
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development 2012-2013 by School of Earth Sciences
13 * Development from 2014 by Centre for Geoscience Computing (GeoComp)
14 *
15 *****************************************************************************/
16
17 /************************************************************************************/
18
19 /* Some utility routines: */
20
21 /************************************************************************************/
22 #include "esysUtils/maths.h"
23 #include "Util.h"
24
25 #ifdef _OPENMP
26 #include <omp.h>
27 #endif
28
29 #include "esysUtils/index.h"
30 #include "esysUtils/mem.h"
31 #include <limits.h>
32 #include "string.h" /* for memcpy*/
33
34 /************************************************************************************/
35
36 /* returns true if any of the values in the short array values is not equal to zero */
37
38 bool Dudley_Util_anyNonZeroDouble(dim_t N, double *values)
39 {
40 dim_t q;
41 for (q = 0; q < N; ++q)
42 if (ABS(values[q]) > 0)
43 return TRUE;
44 return FALSE;
45 }
46
47 /************************************************************************************/
48
49 /* gathers double values out from in by index: */
50
51 /* out(1:numData,1:len)=in(1:numData,index(1:len)) */
52
53 void Dudley_Util_Gather_double(dim_t len, index_t * index, dim_t numData, double *in, double *out)
54 {
55 dim_t s, i;
56 for (s = 0; s < len; s++)
57 {
58 for (i = 0; i < numData; i++)
59 {
60 out[INDEX2(i, s, numData)] = in[INDEX2(i, index[s], numData)];
61 }
62 }
63 }
64
65 /************************************************************************************/
66
67 /* gathers maybelong values out from in by index: */
68
69 /* out(1:numData,1:len)=in(1:numData,index(1:len)) */
70
71 void Dudley_Util_Gather_int(dim_t len, index_t * index, dim_t numData, index_t * in, index_t * out)
72 {
73 dim_t s, i;
74 for (s = 0; s < len; s++)
75 {
76 for (i = 0; i < numData; i++)
77 {
78 out[INDEX2(i, s, numData)] = in[INDEX2(i, index[s], numData)];
79 }
80 }
81 }
82
83 /************************************************************************************/
84
85 /* adds a vector in into out using and index. */
86
87 /* out(1:numData,index[p])+=in(1:numData,p) where p = {k=1...len , index[k]<upperBound}*/
88
89 void Dudley_Util_AddScatter(const dim_t len, const index_t * index, const dim_t numData, const double *in, double *out, const index_t upperBound)
90 {
91 dim_t i, s;
92 for (s = 0; s < len; s++)
93 {
94 for (i = 0; i < numData; i++)
95 {
96 if (index[s] < upperBound)
97 {
98 out[INDEX2(i, index[s], numData)] += in[INDEX2(i, s, numData)];
99 }
100 }
101 }
102 }
103
104 /* multiplies two matrices */
105
106 /* A(1:A1,1:A2)=B(1:A1,1:B2)*C(1:B2,1:A2) */
107
108 void Dudley_Util_SmallMatMult(dim_t A1, dim_t A2, double *A, dim_t B2, const double *B, const double *C)
109 {
110 dim_t i, j, s;
111 register double rtmp;
112 for (i = 0; i < A1; i++)
113 {
114 for (j = 0; j < A2; j++)
115 {
116 rtmp = 0;
117 for (s = 0; s < B2; s++)
118 {
119 rtmp += B[INDEX2(i, s, A1)] * C[INDEX2(s, j, B2)];
120 }
121 A[INDEX2(i, j, A1)] = rtmp;
122 }
123 }
124 }
125
126 /* multiplies two sets of matrices: */
127
128 /* A(1:A1,1:A2,i)=B(1:A1,1:B2,i)*C(1:B2,1:A2,i) i=1,len */
129
130 void Dudley_Util_SmallMatSetMult(dim_t len, dim_t A1, dim_t A2, double *A, dim_t B2, const double *B, const double *C)
131 {
132 dim_t q, i, j, s;
133 register double rtmp;
134 for (q = 0; q < len; q++)
135 {
136 for (i = 0; i < A1; i++)
137 {
138 for (j = 0; j < A2; j++)
139 {
140 rtmp = 0;
141 for (s = 0; s < B2; s++)
142 rtmp += B[INDEX3(i, s, q, A1, B2)] * C[INDEX3(s, j, q, B2, A2)];
143 A[INDEX3(i, j, q, A1, A2)] = rtmp;
144 }
145 }
146 }
147 }
148
149 /* multiplies a set of matrices with a single matrix: */
150
151 /* A(1:A1,1:A2,i)=B(1:A1,1:B2,i)*C(1:B2,1:A2) i=1,len */
152
153 void Dudley_Util_SmallMatSetMult1(dim_t len, dim_t A1, dim_t A2, double *A, dim_t B2, const double *B, const double *C)
154 {
155 dim_t q, i, j, s;
156 register double rtmp;
157 for (q = 0; q < len; q++)
158 {
159 for (i = 0; i < A1; i++)
160 {
161 for (j = 0; j < A2; j++)
162 {
163 rtmp = 0;
164 for (s = 0; s < B2; s++)
165 rtmp += B[INDEX3(i, s, q, A1, B2)] * C[INDEX2(s, j, B2)];
166 A[INDEX3(i, j, q, A1, A2)] = rtmp;
167 }
168 }
169 }
170 }
171
172 /* inverts the set of dim x dim matrices A(:,:,1:len) with dim=1,2,3 */
173 /* the determinant is returned. */
174
175 void Dudley_Util_InvertSmallMat(dim_t len, dim_t dim, double *A, double *invA, double *det)
176 {
177 dim_t q;
178 register double D, A11, A12, A13, A21, A22, A23, A31, A32, A33;
179
180 switch (dim)
181 {
182 case 1:
183 for (q = 0; q < len; q++)
184 {
185 D = A[q];
186 if (ABS(D) > 0)
187 {
188 det[q] = D;
189 D = 1. / D;
190 invA[q] = D;
191 }
192 else
193 {
194 Dudley_setError(ZERO_DIVISION_ERROR, __FILE__ ": Non-regular matrix");
195 return;
196 }
197 }
198 break;
199
200 case 2:
201 for (q = 0; q < len; q++)
202 {
203 A11 = A[INDEX3(0, 0, q, 2, 2)];
204 A12 = A[INDEX3(0, 1, q, 2, 2)];
205 A21 = A[INDEX3(1, 0, q, 2, 2)];
206 A22 = A[INDEX3(1, 1, q, 2, 2)];
207
208 D = A11 * A22 - A12 * A21;
209 if (ABS(D) > 0)
210 {
211 det[q] = D;
212 D = 1. / D;
213 invA[INDEX3(0, 0, q, 2, 2)] = A22 * D;
214 invA[INDEX3(1, 0, q, 2, 2)] = -A21 * D;
215 invA[INDEX3(0, 1, q, 2, 2)] = -A12 * D;
216 invA[INDEX3(1, 1, q, 2, 2)] = A11 * D;
217 }
218 else
219 {
220 Dudley_setError(ZERO_DIVISION_ERROR, __FILE__ ": Non-regular matrix");
221 return;
222 }
223 }
224 break;
225
226 case 3:
227 for (q = 0; q < len; q++)
228 {
229 A11 = A[INDEX3(0, 0, q, 3, 3)];
230 A21 = A[INDEX3(1, 0, q, 3, 3)];
231 A31 = A[INDEX3(2, 0, q, 3, 3)];
232 A12 = A[INDEX3(0, 1, q, 3, 3)];
233 A22 = A[INDEX3(1, 1, q, 3, 3)];
234 A32 = A[INDEX3(2, 1, q, 3, 3)];
235 A13 = A[INDEX3(0, 2, q, 3, 3)];
236 A23 = A[INDEX3(1, 2, q, 3, 3)];
237 A33 = A[INDEX3(2, 2, q, 3, 3)];
238
239 D = A11 * (A22 * A33 - A23 * A32) + A12 * (A31 * A23 - A21 * A33) + A13 * (A21 * A32 - A31 * A22);
240 if (ABS(D) > 0)
241 {
242 det[q] = D;
243 D = 1. / D;
244 invA[INDEX3(0, 0, q, 3, 3)] = (A22 * A33 - A23 * A32) * D;
245 invA[INDEX3(1, 0, q, 3, 3)] = (A31 * A23 - A21 * A33) * D;
246 invA[INDEX3(2, 0, q, 3, 3)] = (A21 * A32 - A31 * A22) * D;
247 invA[INDEX3(0, 1, q, 3, 3)] = (A13 * A32 - A12 * A33) * D;
248 invA[INDEX3(1, 1, q, 3, 3)] = (A11 * A33 - A31 * A13) * D;
249 invA[INDEX3(2, 1, q, 3, 3)] = (A12 * A31 - A11 * A32) * D;
250 invA[INDEX3(0, 2, q, 3, 3)] = (A12 * A23 - A13 * A22) * D;
251 invA[INDEX3(1, 2, q, 3, 3)] = (A13 * A21 - A11 * A23) * D;
252 invA[INDEX3(2, 2, q, 3, 3)] = (A11 * A22 - A12 * A21) * D;
253 }
254 else
255 {
256 Dudley_setError(ZERO_DIVISION_ERROR, __FILE__ ": Non-regular matrix");
257 return;
258 }
259 }
260 break;
261
262 }
263 return;
264 }
265
266 /* sets the derterminat of a set of dim x dim matrices A(:,:,1:len) with dim=1,2,3 */
267
268 void Dudley_Util_DetOfSmallMat(dim_t len, dim_t dim, double *A, double *det)
269 {
270 dim_t q;
271 register double A11, A12, A13, A21, A22, A23, A31, A32, A33;
272
273 switch (dim)
274 {
275 case 1:
276 for (q = 0; q < len; q++)
277 {
278 det[q] = A[q];
279 }
280 break;
281
282 case 2:
283 for (q = 0; q < len; q++)
284 {
285 A11 = A[INDEX3(0, 0, q, 2, 2)];
286 A12 = A[INDEX3(0, 1, q, 2, 2)];
287 A21 = A[INDEX3(1, 0, q, 2, 2)];
288 A22 = A[INDEX3(1, 1, q, 2, 2)];
289
290 det[q] = A11 * A22 - A12 * A21;
291 }
292 break;
293
294 case 3:
295 for (q = 0; q < len; q++)
296 {
297 A11 = A[INDEX3(0, 0, q, 3, 3)];
298 A21 = A[INDEX3(1, 0, q, 3, 3)];
299 A31 = A[INDEX3(2, 0, q, 3, 3)];
300 A12 = A[INDEX3(0, 1, q, 3, 3)];
301 A22 = A[INDEX3(1, 1, q, 3, 3)];
302 A32 = A[INDEX3(2, 1, q, 3, 3)];
303 A13 = A[INDEX3(0, 2, q, 3, 3)];
304 A23 = A[INDEX3(1, 2, q, 3, 3)];
305 A33 = A[INDEX3(2, 2, q, 3, 3)];
306
307 det[q] = A11 * (A22 * A33 - A23 * A32) + A12 * (A31 * A23 - A21 * A33) + A13 * (A21 * A32 - A31 * A22);
308 }
309 break;
310
311 }
312 return;
313 }
314
315 /* returns the normalized vector Normal[dim,len] orthogonal to A(:,0,q) and A(:,1,q) in the case of dim=3 */
316 /* or the vector A(:,0,q) in the case of dim=2 */
317
318 void Dudley_NormalVector(dim_t len, dim_t dim, dim_t dim1, double *A, double *Normal)
319 {
320 dim_t q;
321 register double A11, A12, CO_A13, A21, A22, CO_A23, A31, A32, CO_A33, length, invlength;
322
323 switch (dim)
324 {
325 case 1:
326 for (q = 0; q < len; q++)
327 Normal[q] = 1;
328 break;
329 case 2:
330 for (q = 0; q < len; q++)
331 {
332 A11 = A[INDEX3(0, 0, q, 2, dim1)];
333 A21 = A[INDEX3(1, 0, q, 2, dim1)];
334 length = sqrt(A11 * A11 + A21 * A21);
335 if (length <= 0)
336 {
337 Dudley_setError(ZERO_DIVISION_ERROR, __FILE__ ": area equals zero.");
338 return;
339 }
340 else
341 {
342 invlength = 1. / length;
343 Normal[INDEX2(0, q, 2)] = A21 * invlength;
344 Normal[INDEX2(1, q, 2)] = -A11 * invlength;
345 }
346 }
347 break;
348 case 3:
349 for (q = 0; q < len; q++)
350 {
351 A11 = A[INDEX3(0, 0, q, 3, dim1)];
352 A21 = A[INDEX3(1, 0, q, 3, dim1)];
353 A31 = A[INDEX3(2, 0, q, 3, dim1)];
354 A12 = A[INDEX3(0, 1, q, 3, dim1)];
355 A22 = A[INDEX3(1, 1, q, 3, dim1)];
356 A32 = A[INDEX3(2, 1, q, 3, dim1)];
357 CO_A13 = A21 * A32 - A31 * A22;
358 CO_A23 = A31 * A12 - A11 * A32;
359 CO_A33 = A11 * A22 - A21 * A12;
360 length = sqrt(CO_A13 * CO_A13 + CO_A23 * CO_A23 + CO_A33 * CO_A33);
361 if (length <= 0)
362 {
363 Dudley_setError(ZERO_DIVISION_ERROR, __FILE__ ": area equals zero.");
364 return;
365 }
366 else
367 {
368 invlength = 1. / length;
369 Normal[INDEX2(0, q, 3)] = CO_A13 * invlength;
370 Normal[INDEX2(1, q, 3)] = CO_A23 * invlength;
371 Normal[INDEX2(2, q, 3)] = CO_A33 * invlength;
372 }
373
374 }
375 break;
376
377 }
378 return;
379 }
380
381 /* return the length of the vector which is orthogonal to the vectors A(:,0,q) and A(:,1,q) in the case of dim=3 */
382 /* or the vector A(:,0,q) in the case of dim=2 */
383
384 void Dudley_LengthOfNormalVector(dim_t len, dim_t dim, dim_t dim1, double *A, double *length)
385 {
386 dim_t q;
387 double A11, A12, CO_A13, A21, A22, CO_A23, A31, A32, CO_A33;
388
389 switch (dim)
390 {
391 case 1:
392 for (q = 0; q < len; q++)
393 length[q] = 1;
394 break;
395 case 2:
396 for (q = 0; q < len; q++)
397 {
398 A11 = A[INDEX3(0, 0, q, 2, dim1)];
399 A21 = A[INDEX3(1, 0, q, 2, dim1)];
400 length[q] = sqrt(A11 * A11 + A21 * A21);
401 }
402 break;
403 case 3:
404 for (q = 0; q < len; q++)
405 {
406 A11 = A[INDEX3(0, 0, q, 3, dim1)];
407 A21 = A[INDEX3(1, 0, q, 3, dim1)];
408 A31 = A[INDEX3(2, 0, q, 3, dim1)];
409 A12 = A[INDEX3(0, 1, q, 3, dim1)];
410 A22 = A[INDEX3(1, 1, q, 3, dim1)];
411 A32 = A[INDEX3(2, 1, q, 3, dim1)];
412 CO_A13 = A21 * A32 - A31 * A22;
413 CO_A23 = A31 * A12 - A11 * A32;
414 CO_A33 = A11 * A22 - A21 * A12;
415 length[q] = sqrt(CO_A13 * CO_A13 + CO_A23 * CO_A23 + CO_A33 * CO_A33);
416 }
417 break;
418
419 }
420 return;
421 }
422
423 /* inverts the map map of length len */
424 /* there is no range checking! */
425 /* at output Map[invMap[i]]=i for i=0:lenInvMap */
426
427 void Dudley_Util_InvertMap(dim_t lenInvMap, index_t * invMap, dim_t lenMap, index_t * Map)
428 {
429 dim_t i;
430 for (i = 0; i < lenInvMap; i++)
431 invMap[i] = 0;
432 for (i = 0; i < lenMap; i++)
433 {
434 if (Map[i] >= 0)
435 invMap[Map[i]] = i;
436 }
437 }
438
439 /* orders a Dudley_Util_ValueAndIndex array by value */
440 /* it is assumed that n is large */
441
442 int Dudley_Util_ValueAndIndex_compar(const void *arg1, const void *arg2)
443 {
444 Dudley_Util_ValueAndIndex *e1, *e2;
445 e1 = (Dudley_Util_ValueAndIndex *) arg1;
446 e2 = (Dudley_Util_ValueAndIndex *) arg2;
447 if (e1->value < e2->value)
448 return -1;
449 if (e1->value > e2->value)
450 return 1;
451 if (e1->index < e2->index)
452 return -1;
453 if (e1->index > e2->index)
454 return 1;
455 return 0;
456 }
457
458 void Dudley_Util_sortValueAndIndex(dim_t n, Dudley_Util_ValueAndIndex * array)
459 {
460 /* OMP : needs parallelization ! */
461 qsort(array, n, sizeof(Dudley_Util_ValueAndIndex), Dudley_Util_ValueAndIndex_compar);
462 }
463
464 /************************************************************************************/
465
466 /* calculates the minimum value from a dim X N integer array */
467
468 index_t Dudley_Util_getMinInt(dim_t dim, dim_t N, index_t * values)
469 {
470 dim_t i, j;
471 index_t out, out_local;
472 out = INDEX_T_MAX;
473 if (values != NULL && dim * N > 0)
474 {
475 out = values[0];
476 #pragma omp parallel private(out_local)
477 {
478 out_local = out;
479 #pragma omp for private(i,j) schedule(static)
480 for (j = 0; j < N; j++)
481 {
482 for (i = 0; i < dim; i++)
483 out_local = MIN(out_local, values[INDEX2(i, j, dim)]);
484 }
485 #pragma omp critical
486 out = MIN(out_local, out);
487 }
488 }
489 return out;
490 }
491
492 /* calculates the maximum value from a dim X N integer array */
493
494 index_t Dudley_Util_getMaxInt(dim_t dim, dim_t N, index_t * values)
495 {
496 dim_t i, j;
497 index_t out, out_local;
498 out = -INDEX_T_MAX;
499 if (values != NULL && dim * N > 0)
500 {
501 out = values[0];
502 #pragma omp parallel private(out_local)
503 {
504 out_local = out;
505 #pragma omp for private(i,j) schedule(static)
506 for (j = 0; j < N; j++)
507 {
508 for (i = 0; i < dim; i++)
509 {
510 out_local = MAX(out_local, values[INDEX2(i, j, dim)]);
511
512 }
513 }
514 #pragma omp critical
515 out = MAX(out_local, out);
516 }
517 }
518 return out;
519 }
520
521 /************************************************************************************/
522
523 /* calculates the minimum value from a dim X N integer array */
524
525 index_t Dudley_Util_getFlaggedMinInt(dim_t dim, dim_t N, index_t * values, index_t ignore)
526 {
527 dim_t i, j;
528 index_t out, out_local;
529 out = INDEX_T_MAX;
530 if (values != NULL && dim * N > 0)
531 {
532 out = values[0];
533 #pragma omp parallel private(out_local)
534 {
535 out_local = out;
536 #pragma omp for private(i,j) schedule(static)
537 for (j = 0; j < N; j++)
538 {
539 for (i = 0; i < dim; i++)
540 if (values[INDEX2(i, j, dim)] != ignore)
541 out_local = MIN(out_local, values[INDEX2(i, j, dim)]);
542 }
543 #pragma omp critical
544 out = MIN(out_local, out);
545 }
546 }
547 return out;
548 }
549
550 /* calculates the maximum value from a dim X N integer array */
551
552 index_t Dudley_Util_getFlaggedMaxInt(dim_t dim, dim_t N, index_t * values, index_t ignore)
553 {
554 dim_t i, j;
555 index_t out, out_local;
556 out = -INDEX_T_MAX;
557 if (values != NULL && dim * N > 0)
558 {
559 out = values[0];
560 #pragma omp parallel private(out_local)
561 {
562 out_local = out;
563 #pragma omp for private(i,j) schedule(static)
564 for (j = 0; j < N; j++)
565 {
566 for (i = 0; i < dim; i++)
567 if (values[INDEX2(i, j, dim)] != ignore)
568 out_local = MAX(out_local, values[INDEX2(i, j, dim)]);
569 }
570 #pragma omp critical
571 out = MAX(out_local, out);
572 }
573 }
574 return out;
575 }
576
577 /* set the index of the positive entries in mask. The length of index is returned. */
578
579 dim_t Dudley_Util_packMask(dim_t N, index_t * mask, index_t * index)
580 {
581 dim_t out, k;
582 out = 0;
583 /*OMP */
584 for (k = 0; k < N; k++)
585 {
586 if (mask[k] >= 0)
587 {
588 index[out] = k;
589 out++;
590 }
591 }
592 return out;
593 }
594
595 /* returns true if array contains value */
596 bool Dudley_Util_isAny(dim_t N, index_t * array, index_t value)
597 {
598 bool out = FALSE;
599 dim_t i;
600 #pragma omp parallel for private(i) schedule(static) reduction(||:out)
601 for (i = 0; i < N; i++)
602 out = out || (array[i] == value);
603 return out;
604 }
605
606 /* calculates the cummulative sum in array and returns the total sum */
607 index_t Dudley_Util_cumsum(dim_t N, index_t * array)
608 {
609 index_t out = 0, tmp;
610 dim_t i;
611 #ifdef _OPENMP
612 index_t *partial_sums = NULL, sum;
613 partial_sums = new index_t[omp_get_max_threads()];
614 #pragma omp parallel private(sum,i,tmp)
615 {
616 sum = 0;
617 #pragma omp for schedule(static)
618 for (i = 0; i < N; ++i)
619 sum += array[i];
620 partial_sums[omp_get_thread_num()] = sum;
621 #pragma omp barrier
622 #pragma omp master
623 {
624 out = 0;
625 for (i = 0; i < omp_get_max_threads(); ++i)
626 {
627 tmp = out;
628 out += partial_sums[i];
629 partial_sums[i] = tmp;
630 }
631 }
632 #pragma omp barrier
633 sum = partial_sums[omp_get_thread_num()];
634 #pragma omp for schedule(static)
635 for (i = 0; i < N; ++i)
636 {
637 tmp = sum;
638 sum += array[i];
639 array[i] = tmp;
640 }
641 }
642 delete[] partial_sums;
643 #else
644 for (i = 0; i < N; ++i)
645 {
646 tmp = out;
647 out += array[i];
648 array[i] = tmp;
649 }
650 #endif
651 return out;
652 }
653
654 void Dudley_Util_setValuesInUse(const index_t * values, const dim_t numValues, dim_t * numValuesInUse,
655 index_t ** valuesInUse, Esys_MPIInfo * mpiinfo)
656 {
657 dim_t i;
658 index_t lastFoundValue = INDEX_T_MIN, minFoundValue, local_minFoundValue, *newValuesInUse = NULL;
659 register index_t itmp;
660 bool allFound = FALSE;
661 dim_t nv = 0;
662
663 while (!allFound)
664 {
665 /*
666 * find smallest value bigger than lastFoundValue
667 */
668 minFoundValue = INDEX_T_MAX;
669 #pragma omp parallel private(local_minFoundValue)
670 {
671 local_minFoundValue = minFoundValue;
672 #pragma omp for private(i,itmp) schedule(static)
673 for (i = 0; i < numValues; i++)
674 {
675 itmp = values[i];
676 if ((itmp > lastFoundValue) && (itmp < local_minFoundValue))
677 local_minFoundValue = itmp;
678 }
679 #pragma omp critical
680 {
681 if (local_minFoundValue < minFoundValue)
682 minFoundValue = local_minFoundValue;
683 }
684
685 }
686 #ifdef ESYS_MPI
687 local_minFoundValue = minFoundValue;
688 MPI_Allreduce(&local_minFoundValue, &minFoundValue, 1, MPI_INT, MPI_MIN, mpiinfo->comm);
689 #endif
690 /* if we found a new tag we need to add this too the valuesInUseList */
691
692 if (minFoundValue < INDEX_T_MAX)
693 {
694 newValuesInUse = new index_t[nv + 1];
695 if (*valuesInUse != NULL)
696 {
697 memcpy(newValuesInUse, *valuesInUse, sizeof(index_t) * nv);
698 delete[] *valuesInUse;
699 }
700 newValuesInUse[nv] = minFoundValue;
701 *valuesInUse = newValuesInUse;
702 newValuesInUse = NULL;
703 nv++;
704 lastFoundValue = minFoundValue;
705 }
706 else
707 {
708 allFound = TRUE;
709 }
710 }
711 *numValuesInUse = nv;
712 }
713
714 #ifdef ESYS_MPI
715 void Dudley_printDoubleArray(FILE * fid, dim_t n, double *array, char *name)
716 {
717 index_t i;
718
719 if (name)
720 fprintf(fid, "%s [ ", name);
721 else
722 fprintf(fid, "[ ");
723 for (i = 0; i < (n < 60 ? n : 60); i++)
724 fprintf(fid, "%g ", array[i]);
725 if (n >= 30)
726 fprintf(fid, "... ");
727 fprintf(fid, "]\n");
728 }
729
730 void Dudley_printIntArray(FILE * fid, dim_t n, int *array, char *name)
731 {
732 index_t i;
733
734 if (name)
735 fprintf(fid, "%s [ ", name);
736 else
737 fprintf(fid, "[ ");
738 for (i = 0; i < (n < 60 ? n : 60); i++)
739 fprintf(fid, "%d ", array[i]);
740 if (n >= 30)
741 fprintf(fid, "... ");
742 fprintf(fid, "]\n");
743 }
744
745 void Dudley_printMaskArray(FILE * fid, dim_t n, int *array, char *name)
746 {
747 index_t i;
748
749 if (name)
750 fprintf(fid, "%s [ ", name);
751 else
752 fprintf(fid, "[ ");
753 for (i = 0; i < (n < 60 ? n : 60); i++)
754 if (array[i] != -1)
755 fprintf(fid, "%3d ", array[i]);
756 else
757 fprintf(fid, " * ");
758 if (n >= 30)
759 fprintf(fid, "... ");
760 fprintf(fid, "]\n");
761 }
762 #endif

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision
svn:mergeinfo /branches/lapack2681/finley/src/Util.cpp:2682-2741 /branches/pasowrap/dudley/src/Util.cpp:3661-3674 /branches/py3_attempt2/dudley/src/Util.cpp:3871-3891 /branches/restext/finley/src/Util.cpp:2610-2624 /branches/ripleygmg_from_3668/dudley/src/Util.cpp:3669-3791 /branches/stage3.0/finley/src/Util.cpp:2569-2590 /branches/symbolic_from_3470/dudley/src/Util.cpp:3471-3974 /branches/symbolic_from_3470/ripley/test/python/dudley/src/Util.cpp:3517-3974 /release/3.0/finley/src/Util.cpp:2591-2601 /trunk/dudley/src/Util.cpp:4257-4344 /trunk/ripley/test/python/dudley/src/Util.cpp:3480-3515

  ViewVC Help
Powered by ViewVC 1.1.26