springRenumber.C
Go to the documentation of this file.
1/*---------------------------------------------------------------------------*\
2 ========= |
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4 \\ / O peration |
5 \\ / A nd | www.openfoam.com
6 \\/ M anipulation |
7-------------------------------------------------------------------------------
8 Copyright (C) 2011-2017 OpenFOAM Foundation
9 Copyright (C) 2019-2022 OpenCFD Ltd.
10-------------------------------------------------------------------------------
11License
12 This file is part of OpenFOAM.
13
14 OpenFOAM is free software: you can redistribute it and/or modify it
15 under the terms of the GNU General Public License as published by
16 the Free Software Foundation, either version 3 of the License, or
17 (at your option) any later version.
18
19 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22 for more details.
23
24 You should have received a copy of the GNU General Public License
25 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26
27\*---------------------------------------------------------------------------*/
28
29#include "springRenumber.H"
31#include "decompositionMethod.H"
32
33// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34
35namespace Foam
36{
38
40 (
44 );
45}
46
47
48// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
49
51:
53 coeffsDict_(dict.optionalSubDict(typeName+"Coeffs")),
54 maxIter_(coeffsDict_.get<label>("maxIter")),
55 maxCo_(coeffsDict_.get<scalar>("maxCo")),
56 freezeFraction_(coeffsDict_.get<scalar>("freezeFraction"))
57{}
58
59
60// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
61
62template<class ConnectionListListType>
63Foam::labelList Foam::springRenumber::renumberImpl
64(
65 const ConnectionListListType& cellCells
66) const
67{
68 const label nOldCells(cellCells.size());
69
70 // Look at cell index as a 1D position parameter.
71 // Move cells to the average 'position' of their neighbour.
72
73 scalarField position(nOldCells);
74 forAll(position, celli)
75 {
76 position[celli] = celli;
77 }
78
79 // Sum force per cell. Also reused for the displacement
80 scalarField sumForce(nOldCells);
81
82 labelList oldToNew(identity(nOldCells));
83
84 scalar maxCo = (maxCo_ * nOldCells);
85
86 for (label iter = 0; iter < maxIter_; ++iter)
87 {
88 //Pout<< "Iteration : " << iter << nl
89 // << "------------"
90 // << endl;
91
92 //Pout<< "Position :" << nl
93 // << " min : " << min(position) << nl
94 // << " max : " << max(position) << nl
95 // << " avg : " << average(position) << nl
96 // << endl;
97
98 // Sum force per cell.
99 sumForce = Zero;
100 for (label oldCelli = 0; oldCelli < nOldCells; ++oldCelli)
101 {
102 const label celli = oldToNew[oldCelli];
103 const auto& neighbours = cellCells[oldCelli];
104
105 for (const label nbr : neighbours)
106 {
107 const label nbrCelli = oldToNew[nbr];
108
109 sumForce[celli] += (position[nbrCelli]-position[celli]);
110 }
111 }
112
113 //Pout<< "Force :" << nl
114 // << " min : " << min(sumForce) << nl
115 // << " max : " << max(sumForce) << nl
116 // << " avgMag : " << average(mag(sumForce)) << nl
117 // << "DeltaT : " << deltaT << nl
118 // << endl;
119
120 // Limit displacement
121 scalar deltaT = maxCo / max(mag(sumForce));
122
123 Info<< "Iter:" << iter
124 << " maxCo:" << maxCo
125 << " deltaT:" << deltaT
126 << " average force:" << average(mag(sumForce)) << endl;
127
128 // Determine displacement
129 scalarField& displacement = sumForce;
130 displacement *= deltaT;
131
132 //Pout<< "Displacement :" << nl
133 // << " min : " << min(displacement) << nl
134 // << " max : " << max(displacement) << nl
135 // << " avgMag : " << average(mag(displacement)) << nl
136 // << endl;
137
138 // Calculate new position and scale to be within original range
139 // (0..nCells-1) for ease of postprocessing.
140 position += displacement;
141 position -= min(position);
142 position *= (position.size()-1)/max(position);
143
144 // Slowly freeze.
145 maxCo *= freezeFraction_;
146 }
147
148 //writeOBJ("endPosition.obj", cellCells, position);
149
150 // Move cells to new position
151 labelList shuffle(sortedOrder(position));
152
153 // Reorder oldToNew
154 inplaceReorder(shuffle, oldToNew);
155
156 return invert(nOldCells, oldToNew);
157}
158
159
161(
162 const polyMesh& mesh,
163 const pointField&
164) const
165{
166 CompactListList<label> cellCells;
167 decompositionMethod::calcCellCells
168 (
169 mesh,
171 mesh.nCells(),
172 false, // local only
173 cellCells
174 );
175
176 return renumberImpl(cellCells);
177}
178
179
181(
182 const CompactListList<label>& cellCells,
183 const pointField&
184) const
185{
186 return renumberImpl(cellCells);
187}
188
189
191(
192 const labelListList& cellCells,
193 const pointField&
194) const
195{
196 return renumberImpl(cellCells);
197}
198
199
200// ************************************************************************* //
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
A packed storage unstructured matrix of objects of type <T> using an offset table for access.
void size(const label n)
Older name for setAddressableSize.
Definition: UList.H:114
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:81
label nCells() const noexcept
Number of mesh cells.
Abstract base class for renumbering.
Use spring analogy - attract neighbouring cells according to the distance of their cell indices.
virtual labelList renumber(const pointField &) const
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
dynamicFvMesh & mesh
scalar maxCo
Namespace for OpenFOAM.
void shuffle(UList< T > &list)
Randomise the list order.
Definition: UList.C:370
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
labelList identity(const label len, label start=0)
Return an identity map of the given length with (map[i] == i)
Definition: labelList.C:38
messageStream Info
Information stream (stdout output on master, null elsewhere)
dimensioned< Type > average(const DimensionedField< Type, GeoMesh > &df)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
void inplaceReorder(const labelUList &oldToNew, ListType &input, const bool prune=false)
Inplace reorder the elements of a list.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
labelList sortedOrder(const UList< T > &input)
Return the (stable) sort order for the list.
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
labelList invert(const label len, const labelUList &map)
Create an inverse one-to-one mapping.
Definition: ListOps.C:36
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:333