FreeFOAM The Cross-Platform CFD Toolkit
BSpline.H
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2009-2010 OpenCFD Ltd.
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
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.
15 
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
19  for more details.
20 
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/>.
23 
24 Class
25  Foam::BSpline
26 
27 Description
28  An implementation of B-splines.
29 
30  In this implementation, the end tangents are created automatically
31  by reflection.
32 
33  In matrix form, the @e local interpolation on the interval t=[0..1] is
34  described as follows:
35  @verbatim
36  P(t) = 1/6 * [ t^3 t^2 t 1 ] * [ -1 3 -3 1 ] * [ P-1 ]
37  [ 3 -6 3 0 ] [ P0 ]
38  [ -3 0 3 0 ] [ P1 ]
39  [ 1 4 1 0 ] [ P2 ]
40  @endverbatim
41 
42  Where P-1 and P2 represent the neighbouring points or the extrapolated
43  end points. Simple reflection is used to automatically create the end
44  points.
45 
46  The spline is discretized based on the chord length of the individual
47  segments. In rare cases (sections with very high curvatures), the
48  resulting distribution may be sub-optimal.
49 
50 SeeAlso
51  CatmullRomSpline
52 
53 ToDo
54  A future implementation could also handle closed splines.
55 
56 SourceFiles
57  BSpline.C
58 
59 \*---------------------------------------------------------------------------*/
60 
61 #ifndef BSpline_H
62 #define BSpline_H
63 
64 #include "polyLine.H"
65 
66 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
67 
68 namespace Foam
69 {
70 
71 /*---------------------------------------------------------------------------*\
72  Class BSpline Declaration
73 \*---------------------------------------------------------------------------*/
74 
75 class BSpline
76 :
77  public polyLine
78 {
79  // Private Member Functions
80 
81  //- Disallow default bitwise copy construct
82  BSpline(const BSpline&);
83 
84  //- Disallow default bitwise assignment
85  void operator=(const BSpline&);
86 
87 
88 public:
89 
90  // Constructors
91 
92  //- Construct from components
93  BSpline
94  (
95  const pointField& knots,
96  const bool notImplementedClosed = false
97  );
98 
99 
100  // Member Functions
101 
102  //- Return the point position corresponding to the curve parameter
103  // 0 <= lambda <= 1
104  point position(const scalar lambda) const;
105 
106  //- Return the point position corresponding to the local parameter
107  // 0 <= lambda <= 1 on the given segment
108  point position(const label segment, const scalar lambda) const;
109 
110  //- Return the length of the curve
111  scalar length() const;
112 };
113 
114 
115 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
116 
117 } // End namespace Foam
118 
119 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
120 
121 #endif
122 
123 // ************************ vim: set sw=4 sts=4 et: ************************ //