mirrorFvMesh.C 12.2 KB
Newer Older
1
2
3
4
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
5
    \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
6
7
8
9
10
     \\/     M anipulation  |
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

11
12
13
14
    OpenFOAM is free software: you can redistribute it and/or modify it
    under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.
15
16
17
18
19
20
21

    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
    for more details.

    You should have received a copy of the GNU General Public License
22
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41

\*---------------------------------------------------------------------------*/

#include "mirrorFvMesh.H"
#include "Time.H"
#include "plane.H"

// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //

Foam::mirrorFvMesh::mirrorFvMesh(const IOobject& io)
:
    fvMesh(io),
    mirrorMeshDict_
    (
        IOobject
        (
            "mirrorMeshDict",
            time().system(),
            *this,
42
            IOobject::MUST_READ_IF_MODIFIED,
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
            IOobject::NO_WRITE
        )
    ),
    mirrorMeshPtr_(NULL)
{
    plane mirrorPlane(mirrorMeshDict_);

    scalar planeTolerance
    (
        readScalar(mirrorMeshDict_.lookup("planeTolerance"))
    );

    const pointField& oldPoints = points();
    const faceList& oldFaces = faces();
    const cellList& oldCells = cells();
    const label nOldInternalFaces = nInternalFaces();
    const polyPatchList& oldPatches = boundaryMesh();

    // Mirror the points
Mark Olesen's avatar
Mark Olesen committed
62
    Info<< "Mirroring points. Old points: " << oldPoints.size();
63
64
65
66
67
68
69

    pointField newPoints(2*oldPoints.size());
    label nNewPoints = 0;

    labelList mirrorPointLookup(oldPoints.size(), -1);

    // Grab the old points
70
    forAll(oldPoints, pointI)
71
72
73
74
75
    {
        newPoints[nNewPoints] = oldPoints[pointI];
        nNewPoints++;
    }

76
    forAll(oldPoints, pointI)
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
    {
        scalar alpha =
            mirrorPlane.normalIntersect
            (
                oldPoints[pointI],
                mirrorPlane.normal()
            );

        // Check plane on tolerance
        if (mag(alpha) > planeTolerance)
        {
            // The point gets mirrored
            newPoints[nNewPoints] =
                oldPoints[pointI] + 2.0*alpha*mirrorPlane.normal();

            // remember the point correspondence
            mirrorPointLookup[pointI] = nNewPoints;
            nNewPoints++;
        }
        else
        {
            // The point is on the plane and does not get mirrored
            // Adjust plane location
            newPoints[nNewPoints] =
                oldPoints[pointI] + alpha*mirrorPlane.normal();

            mirrorPointLookup[pointI] = pointI;
        }
    }

    // Reset the size of the point list
Mark Olesen's avatar
Mark Olesen committed
108
    Info<< " New points: " << nNewPoints << endl;
109
110
    newPoints.setSize(nNewPoints);

Mark Olesen's avatar
Mark Olesen committed
111
    Info<< "Mirroring faces. Old faces: " << oldFaces.size();
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132

    // Algorithm:
    // During mirroring, the faces that were previously boundary faces
    // in the mirror plane may become ineternal faces. In order to
    // deal with the ordering of the faces, the algorithm is split
    // into two parts.  For original faces, the internal faces are
    // distributed to their owner cells.  Once all internal faces are
    // distributed, the boundary faces are visited and if they are in
    // the mirror plane they are added to the master cells (the future
    // boundary faces are not touched).  After the first phase, the
    // internal faces are collected in the cell order and numbering
    // information is added.  Then, the internal faces are mirrored
    // and the face numbering data is stored for the mirrored section.
    // Once all the internal faces are mirrored, the boundary faces
    // are added by mirroring the faces patch by patch.

    // Distribute internal faces
    labelListList newCellFaces(oldCells.size());

    const unallocLabelList& oldOwnerStart = lduAddr().ownerStartAddr();

133
    forAll(newCellFaces, cellI)
134
135
136
137
138
139
140
141
    {
        labelList& curFaces = newCellFaces[cellI];

        const label s = oldOwnerStart[cellI];
        const label e = oldOwnerStart[cellI + 1];

        curFaces.setSize(e - s);

142
        forAll(curFaces, i)
143
144
145
146
147
148
149
150
151
        {
            curFaces[i] = s + i;
        }
    }

    // Distribute boundary faces.  Remember the faces that have been inserted
    // as internal
    boolListList insertedBouFace(oldPatches.size());

152
    forAll(oldPatches, patchI)
153
154
    {
        const polyPatch& curPatch = oldPatches[patchI];
155
156
157
158
159
160
161
162
163
164

        if (curPatch.coupled())
        {
            WarningIn("mirrorFvMesh::mirrorFvMesh(const IOobject&)")
                << "Found coupled patch " << curPatch.name() << endl
                << "    Mirroring faces on coupled patches destroys"
                << " the ordering. This might be fixed by running a dummy"
                << " createPatch afterwards." << endl;
        }

165
166
167
168
169
170
171
172
173
        boolList& curInsBouFace = insertedBouFace[patchI];

        curInsBouFace.setSize(curPatch.size());
        curInsBouFace = false;

        // Get faceCells for face insertion
        const unallocLabelList& curFaceCells = curPatch.faceCells();
        const label curStart = curPatch.start();

174
        forAll(curPatch, faceI)
175
176
177
178
179
180
181
        {
            // Find out if the mirrored face is identical to the
            // original.  If so, the face needs to become internal and
            // added to its owner cell
            const face& origFace = curPatch[faceI];

            face mirrorFace(origFace.size());
182
            forAll(mirrorFace, pointI)
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
            {
                mirrorFace[pointI] = mirrorPointLookup[origFace[pointI]];
            }

            if (origFace == mirrorFace)
            {
                // The mirror is identical to current face.  This will
                // become an internal face
                const label oldSize = newCellFaces[curFaceCells[faceI]].size();

                newCellFaces[curFaceCells[faceI]].setSize(oldSize + 1);
                newCellFaces[curFaceCells[faceI]][oldSize] = curStart + faceI;

                curInsBouFace[faceI] = true;
            }
        }
    }

    // Construct the new list of faces.  Boundary faces are added
    // last, cush that each patch is mirrored separately.  The
    // addressing is stored in two separate arrays: first for the
    // original cells (face order has changed) and then for the
    // mirrored cells.
    labelList masterFaceLookup(oldFaces.size(), -1);
    labelList mirrorFaceLookup(oldFaces.size(), -1);

    faceList newFaces(2*oldFaces.size());
    label nNewFaces = 0;

    // Insert original (internal) faces
213
    forAll(newCellFaces, cellI)
214
215
216
    {
        const labelList& curCellFaces = newCellFaces[cellI];

217
        forAll(curCellFaces, cfI)
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
        {
            newFaces[nNewFaces] = oldFaces[curCellFaces[cfI]];
            masterFaceLookup[curCellFaces[cfI]] = nNewFaces;

            nNewFaces++;
        }
    }

    // Mirror internal faces
    for (label faceI = 0; faceI < nOldInternalFaces; faceI++)
    {
        const face& oldFace = oldFaces[faceI];
        face& nf = newFaces[nNewFaces];
        nf.setSize(oldFace.size());

        nf[0] = mirrorPointLookup[oldFace[0]];

        for (label i = 1; i < oldFace.size(); i++)
        {
            nf[i] = mirrorPointLookup[oldFace[oldFace.size() - i]];
        }

        mirrorFaceLookup[faceI] = nNewFaces;
        nNewFaces++;
    }

    // Mirror boundary faces patch by patch

    wordList newPatchTypes(boundary().size());
    wordList newPatchNames(boundary().size());
    labelList newPatchSizes(boundary().size(), -1);
    labelList newPatchStarts(boundary().size(), -1);
    label nNewPatches = 0;

252
    forAll(boundaryMesh(), patchI)
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
    {
        const label curPatchSize = boundaryMesh()[patchI].size();
        const label curPatchStart = boundaryMesh()[patchI].start();
        const boolList& curInserted = insertedBouFace[patchI];

        newPatchStarts[nNewPatches] = nNewFaces;

        // Master side
        for (label faceI = 0; faceI < curPatchSize; faceI++)
        {
            // Check if the face has already been added.  If not, add it and
            // insert the numbering details.
            if (!curInserted[faceI])
            {
                newFaces[nNewFaces] = oldFaces[curPatchStart + faceI];

                masterFaceLookup[curPatchStart + faceI] = nNewFaces;
                nNewFaces++;
            }
        }

        // Mirror side
        for (label faceI = 0; faceI < curPatchSize; faceI++)
        {
            // Check if the face has already been added.  If not, add it and
            // insert the numbering details.
            if (!curInserted[faceI])
            {
                const face& oldFace = oldFaces[curPatchStart + faceI];
                face& nf = newFaces[nNewFaces];
                nf.setSize(oldFace.size());

                nf[0] = mirrorPointLookup[oldFace[0]];

                for (label i = 1; i < oldFace.size(); i++)
                {
                    nf[i] = mirrorPointLookup[oldFace[oldFace.size() - i]];
                }

                mirrorFaceLookup[curPatchStart + faceI] = nNewFaces;
                nNewFaces++;
            }
            else
            {
                // Grab the index of the master face for the mirror side
                mirrorFaceLookup[curPatchStart + faceI] =
                    masterFaceLookup[curPatchStart + faceI];
            }
        }

        // If patch exists, grab the name and type of the original patch
        if (nNewFaces > newPatchStarts[nNewPatches])
        {
            newPatchTypes[nNewPatches] = boundaryMesh()[patchI].type();
            newPatchNames[nNewPatches] = boundaryMesh()[patchI].name();
            newPatchSizes[nNewPatches] =
                nNewFaces - newPatchStarts[nNewPatches];

            nNewPatches++;
        }
    }

    // Tidy up the lists
    newFaces.setSize(nNewFaces);
Mark Olesen's avatar
Mark Olesen committed
317
    Info<< " New faces: " << nNewFaces << endl;
318
319
320
321
322
323

    newPatchTypes.setSize(nNewPatches);
    newPatchNames.setSize(nNewPatches);
    newPatchSizes.setSize(nNewPatches);
    newPatchStarts.setSize(nNewPatches);

Mark Olesen's avatar
Mark Olesen committed
324
    Info<< "Mirroring patches. Old patches: " << boundary().size()
325
326
327
328
329
330
331
332
333
        << " New patches: " << nNewPatches << endl;

    Info<< "Mirroring cells.  Old cells: " << oldCells.size()
        << " New cells: " << 2*oldCells.size() << endl;

    cellList newCells(2*oldCells.size());
    label nNewCells = 0;

    // Grab the original cells.  Take care of face renumbering.
334
    forAll(oldCells, cellI)
335
336
337
338
339
340
    {
        const cell& oc = oldCells[cellI];

        cell& nc = newCells[nNewCells];
        nc.setSize(oc.size());

341
        forAll(oc, i)
342
343
344
345
346
347
348
349
        {
            nc[i] = masterFaceLookup[oc[i]];
        }

        nNewCells++;
    }

    // Mirror the cells
350
    forAll(oldCells, cellI)
351
352
353
354
355
356
    {
        const cell& oc = oldCells[cellI];

        cell& nc = newCells[nNewCells];
        nc.setSize(oc.size());

357
        forAll(oc, i)
358
359
360
361
362
363
364
365
        {
            nc[i] = mirrorFaceLookup[oc[i]];
        }

        nNewCells++;
    }

    // Mirror the cell shapes
Mark Olesen's avatar
Mark Olesen committed
366
    Info<< "Mirroring cell shapes." << endl;
367

Mark Olesen's avatar
Mark Olesen committed
368
    Info<< nl << "Creating new mesh" << endl;
369
370
371
    mirrorMeshPtr_ = new fvMesh
    (
        io,
372
373
374
        xferMove(newPoints),
        xferMove(newFaces),
        xferMove(newCells)
375
376
377
378
379
380
381
    );

    fvMesh& pMesh = *mirrorMeshPtr_;

    // Add the boundary patches
    List<polyPatch*> p(newPatchTypes.size());

382
    forAll(p, patchI)
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
    {
        p[patchI] = polyPatch::New
        (
            newPatchTypes[patchI],
            newPatchNames[patchI],
            newPatchSizes[patchI],
            newPatchStarts[patchI],
            patchI,
            pMesh.boundaryMesh()
        ).ptr();
    }

    pMesh.addPatches(p);
}


// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //

Foam::mirrorFvMesh::~mirrorFvMesh()
{}


// ************************************************************************* //