58 void Foam::hierarchGeomDecomp::setOrder()
66 else if (order.size() != 3)
69 <<
"Number of characters in order (" << order <<
") != 3"
73 for (
int i = 0; i < 3; ++i)
79 case 'x': order_[i] = 0;
break;
80 case 'y': order_[i] = 1;
break;
81 case 'z': order_[i] = 2;
break;
85 <<
"Illegal decomposition order " << order <<
nl
86 <<
"It should only contain x, y or z"
94 Foam::label Foam::hierarchGeomDecomp::findLower
96 const UList<scalar>& list,
110 while ((high - low) > 1)
112 const label mid = (low + high)/2;
126 if (list[high-1] < val)
140 void Foam::hierarchGeomDecomp::calculateSortedWeightedSizes
145 const label globalCurrentSize,
151 sortedWeightedSizes[0] = 0;
154 const label pointi = current[indices[i]];
155 sortedWeightedSizes[i + 1] = sortedWeightedSizes[i] + weights[pointi];
160 sortedWeightedSizes[current.size()],
165 sortedWeightedSizes *= (globalCurrentSize/globalCurrentLength);
171 bool Foam::hierarchGeomDecomp::findBinary
174 const List<scalar>&
values,
175 const label minIndex,
178 const scalar wantedSize,
188 label low = minIndex;
189 label high =
values.size();
192 scalar midValuePrev = VGREAT;
196 label size =
returnReduce(mid-minIndex, sumOp<label>());
200 Pout<<
" low:" << low <<
" lowValue:" << lowValue
201 <<
" high:" << high <<
" highValue:" << highValue
202 <<
" mid:" << mid <<
" midValue:" << midValue <<
endl
203 <<
" globalSize:" << size <<
" wantedSize:" << wantedSize
204 <<
" sizeTol:" << sizeTol <<
endl;
207 if (wantedSize < size - sizeTol)
210 highValue = midValue;
212 else if (wantedSize > size + sizeTol)
223 midValue = 0.5*(lowValue+highValue);
227 bool hasNotChanged = (
mag(midValue-midValuePrev) < SMALL);
234 <<
"unable to find desired decomposition split, making do!"
241 midValuePrev = midValue;
250 bool Foam::hierarchGeomDecomp::findBinary
253 const List<scalar>& sortedWeightedSizes,
254 const List<scalar>&
values,
255 const label minIndex,
258 const scalar wantedSize,
265 label low = minIndex;
270 label high =
values.size();
273 scalar midValuePrev = VGREAT;
279 sortedWeightedSizes[mid] - sortedWeightedSizes[minIndex],
285 Pout<<
" low:" << low <<
" lowValue:" << lowValue
286 <<
" high:" << high <<
" highValue:" << highValue
287 <<
" mid:" << mid <<
" midValue:" << midValue <<
endl
288 <<
" globalSize:" << weightedSize
289 <<
" wantedSize:" << wantedSize
290 <<
" sizeTol:" << sizeTol <<
endl;
293 if (wantedSize < weightedSize - sizeTol)
296 highValue = midValue;
298 else if (wantedSize > weightedSize + sizeTol)
309 midValue = 0.5*(lowValue+highValue);
313 bool hasNotChanged = (
mag(midValue-midValuePrev) < SMALL);
320 <<
"Unable to find desired decomposition split, making do!"
327 midValuePrev = midValue;
335 Foam::label Foam::hierarchGeomDecomp::sortComponent
346 const label compI = order_[componentIndex];
353 Pout<<
"sortComponent : Sorting slice of size " << current.size()
354 <<
" in component " << compI <<
endl;
358 SortableList<scalar> sortedCoord(current.size());
362 const label pointi = current[i];
364 sortedCoord[i] =
points[pointi][compI];
368 label globalCurrentSize =
returnReduce(current.size(), sumOp<label>());
392 Pout<<
"sortComponent : minCoord:" << minCoord
393 <<
" maxCoord:" << maxCoord <<
endl;
399 scalar leftCoord = minCoord;
402 for (label bin = 0; bin < n_[compI]; bin++)
410 label localSize = -1;
413 scalar rightCoord = -GREAT;
415 if (bin == n_[compI]-1)
418 localSize = current.size()-leftIndex;
419 rightCoord = maxCoord;
424 localSize = label(current.size()/n_[compI]);
425 if (leftIndex+localSize < sortedCoord.size())
427 rightCoord = sortedCoord[leftIndex+localSize];
431 rightCoord = maxCoord;
440 label rightIndex = current.size();
441 rightCoord = maxCoord;
451 globalCurrentSize/n_[compI],
456 localSize = rightIndex - leftIndex;
466 Pout<<
"For component " << compI <<
", bin " << bin
467 <<
" copying" <<
endl
468 <<
"from " << leftCoord <<
" at local index "
470 <<
"to " << rightCoord <<
" localSize:"
481 label pointi = current[sortedCoord.indices()[leftIndex+i]];
484 finalDecomp[pointi] += bin*mult;
491 if (componentIndex < 2)
500 nWarnings += sortComponent
517 leftIndex += localSize;
518 leftCoord = rightCoord;
526 Foam::label Foam::hierarchGeomDecomp::sortComponent
538 const label compI = order_[componentIndex];
545 Pout<<
"sortComponent : Sorting slice of size " << current.size()
546 <<
" in component " << compI <<
endl;
550 SortableList<scalar> sortedCoord(current.size());
554 label pointi = current[i];
556 sortedCoord[i] =
points[pointi][compI];
560 label globalCurrentSize =
returnReduce(current.size(), sumOp<label>());
565 calculateSortedWeightedSizes
568 sortedCoord.indices(),
596 Pout<<
"sortComponent : minCoord:" << minCoord
597 <<
" maxCoord:" << maxCoord <<
endl;
603 scalar leftCoord = minCoord;
606 for (label bin = 0; bin < n_[compI]; bin++)
614 label localSize = -1;
617 scalar rightCoord = -GREAT;
619 if (bin == n_[compI]-1)
622 localSize = current.size()-leftIndex;
623 rightCoord = maxCoord;
632 label rightIndex = current.size();
633 rightCoord = maxCoord;
644 globalCurrentSize/n_[compI],
649 localSize = rightIndex - leftIndex;
659 Pout<<
"For component " << compI <<
", bin " << bin
660 <<
" copying" <<
endl
661 <<
"from " << leftCoord <<
" at local index "
663 <<
"to " << rightCoord <<
" localSize:"
674 label pointi = current[sortedCoord.indices()[leftIndex+i]];
677 finalDecomp[pointi] += bin*mult;
684 if (componentIndex < 2)
693 nWarnings += sortComponent
711 leftIndex += localSize;
712 leftCoord = rightCoord;
721 Foam::hierarchGeomDecomp::hierarchGeomDecomp
733 Foam::hierarchGeomDecomp::hierarchGeomDecomp
765 const label sizeTol =
max(1, label(1
e-3*allSize/nDomains_));
768 const label nWarnings = sortComponent
781 <<
"\nEncountered " << nWarnings <<
" occurrences where the desired"
782 " decomposition split could not be properly satisfied" <<
endl;
807 const label sizeTol =
max(1, label(1
e-3*allSize/nDomains_));
810 const label nWarnings = sortComponent
824 <<
"\nEncountered " << nWarnings <<
" occurrences where the desired"
825 " decomposition split could not be properly satisfied" <<
endl;