patchProbes.C 7.54 KB
Newer Older
Henry's avatar
Henry committed
1
2
3
4
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
5
    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
Henry's avatar
Henry committed
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
     \\/     M anipulation  |
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

    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.

    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
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.

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

#include "patchProbes.H"
#include "volFields.H"
#include "IOmanip.H"
#include "mappedPatchBase.H"
#include "treeBoundBox.H"
#include "treeDataFace.H"
32
33
#include "addToRunTimeSelectionTable.H"

Henry's avatar
Henry committed
34
35
36
37
38
39

// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //

namespace Foam
{
    defineTypeNameAndDebug(patchProbes, 0);
40
41
42
43
44
45
46

    addToRunTimeSelectionTable
    (
        functionObject,
        patchProbes,
        dictionary
    );
Henry's avatar
Henry committed
47
48
49
50
51
52
53
54
55
56
}

// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //

void Foam::patchProbes::findElements(const fvMesh& mesh)
{
    (void)mesh.tetBasePtIs();

    const polyBoundaryMesh& bm = mesh.boundaryMesh();

57
    label patchi = bm.findPatchID(patchName_);
Henry's avatar
Henry committed
58

59
    if (patchi == -1)
Henry's avatar
Henry committed
60
    {
61
62
        FatalErrorInFunction
            << " Unknown patch name "
Henry's avatar
Henry committed
63
64
65
66
67
68
69
            << patchName_ << endl
            << exit(FatalError);
    }

     // All the info for nearest. Construct to miss
    List<mappedPatchBase::nearInfo> nearest(this->size());

70
    const polyPatch& pp = bm[patchi];
Henry's avatar
Henry committed
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100

    if (pp.size() > 0)
    {
        labelList bndFaces(pp.size());
        forAll(bndFaces, i)
        {
            bndFaces[i] =  pp.start() + i;
        }

        treeBoundBox overallBb(pp.points());
        Random rndGen(123456);
        overallBb = overallBb.extend(rndGen, 1e-4);
        overallBb.min() -= point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);
        overallBb.max() += point(ROOTVSMALL, ROOTVSMALL, ROOTVSMALL);

        const indexedOctree<treeDataFace> boundaryTree
        (
            treeDataFace    // all information needed to search faces
            (
                false,                      // do not cache bb
                mesh,
                bndFaces                    // patch faces only
            ),
            overallBb,                      // overall search domain
            8,                              // maxLevel
            10,                             // leafsize
            3.0                             // duplicity
        );


101
        forAll(probeLocations(), probei)
Henry's avatar
Henry committed
102
        {
103
            const point sample = probeLocations()[probei];
Henry's avatar
Henry committed
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121

            scalar span = boundaryTree.bb().mag();

            pointIndexHit info = boundaryTree.findNearest
            (
                sample,
                Foam::sqr(span)
            );

            if (!info.hit())
            {
                info = boundaryTree.findNearest
                (
                    sample,
                    Foam::sqr(GREAT)
                );
            }

122
            label facei = boundaryTree.shapes().faceLabels()[info.index()];
Henry's avatar
Henry committed
123

124
            const label patchi = bm.whichPatch(facei);
Henry's avatar
Henry committed
125
126
127

            if (isA<emptyPolyPatch>(bm[patchi]))
            {
128
                WarningInFunction
Henry's avatar
Henry committed
129
130
131
132
133
134
135
136
                << " The sample point: " << sample
                << " belongs to " << patchi
                << " which is an empty patch. This is not permitted. "
                << " This sample will not be included "
                << endl;
            }
            else
            {
137
                const point& fc = mesh.faceCentres()[facei];
Henry's avatar
Henry committed
138
139
140
141
142
143
144

                mappedPatchBase::nearInfo sampleInfo;

                sampleInfo.first() = pointIndexHit
                (
                    true,
                    fc,
145
                    facei
Henry's avatar
Henry committed
146
147
148
149
150
                );

                sampleInfo.second().first() = magSqr(fc-sample);
                sampleInfo.second().second() = Pstream::myProcNo();

151
                nearest[probei]= sampleInfo;
Henry's avatar
Henry committed
152
153
154
155
156
157
158
159
160
161
162
            }
        }
    }


    // Find nearest.
    Pstream::listCombineGather(nearest, mappedPatchBase::nearestEqOp());
    Pstream::listCombineScatter(nearest);

    if (debug)
    {
163
        InfoInFunction << endl;
Henry's avatar
Henry committed
164
165
        forAll(nearest, sampleI)
        {
166
            label proci = nearest[sampleI].second().second();
Henry's avatar
Henry committed
167
168
169
            label localI = nearest[sampleI].first().index();

            Info<< "    " << sampleI << " coord:"<< operator[](sampleI)
170
                << " found on processor:" << proci
Henry's avatar
Henry committed
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
                << " in local cell/face:" << localI
                << " with fc:" << nearest[sampleI].first().rawPoint() << endl;
        }
    }


    // Extract any local faces to sample
    elementList_.setSize(nearest.size(), -1);

    forAll(nearest, sampleI)
    {
        if (nearest[sampleI].second().second() == Pstream::myProcNo())
        {
            // Store the face to sample
            elementList_[sampleI] = nearest[sampleI].first().index();
        }
    }
}


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

193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
Foam::patchProbes::patchProbes
(
    const word& name,
    const Time& t,
    const dictionary& dict
)
:
    probes(name, t, dict)
{
    // When constructing probes above it will have called the
    // probes::findElements (since the virtual mechanism not yet operating).
    // Not easy to workaround (apart from feeding through flag into constructor)
    // so clear out any cells found for now.
    elementList_.clear();
    faceList_.clear();

    read(dict);
}


Henry's avatar
Henry committed
213
214
215
216
217
218
219
220
Foam::patchProbes::patchProbes
(
    const word& name,
    const objectRegistry& obr,
    const dictionary& dict,
    const bool loadFromFiles
)
:
221
    probes(name, obr, dict)
Henry's avatar
Henry committed
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
{
    // When constructing probes above it will have called the
    // probes::findElements (since the virtual mechanism not yet operating).
    // Not easy to workaround (apart from feeding through flag into constructor)
    // so clear out any cells found for now.
    elementList_.clear();
    faceList_.clear();

    read(dict);
}


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

Foam::patchProbes::~patchProbes()
{}


240
bool Foam::patchProbes::write(const bool postProcess)
Henry's avatar
Henry committed
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
{
    if (this->size() && prepare())
    {
        sampleAndWrite(scalarFields_);
        sampleAndWrite(vectorFields_);
        sampleAndWrite(sphericalTensorFields_);
        sampleAndWrite(symmTensorFields_);
        sampleAndWrite(tensorFields_);

        sampleAndWriteSurfaceFields(surfaceScalarFields_);
        sampleAndWriteSurfaceFields(surfaceVectorFields_);
        sampleAndWriteSurfaceFields(surfaceSphericalTensorFields_);
        sampleAndWriteSurfaceFields(surfaceSymmTensorFields_);
        sampleAndWriteSurfaceFields(surfaceTensorFields_);
    }
256
257

    return true;
Henry's avatar
Henry committed
258
259
}

260
261

bool Foam::patchProbes::read(const dictionary& dict)
Henry's avatar
Henry committed
262
263
{
    dict.lookup("patchName") >> patchName_;
264
    return probes::read(dict);
Henry's avatar
Henry committed
265
266
267
268
}


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