From 066b3158e7e1e46dcb03deed16ab70f9157dfdc7 Mon Sep 17 00:00:00 2001
From: Mark Olesen <Mark.Olesen@Germany>
Date: Fri, 4 Dec 2009 15:34:42 +0100
Subject: [PATCH] quickly implemented BSpline2 as an alternative B-Spline
implementation
- also looks reasonable and doesn't deviate much from Catmull-Rom
---
applications/test/spline/splineTest.C | 35 +++--
src/mesh/blockMesh/Make/files | 1 +
src/mesh/blockMesh/curvedEdges/BSpline2.C | 131 ++++++++++++++++++
src/mesh/blockMesh/curvedEdges/BSpline2.H | 126 +++++++++++++++++
.../blockMesh/curvedEdges/CatmullRomSpline.C | 2 +-
5 files changed, 286 insertions(+), 9 deletions(-)
create mode 100644 src/mesh/blockMesh/curvedEdges/BSpline2.C
create mode 100644 src/mesh/blockMesh/curvedEdges/BSpline2.H
diff --git a/applications/test/spline/splineTest.C b/applications/test/spline/splineTest.C
index c4c8f1434db..4c793c21361 100644
--- a/applications/test/spline/splineTest.C
+++ b/applications/test/spline/splineTest.C
@@ -27,7 +27,9 @@ License
#include "vector.H"
#include "IFstream.H"
+
#include "BSpline.H"
+#include "BSpline2.H"
#include "CatmullRomSpline.H"
using namespace Foam;
@@ -39,8 +41,9 @@ int main(int argc, char *argv[])
{
argList::noParallel();
argList::validArgs.insert("file .. fileN");
- argList::addBoolOption("B", "B-Spline");
- argList::addBoolOption("cmr", "catmull-rom spline (default)");
+ argList::addBoolOption("Bold", "B-Spline implementation OLD");
+ argList::addBoolOption("Bnew", "B-Spline implementation NEW");
+ argList::addBoolOption("CMR", "catmull-rom spline (default)");
argList::addOption
(
"n",
@@ -55,11 +58,12 @@ int main(int argc, char *argv[])
args.printUsage();
}
- bool useBSpline = args.optionFound("B");
- bool useCatmullRom = args.optionFound("cmr");
+ bool useBSplineOld = args.optionFound("Bold");
+ bool useBSplineNew = args.optionFound("Bnew");
+ bool useCatmullRom = args.optionFound("CMR");
label nSeg = args.optionLookupOrDefault<label>("n", 20);
- if (!useBSpline && !useCatmullRom)
+ if (!useCatmullRom && !useBSplineOld && !useBSplineNew)
{
Info<<"defaulting to Catmull-Rom spline" << endl;
useCatmullRom = true;
@@ -77,12 +81,27 @@ int main(int argc, char *argv[])
{
Info<<"\noriginal points: " << pointFields[splineI] << nl;
- if (useBSpline)
+ if (useBSplineOld)
+ {
+ BSpline spl(pointFields[splineI]);
+
+ Info<< nl
+ << "B-Spline interpolation: OLD" << nl
+ << "----------------------" << endl;
+
+ for (label segI = 0; segI <= nSeg; ++segI)
+ {
+ scalar lambda = scalar(segI)/scalar(nSeg);
+ Info<< spl.position(lambda) << " // " << lambda << endl;
+ }
+ }
+
+ if (useBSplineNew)
{
- BSpline spl(pointFields[splineI], vector::zero, vector::zero);
+ BSpline2 spl(pointFields[splineI]);
Info<< nl
- << "B-Spline interpolation:" << nl
+ << "B-Spline interpolation: NEW" << nl
<< "----------------------" << endl;
for (label segI = 0; segI <= nSeg; ++segI)
diff --git a/src/mesh/blockMesh/Make/files b/src/mesh/blockMesh/Make/files
index d68caa1c244..3b08f6062bc 100644
--- a/src/mesh/blockMesh/Make/files
+++ b/src/mesh/blockMesh/Make/files
@@ -1,3 +1,4 @@
+curvedEdges/BSpline2.C
curvedEdges/CatmullRomSpline.C
curvedEdges/polyLine.C
diff --git a/src/mesh/blockMesh/curvedEdges/BSpline2.C b/src/mesh/blockMesh/curvedEdges/BSpline2.C
new file mode 100644
index 00000000000..8208a597002
--- /dev/null
+++ b/src/mesh/blockMesh/curvedEdges/BSpline2.C
@@ -0,0 +1,131 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd.
+ \\/ 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 2 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, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "error.H"
+#include "BSpline2.H"
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+Foam::BSpline2::BSpline2
+(
+ const pointField& Knots,
+ const vector&,
+ const vector&
+)
+:
+ polyLine(Knots)
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
+
+Foam::point Foam::BSpline2::position(const scalar mu) const
+{
+ // endpoints
+ if (mu < SMALL)
+ {
+ return points().first();
+ }
+ else if (mu > 1 - SMALL)
+ {
+ return points().last();
+ }
+
+ scalar lambda = mu;
+ label segment = localParameter(lambda);
+ return position(segment, lambda);
+}
+
+
+Foam::point Foam::BSpline2::position
+(
+ const label segment,
+ const scalar mu
+) const
+{
+ const point& p0 = points()[segment];
+ const point& p1 = points()[segment+1];
+
+ // special cases - no calculation needed
+ if (segment < 0 || mu < 0.0)
+ {
+ return p0;
+ }
+ else if (segment > nSegments() || mu >= 1.0)
+ {
+ return p1;
+ }
+
+ // determine the end points
+ point e0;
+ point e1;
+
+ if (segment == 0)
+ {
+ // end: simple reflection
+ e0 = 2.0 * p0 - p1;
+ }
+ else
+ {
+ e0 = points()[segment-1];
+ }
+
+ if (segment+1 == nSegments())
+ {
+ // end: simple reflection
+ e1 = 2.0 * p1 - p0;
+ }
+ else
+ {
+ e1 = points()[segment+2];
+ }
+
+
+ return 1.0/6.0 *
+ (
+ ( e0 + 4*p0 + p1 )
+ + mu *
+ (
+ ( -3*e0 + 3*p1 )
+ + mu *
+ (
+ ( 3*e0 - 6*p0 + 3*p1 )
+ + mu *
+ ( -e0 + 3*p0 - 3*p1 + e1 )
+ )
+ )
+ );
+}
+
+
+Foam::scalar Foam::BSpline2::length() const
+{
+ notImplemented("BSpline2::length() const");
+ return 1.0;
+}
+
+
+// ************************************************************************* //
diff --git a/src/mesh/blockMesh/curvedEdges/BSpline2.H b/src/mesh/blockMesh/curvedEdges/BSpline2.H
new file mode 100644
index 00000000000..9aa5d1760e2
--- /dev/null
+++ b/src/mesh/blockMesh/curvedEdges/BSpline2.H
@@ -0,0 +1,126 @@
+/*---------------------------------------------------------------------------*\
+ ========= |
+ \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
+ \\ / O peration |
+ \\ / A nd | Copyright (C) 2009-2009 OpenCFD Ltd.
+ \\/ 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 2 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, write to the Free Software Foundation,
+ Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+ Foam::BSpline2
+
+Description
+ An implementation of B-splines.
+
+ In this implementation, the end tangents are created
+ automatically by reflection.
+
+ In matrix form, the @e local interpolation on the interval t=[0..1] is
+ described as follows:
+ @verbatim
+ P(t) = 1/6 * [ t^3 t^2 t 1 ] * [ -1 3 -3 1 ] * [ P-1 ]
+ [ 3 -6 3 0 ] [ P0 ]
+ [ -3 0 3 0 ] [ P1 ]
+ [ 1 4 1 0 ] [ P2 ]
+ @endverbatim
+
+ Where P-1 and P2 represent the neighbouring points or the
+ extrapolated end points. Simple reflection is used to
+ automatically create the end points.
+
+ The spline is discretized based on the chord length of the
+ individual segments. In rare cases (sections with very high
+ curvatures), the resulting distribution may be sub-optimal.
+
+SeeAlso
+ CatmullRomSpline
+
+ToDo
+ A future implementation could also handle closed splines - either
+ when the start/end points are identically or when specified.
+
+SourceFiles
+ BSpline2.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef BSpline2_H
+#define BSpline2_H
+
+#include "polyLine.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+ Class BSpline2 Declaration
+\*---------------------------------------------------------------------------*/
+
+class BSpline2
+:
+ public polyLine
+{
+ // Private Member Functions
+
+ //- Disallow default bitwise copy construct
+ BSpline2(const BSpline2&);
+
+ //- Disallow default bitwise assignment
+ void operator=(const BSpline2&);
+
+
+public:
+
+ // Constructors
+
+ //- Construct from components
+ BSpline2
+ (
+ const pointField& knots,
+ const vector& begTangentNotImplemented = vector::zero,
+ const vector& endTangentNotImplemented = vector::zero
+ );
+
+
+ // Member Functions
+
+ //- Return the point position corresponding to the curve parameter
+ // 0 <= lambda <= 1
+ point position(const scalar lambda) const;
+
+ //- Return the point position corresponding to the local parameter
+ // 0 <= lambda <= 1 on the given segment
+ point position(const label segment, const scalar lambda) const;
+
+ //- Return the length of the curve
+ scalar length() const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/mesh/blockMesh/curvedEdges/CatmullRomSpline.C b/src/mesh/blockMesh/curvedEdges/CatmullRomSpline.C
index d1647de07f7..74dbe69f85a 100644
--- a/src/mesh/blockMesh/curvedEdges/CatmullRomSpline.C
+++ b/src/mesh/blockMesh/curvedEdges/CatmullRomSpline.C
@@ -106,7 +106,7 @@ Foam::point Foam::CatmullRomSpline::position
return 0.5 *
(
- ( 2 * p0 )
+ ( 2*p0 )
+ mu *
(
( -e0 + p1 )
--
GitLab