1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
11 OpenFOAM is free software: you can redistribute it and/or modify it
12 under the terms of the GNU General Public License as published by
13 the Free Software Foundation, either version 3 of the License, or
14 (at your option) any later version.
16 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
24 \*---------------------------------------------------------------------------*/
26 #include "distanceSurface.H"
27 #include "dictionary.H"
28 #include "volFields.H"
29 #include "volPointInterpolation.H"
30 #include "addToRunTimeSelectionTable.H"
32 #include "isoSurface.H"
33 // #include "isoSurfaceCell.H"
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
39 defineTypeNameAndDebug(distanceSurface, 0);
40 addToRunTimeSelectionTable(sampledSurface, distanceSurface, word);
43 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
45 void Foam::distanceSurface::createGeometry()
49 Pout<< "distanceSurface::createGeometry :updating geometry." << endl;
52 // Clear any stored topologies
58 const fvMesh& fvm = static_cast<const fvMesh&>(mesh());
60 // Distance to cell centres
61 // ~~~~~~~~~~~~~~~~~~~~~~~~
63 cellDistancePtr_.reset
70 fvm.time().timeName(),
77 dimensionedScalar("zero", dimLength, 0)
80 volScalarField& cellDistance = cellDistancePtr_();
84 const pointField& cc = fvm.C();
85 scalarField& fld = cellDistance.internalField();
87 List<pointIndexHit> nearest;
88 surfPtr_().findNearest
91 scalarField(cc.size(), GREAT),
98 surfPtr_().getNormal(nearest, normal);
102 vector d(cc[i]-nearest[i].hitPoint());
104 if ((d&normal[i]) > 0)
106 fld[i] = Foam::mag(d);
110 fld[i] = -Foam::mag(d);
118 fld[i] = Foam::mag(cc[i] - nearest[i].hitPoint());
125 forAll(fvm.C().boundaryField(), patchI)
127 const pointField& cc = fvm.C().boundaryField()[patchI];
128 fvPatchScalarField& fld = cellDistance.boundaryField()[patchI];
130 List<pointIndexHit> nearest;
131 surfPtr_().findNearest
134 scalarField(cc.size(), GREAT),
141 surfPtr_().getNormal(nearest, normal);
145 vector d(cc[i]-nearest[i].hitPoint());
147 if ((d&normal[i]) > 0)
149 fld[i] = Foam::mag(d);
153 fld[i] = -Foam::mag(d);
161 fld[i] = Foam::mag(cc[i] - nearest[i].hitPoint());
168 // On processor patches the mesh.C() will already be the cell centre
169 // on the opposite side so no need to swap cellDistance.
172 // Distance to points
173 pointDistance_.setSize(fvm.nPoints());
175 const pointField& pts = fvm.points();
177 List<pointIndexHit> nearest;
178 surfPtr_().findNearest
181 scalarField(pts.size(), GREAT),
188 surfPtr_().getNormal(nearest, normal);
192 vector d(pts[i]-nearest[i].hitPoint());
194 if ((d&normal[i]) > 0)
196 pointDistance_[i] = Foam::mag(d);
200 pointDistance_[i] = -Foam::mag(d);
208 pointDistance_[i] = Foam::mag(pts[i]-nearest[i].hitPoint());
216 Pout<< "Writing cell distance:" << cellDistance.objectPath() << endl;
217 cellDistance.write();
218 pointScalarField pDist
223 fvm.time().timeName(),
230 dimensionedScalar("zero", dimLength, 0)
232 pDist.internalField() = pointDistance_;
234 Pout<< "Writing point distance:" << pDist.objectPath() << endl;
239 //- Direct from cell field and point field.
259 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
261 Foam::distanceSurface::distanceSurface
264 const polyMesh& mesh,
265 const dictionary& dict
268 sampledSurface(name, mesh, dict),
271 searchableSurface::New
273 dict.lookup("surfaceType"),
276 dict.lookupOrDefault("surfaceName", name), // name
277 mesh.time().constant(), // directory
278 "triSurface", // instance
279 mesh.time(), // registry
286 distance_(readScalar(dict.lookup("distance"))),
287 signed_(readBool(dict.lookup("signed"))),
288 regularise_(dict.lookupOrDefault("regularise", true)),
289 zoneName_(word::null),
294 // dict.readIfPresent("zone", zoneName_);
296 // if (debug && zoneName_.size())
298 // if (mesh.cellZones().findZoneID(zoneName_) < 0)
300 // Info<< "cellZone \"" << zoneName_
301 // << "\" not found - using entire mesh" << endl;
307 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
309 Foam::distanceSurface::~distanceSurface()
313 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
315 bool Foam::distanceSurface::needsUpdate() const
321 bool Foam::distanceSurface::expire()
325 Pout<< "distanceSurface::expire :"
326 << " have-facesPtr_:" << facesPtr_.valid()
327 << " needsUpdate_:" << needsUpdate_ << endl;
330 // Clear any stored topologies
333 // Clear derived data
336 // already marked as expired
347 bool Foam::distanceSurface::update()
351 Pout<< "distanceSurface::update :"
352 << " have-facesPtr_:" << facesPtr_.valid()
353 << " needsUpdate_:" << needsUpdate_ << endl;
363 needsUpdate_ = false;
368 Foam::tmp<Foam::scalarField>
369 Foam::distanceSurface::sample
371 const volScalarField& vField
374 return sampleField(vField);
378 Foam::tmp<Foam::vectorField>
379 Foam::distanceSurface::sample
381 const volVectorField& vField
384 return sampleField(vField);
388 Foam::tmp<Foam::sphericalTensorField>
389 Foam::distanceSurface::sample
391 const volSphericalTensorField& vField
394 return sampleField(vField);
398 Foam::tmp<Foam::symmTensorField>
399 Foam::distanceSurface::sample
401 const volSymmTensorField& vField
404 return sampleField(vField);
408 Foam::tmp<Foam::tensorField>
409 Foam::distanceSurface::sample
411 const volTensorField& vField
414 return sampleField(vField);
418 Foam::tmp<Foam::scalarField>
419 Foam::distanceSurface::interpolate
421 const interpolation<scalar>& interpolator
424 return interpolateField(interpolator);
428 Foam::tmp<Foam::vectorField>
429 Foam::distanceSurface::interpolate
431 const interpolation<vector>& interpolator
434 return interpolateField(interpolator);
437 Foam::tmp<Foam::sphericalTensorField>
438 Foam::distanceSurface::interpolate
440 const interpolation<sphericalTensor>& interpolator
443 return interpolateField(interpolator);
447 Foam::tmp<Foam::symmTensorField>
448 Foam::distanceSurface::interpolate
450 const interpolation<symmTensor>& interpolator
453 return interpolateField(interpolator);
457 Foam::tmp<Foam::tensorField>
458 Foam::distanceSurface::interpolate
460 const interpolation<tensor>& interpolator
463 return interpolateField(interpolator);
467 void Foam::distanceSurface::print(Ostream& os) const
469 os << "distanceSurface: " << name() << " :"
470 << " surface:" << surfPtr_().name()
471 << " distance:" << distance_
472 << " faces:" << faces().size()
473 << " points:" << points().size();
477 // ************************************************************************* //