From 54b5b660231607c5c4647e2cc1d3eeee2425931d Mon Sep 17 00:00:00 2001
From: mattijs <mattijs>
Date: Tue, 7 Feb 2012 08:47:50 +0000
Subject: [PATCH] ENH: wallBoundedStreamLine: updated motorBike tutorial.

---
 .../field/nearWallFields/controlDict          |   1 +
 .../field/nearWallFields/nearWallFields.C     |  25 ++--
 .../field/streamLine/streamLine.C             | 131 +-----------------
 .../field/streamLine/streamLine.H             |   7 +-
 .../field/streamLine/streamLineParticle.C     |  51 +------
 .../field/streamLine/streamLineParticle.H     |  16 +--
 .../wallBoundedStreamLine.C                   |   4 +-
 .../simpleFoam/motorBike/system/controlDict   |   1 +
 .../motorBike/system/wallBoundedStreamLines   |  96 +++++++++++++
 9 files changed, 122 insertions(+), 210 deletions(-)
 create mode 100644 tutorials/incompressible/simpleFoam/motorBike/system/wallBoundedStreamLines

diff --git a/src/postProcessing/functionObjects/field/nearWallFields/controlDict b/src/postProcessing/functionObjects/field/nearWallFields/controlDict
index 3c8d512f12..a4f7bec09c 100644
--- a/src/postProcessing/functionObjects/field/nearWallFields/controlDict
+++ b/src/postProcessing/functionObjects/field/nearWallFields/controlDict
@@ -59,6 +59,7 @@ functions
 
         // Fields to be sampled. Per field original name and mapped field to
         // create.
+        // Note: fields only get updated when writing!
         fields
         (
             (p pNear)
diff --git a/src/postProcessing/functionObjects/field/nearWallFields/nearWallFields.C b/src/postProcessing/functionObjects/field/nearWallFields/nearWallFields.C
index d586dfba0c..ee8eb8e334 100644
--- a/src/postProcessing/functionObjects/field/nearWallFields/nearWallFields.C
+++ b/src/postProcessing/functionObjects/field/nearWallFields/nearWallFields.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -121,6 +121,7 @@ void Foam::nearWallFields::read(const dictionary& dict)
             reverseFieldMap_.insert(sampleFldName, fldName);
         }
 
+        Info<< "Creating " << fieldMap_.size() << " fields" << endl;
         createFields(vsf_);
         createFields(vvf_);
         createFields(vSpheretf_);
@@ -137,14 +138,14 @@ void Foam::nearWallFields::execute()
         Info<< "nearWallFields:execute()" << endl;
     }
 
-    if (active_)
-    {
-        sampleFields(vsf_);
-        sampleFields(vvf_);
-        sampleFields(vSpheretf_);
-        sampleFields(vSymmtf_);
-        sampleFields(vtf_);
-    }
+    //if (active_)
+    //{
+    //    sampleFields(vsf_);
+    //    sampleFields(vvf_);
+    //    sampleFields(vSpheretf_);
+    //    sampleFields(vSymmtf_);
+    //    sampleFields(vtf_);
+    //}
 }
 
 
@@ -172,6 +173,12 @@ void Foam::nearWallFields::write()
         Info<< "Writing sampled fields to " << obr_.time().timeName()
             << endl;
 
+        sampleFields(vsf_);
+        sampleFields(vvf_);
+        sampleFields(vSpheretf_);
+        sampleFields(vSymmtf_);
+        sampleFields(vtf_);
+
         // Write fields
         forAll(vsf_, i)
         {
diff --git a/src/postProcessing/functionObjects/field/streamLine/streamLine.C b/src/postProcessing/functionObjects/field/streamLine/streamLine.C
index ab815a8df8..13400cbcef 100644
--- a/src/postProcessing/functionObjects/field/streamLine/streamLine.C
+++ b/src/postProcessing/functionObjects/field/streamLine/streamLine.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -292,119 +292,6 @@ void Foam::streamLine::track()
     }
 
 
-    //- For surface following: per face the normal to constrain the velocity
-    //  with.
-    pointField patchNormals;
-
-
-    if (followSurface_)
-    {
-        // Determine geometric data on the non-constraint patches
-        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-        // - point normals
-        // - neighbouring cells
-        autoPtr<indirectPrimitivePatch> boundaryPatch(wallPatch());
-
-        // Calculate parallel consistent normal (on boundaryPatch points only)
-        pointField boundaryNormals
-        (
-            PatchTools::pointNormals
-            (
-                mesh,
-                boundaryPatch(),
-                boundaryPatch().addressing()
-            )
-        );
-        // Make sure all points know about point normal. There might be points
-        // which are not on boundaryPatch() but coupled to points that are.
-        pointField pointNormals(mesh.nPoints(), vector::zero);
-        UIndirectList<point>
-        (
-            pointNormals,
-            boundaryPatch().meshPoints()
-        ) = boundaryNormals;
-        boundaryNormals.clear();
-
-        mesh.globalData().syncPointData
-        (
-            pointNormals,
-            maxMagSqrEqOp<point>(),
-            mapDistribute::transform()
-        );
-
-        // Calculate per-face a single constraint normal:
-        // - faces on patch: face-normal
-        // - faces on cell with point on patch: point-normal
-        // - rest: vector::zero
-        //
-        // Note: the disadvantage is that we can only handle faces with one
-        //       point on the patch - otherwise we might pick up the wrong
-        //       velocity.
-
-
-
-        patchNormals.setSize(mesh.nFaces(), vector::zero);
-        label nPointNormals = 0;
-
-        // 1. faces on patch
-        UIndirectList<point>(patchNormals, boundaryPatch().addressing()) =
-             boundaryPatch().faceNormals();
-
-        // 2. faces with a point on the patch
-        forAll(mesh.faces(), faceI)
-        {
-            if (patchNormals[faceI] == vector::zero)
-            {
-                const face& f = mesh.faces()[faceI];
-                forAll(f, fp)
-                {
-                    label pointI = f[fp];
-                    if (pointNormals[pointI] != vector::zero)
-                    {
-                        patchNormals[faceI] = pointNormals[pointI];
-                        nPointNormals++;
-                    }
-                }
-            }
-        }
-
-        // 3. cells on faces with a point on the patch
-        forAll(mesh.points(), pointI)
-        {
-            if (pointNormals[pointI] != vector::zero)
-            {
-                const labelList& pCells = mesh.pointCells(pointI);
-                forAll(pCells, i)
-                {
-                    const cell& cFaces = mesh.cells()[pCells[i]];
-                    forAll(cFaces, j)
-                    {
-                        label faceI = cFaces[j];
-
-                        if (patchNormals[faceI] == vector::zero)
-                        {
-                            patchNormals[faceI] = pointNormals[pointI];
-                            nPointNormals++;
-                        }
-                    }
-                }
-            }
-        }
-
-        Info<< type() << " : calculated constrained-tracking normals" << nl
-            << "Number of faces in mesh                     : "
-            << returnReduce(mesh.nFaces(), sumOp<label>()) << nl
-            << "    of which on walls                       : "
-            <<  returnReduce
-                (
-                    boundaryPatch().addressing().size(),
-                    sumOp<label>()
-                ) << nl
-            << "    of which on cell with point on boundary : "
-            << returnReduce(nPointNormals, sumOp<label>()) << nl
-            << endl;
-    }
-
     // additional particle info
     streamLineParticle::trackingData td
     (
@@ -416,12 +303,6 @@ void Foam::streamLine::track()
         nSubCycle_,     // automatic track control:step through cells in steps?
         trackLength_,   // fixed track length
 
-        //- For surface-constrained streamlines
-        followSurface_, // stay on surface?
-        patchNormals,   // per face normal to constrain tracking to
-        //pointNormals,   // per point          ,,
-        trackHeight_,   // track height
-
         allTracks_,
         allScalars_,
         allVectors_
@@ -560,16 +441,6 @@ void Foam::streamLine::read(const dictionary& dict)
         }
 
 
-        followSurface_ = false;
-        trackHeight_ = VGREAT;
-        if (dict.readIfPresent("followSurface", followSurface_))
-        {
-            dict.lookup("trackHeight") >> trackHeight_;
-
-            Info<< type() << " : surface-constrained streamline at "
-                << trackHeight_ << " distance above the surface" << nl << endl;
-        }
-
         interpolationScheme_ = dict.lookupOrDefault
         (
             "interpolationScheme",
diff --git a/src/postProcessing/functionObjects/field/streamLine/streamLine.H b/src/postProcessing/functionObjects/field/streamLine/streamLine.H
index 77984ed2c8..f4fb4e42dc 100644
--- a/src/postProcessing/functionObjects/field/streamLine/streamLine.H
+++ b/src/postProcessing/functionObjects/field/streamLine/streamLine.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -96,11 +96,6 @@ class streamLine
         //- Maximum lifetime (= number of cells) of particle
         label lifeTime_;
 
-        //- Surface-constrained?
-        bool followSurface_;
-
-        scalar trackHeight_;
-
         //- Number of subcycling steps
         label nSubCycle_;
 
diff --git a/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.C b/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.C
index 279e17ef80..cf1f749baf 100644
--- a/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.C
+++ b/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -211,17 +211,6 @@ bool Foam::streamLineParticle::move
                 U = -U;
             }
 
-            if (td.followSurface_)
-            {
-                //- Simple: remove wall-normal component. Should keep particle
-                //  at same height. Does not work well on warped faces.
-                const vector& n = td.patchNormals_[tetFaceI_];
-                if (n != vector::zero)
-                {
-                    U -= (U&n)*n;
-                }
-            }
-
             scalar magU = mag(U);
 
             if (magU < SMALL)
@@ -396,42 +385,8 @@ void Foam::streamLineParticle::hitWallPatch
     const tetIndices&
 )
 {
-    if (td.followSurface_)
-    {
-        // Push slightly in. Since we're pushing to tet centre we should still
-        // be in the same tet.
-        tetPointRef tet = currentTet();
-
-        const point oldPosition(position());
-
-        position() += trackingCorrectionTol*(tet.centre() - position());
-
-        if (debug)
-        {
-            label patchID = patch(face());
-
-            Pout<< "hitWallPatch : on wall "
-                << mesh().boundaryMesh()[patchID].name()
-                << " moved from " << oldPosition
-                << " to " << position()
-                << " due to centre " << tet.centre() << endl;
-        }
-
-        // Check that still in tet (should always be the case
-        if (debug && !tet.inside(position()))
-        {
-            FatalErrorIn("streamLineParticle::hitWallPatch()")
-                << "Pushing particle " << *this
-                << " away from wall " << wpp.name()
-                << " moved it outside of tet."
-                << exit(FatalError);
-        }
-    }
-    else
-    {
-        // Remove particle
-        td.keepParticle = false;
-    }
+    // Remove particle
+    td.keepParticle = false;
 }
 
 
diff --git a/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.H b/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.H
index c5feff4285..02bd609933 100644
--- a/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.H
+++ b/src/postProcessing/functionObjects/field/streamLine/streamLineParticle.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -75,12 +75,6 @@ public:
         const label nSubCycle_;
         const scalar trackLength_;
 
-        //- For surface-constrained tracking
-        const bool followSurface_;
-        const pointField& patchNormals_;
-        const scalar trackHeight_;
-
-
         DynamicList<vectorList>& allPositions_;
         List<DynamicList<scalarList> >& allScalars_;
         List<DynamicList<vectorList> >& allVectors_;
@@ -97,10 +91,6 @@ public:
                 const bool trackForward,
                 const label nSubCycle,
                 const scalar trackLength,
-                const bool followSurface,
-
-                const pointField& patchNormals,
-                const scalar trackHeight,
 
                 DynamicList<List<point> >& allPositions,
                 List<DynamicList<scalarList> >& allScalars,
@@ -114,10 +104,6 @@ public:
                 trackForward_(trackForward),
                 nSubCycle_(nSubCycle),
                 trackLength_(trackLength),
-                followSurface_(followSurface),
-
-                patchNormals_(patchNormals),
-                trackHeight_(trackHeight),
 
                 allPositions_(allPositions),
                 allScalars_(allScalars),
diff --git a/src/postProcessing/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLine.C b/src/postProcessing/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLine.C
index 4b6f6761f7..9089749eef 100644
--- a/src/postProcessing/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLine.C
+++ b/src/postProcessing/functionObjects/field/wallBoundedStreamLine/wallBoundedStreamLine.C
@@ -571,11 +571,11 @@ void Foam::wallBoundedStreamLine::read(const dictionary& dict)
             {
                 label nFaces = returnReduce(faces.size(), sumOp<label>());
 
-                FatalErrorIn("wallBoundedStreamLine::track()")
+                WarningIn("wallBoundedStreamLine::track()")
                     << "Found " << nFaces
                     <<" faces with low quality or negative volume "
                     << "decomposition tets. Writing to faceSet " << faces.name()
-                    << exit(FatalError);
+                    << endl;    //exit(FatalError);
             }
 
             // 2. all edges on a cell having two faces
diff --git a/tutorials/incompressible/simpleFoam/motorBike/system/controlDict b/tutorials/incompressible/simpleFoam/motorBike/system/controlDict
index ed21ad27f0..a740c9702e 100644
--- a/tutorials/incompressible/simpleFoam/motorBike/system/controlDict
+++ b/tutorials/incompressible/simpleFoam/motorBike/system/controlDict
@@ -48,6 +48,7 @@ functions
 {
     #include "readFields"
     #include "streamLines"
+    #include "wallBoundedStreamLines"
     #include "cuttingPlane"
     #include "forceCoeffs"
 }
diff --git a/tutorials/incompressible/simpleFoam/motorBike/system/wallBoundedStreamLines b/tutorials/incompressible/simpleFoam/motorBike/system/wallBoundedStreamLines
new file mode 100644
index 0000000000..09db1038e7
--- /dev/null
+++ b/tutorials/incompressible/simpleFoam/motorBike/system/wallBoundedStreamLines
@@ -0,0 +1,96 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+
+// Interpolate U to create near-wall UNear
+near
+{
+    // Where to load it from
+    functionObjectLibs ("libfieldFunctionObjects.so");
+
+    type nearWallFields;
+
+    // Output every
+    outputControl   outputTime;
+    //outputInterval  1;
+
+    // Fields to be sampled. Per field original name and mapped field to
+    // create.
+    fields
+    (
+        (U UNear)
+    );
+
+    // Patches/groups to sample (regular expressions)
+    patches (motorBike);
+
+    // Distance to sample
+    distance 0.001;
+}
+
+// Use UNear to track along wall
+wallBoundedStreamLines
+{
+    // Where to load it from (if not already in solver)
+    functionObjectLibs ("libfieldFunctionObjects.so");
+    type            wallBoundedStreamLine;
+
+    // Output every
+    outputControl   outputTime;
+    // outputInterval 10;
+
+    setFormat       vtk; //gnuplot; //xmgr; //raw; //jplot;
+
+    // Velocity field to use for tracking.
+    UName UNear;
+
+    // Interpolation method. Default is cellPoint. See sampleDict.
+    //interpolationScheme pointMVC;
+
+    // Tracked forwards (+U) or backwards (-U)
+    trackForward    true;
+
+    interpolationScheme cellPoint;
+
+    // Names of fields to sample. Should contain above velocity field!
+    fields (p U k UNear);
+
+    // Steps particles can travel before being removed
+    lifeTime        100;
+
+    // Cloud name to use
+    cloudName       wallBoundedParticleTracks;
+
+    // Seeding method. See the sampleSets in sampleDict.
+    seedSampleSet   patchSeed;    //cloud;//triSurfaceMeshPointSet;
+
+    uniformCoeffs
+    {
+        type        uniform;
+        axis        x;  //distance;
+
+        start       (0.0035 0.0999 0.0001);
+        end         (0.0035 0.0999 0.0099);
+        nPoints     20;
+    }
+    cloudCoeffs
+    {
+        type        cloud;
+        axis        x;  //distance;
+        points      ((0.351516548679288 -0.0116085375585099 1.24));
+    }
+    patchSeedCoeffs
+    {
+        type        patchSeed;//patchCloud; //cloud;  //uniform;
+        patches     (motorBike);
+        axis        x;  //distance;
+        maxPoints   20000;
+    }
+}
+
+
+// ************************************************************************* //