Interpolation von Durchflussdaten

Numerische Berechnungen einer Füllstandsfunktion durch Interpolation, Ausgleichsrechnung, Min-Max-Approximation sowie Splining

Prof. Dr. Bernhard Zgraggen, HSR

Einführung in dieses Notebook

Dieses Dokument ist ein "selbsterklärendes" Mathematica-Notebook, das interaktiv bearbeitet werden kann ("Hands-On"). Input-Zellen können mittels der Tastenkombination ⇧⌤ (in Worten: Shift-Enter) ausgeführt werden. Falls Textzellen geschlossen sind, so können sie durch Klick auf das Icon ▼ oder auf die blaue Zellklammer rechterhand geöffnet werden.

Wenn Sie in Mathematica auf Befehle oder Begriffe stossen, die Sie nicht verstehen, so steht Ihnen die mächtige Hilfe/Doc Center zur Verfügung; diese ist interaktiv.

Allfällige Aufgaben/Übungen befinden sich am Ende des Dokuments.

Datenimport sowie Datenaufbereitung und -darstellung

Die Daten der Durchflussdaten sind in einer XLS-Tabelle (XLS) gegeben (InputData_Anschnitte.xls). Sie werden importiert und danach die einzelnen Durchflussfunktionen tabelliert. Die Ausgaben können mittels Semikolon unterdrückt werden..

SetDirectory[NotebookDirectory[]]

Die Durchflussdaten werden importiert und danach die Spalten in Listen von {x, y}-Koordinatenpaaren transformiert (Leerwerte "" werden gelöscht).

durchflussdata=Import["InputData_AnschnitteRed.xls", "XLS"];

(* Alternative: Chebyshev Messargumente: *) 
(* durchflussdata=Import["InputData_AnschnitteChebyshev.xls", "XLS"]; *)

Do[durchfluss[i]=DeleteCases[Extract[Transpose[durchflussdata[[1]]], {{1},{i+1}}]//Transpose,{_,""}], {i,1,5}]

Im folgenden werden Listenplots erstellt.

Do[
durchflussplot[i]=ListPlot[durchfluss[i], AxesLabel->{"Füllzeit sec", "Durchfluss (cm^3)/sec"}, PlotLabel->"Durchfluss "<>ToString[i],PlotStyle->Directive[{Hue[i/5],PointSize[0.005]}]],
{i,1,5}
]

Die Listenplots können zum Beispiel gemeinsam angezeigt werden. Die negativen Funktionswerte (ab ca. 11 sec) waren in den Originaldaten (Excel) enthalten; sie wurden unverändert übernommen.

Show[Table[durchflussplot[i],{i,1,5}],PlotLabel->"Füllkurven Anschnitte",PlotRange->Automatic]

Die folgende Code-Zeile ist funktional programmiert und erzeugt summierte Daten. Das Funktional {#[[1]], Sum[#[[i]], {i, 2, 6}]}& wird mittels des Map-Operators /@ auf die einzelnen Bestandteile der Durchflussdaten angewendet.

durchflussTotal={#[[1]],Sum[#[[i]],{i,2,6}]}&/@durchflussdata[[1]];

Show[Table[durchflussplot[i],{i,1,5}],ListPlot[durchflussTotal,PlotStyle->PointSize[0.005]],PlotLabel->"Füllkurven Anschnitte",PlotRange->All];

durchflussplot[5]

Show[durchflussplot[5],ListPlot[durchfluss[5],Joined->True,InterpolationOrder->1]]

Polynom-Interpolationen (Collocation) → {InterpolatingPolynomial}

Das "Collocation Polynomial" kann mittels des mächtigen und flexiblen Operators  InterpolatingPolynomial  gefunden werden. Dieser Operator gibt ein Objekt zurück, das mannigfache Art und Weise in Mathematica weiter verarbeitet werden kann (Differentiation, Integration, Grafik usw.). Er erlaubt symbolische "Daten", mehrdimensionale Interpolation, Ableitungsvorgaben, fehlende Daten u.v.a.m.

Clear[interpoly]
Do[interpoly[i]=InterpolatingPolynomial[durchfluss[i],x],{i,1,5}]

Ein Interpolationspolynom wird in Horner-Form (s. HornerForm) gespeichert (zwecks effizienter und numerisch eher stabiler Auswertung). Die mit Expand entwickelte Form zeigt die Standardform des Polynoms (nur informativ). Das Strichpunkt verhindert die Ausgabe (spart Speicher und Platz bei grossen Ausgaben).

interpoly[5];

interpoly[5]//Expand;

Der folgende kombinierte Plot zeigt den eingeschränkten Nutzen der Polynom-Interpolation bei relativ grossen Datenmengen. Das Schwingverhalten ist stark und deshalb schwankt der Fehler in einer grossen Bandbreite. Die Polynom-Interpolation führt zu hochgradigen Polynomen (diese divergieren absolut immer gegen ∞, wenn x → ±∞, wodurch sie auch für Extrapolationen ungeeignet sind).

i=5;
Show[{Plot[interpoly[i],{x,1.7818,11.14},PlotRange->{0,600}],durchflussplot[i]}]

interpoly=InterpolatingPolynomial[durchflussTotal,x];
interpoly//Expand;

In ähnlicher Art werden auch die summierten Durchflussdaten interpoliert und die Resultate grafisch dargestellt.

Show[Table[durchflussplot[i],{i,1,5}],ListPlot[durchflussTotal,PlotStyle->PointSize[0.005]],Plot[interpoly,{x,1.7818,11.14}],PlotLabel->"Füllkurven Anschnitte und Füllstand",PlotRange->All,AxesLabel->{"Füllzeit sec", "(cm^3)/sec bzw. (cm^3)"}]

*Integrationsversuch mit Polynom-Interpolation (Gesamtfüllstand)

Das Integral einer Durchflussfunktion ergibt theoretisch die Füllstandsfunktion (bis auf additive Verschiebung: Zum Zeitpunkt x = 1.7818 ist der Füllstand jeweils 0). Der Operator /. ist der in Mathematica wichtige Ersetzungsoperator (ReplaceAll). Das Verwenden der HornerForm verbessert die numerische Auswertung hochgradiger Polynome signifikant.

Do[fillpoly1[i]=HornerForm@Integrate[interpoly[i],x],{i,1,5}];
Do[fillpoly[i]=fillpoly1[i]-(fillpoly1[i]/.{x->1.7818}),{i,1,5}];  
(* Normierung so, dass der Wert bei x = 1.7818 Null wird. *)

fill=Sum[fillpoly[i],{i,1,5}];

Show[Table[durchflussplot[i],{i,1,5}],Plot[fill,{x,0,12}],PlotLabel->"Füllkurven Anschnitte und Füllstand",PlotRange->Automatic,AxesLabel->{"Füllzeit sec", "(cm^3)/sec bzw. (cm^3)"}]

*Datenreduktion (Provisorisch)

Durch folgende Code-Zeilen werden Daten generisch ausgespart (reduziert). Natürlich wäre es besser, sorgfältiger Daten auszusparen (bspw. Ausreisser mit erster Priorität zu eliminieren) bzw. die Daten auszuglätten. Weil es Approximationsverfahren gibt, die das Schwingproblem relativ effizient und elegant lösen werden, wird hier auf eine Ausarbeitung der Datenreduktion verzichtet.

step={32,16,16,16,16}; (* Schrittweiten *)
Do[durchflussRed[i]=Take[durchfluss[i],{1,Length[durchfluss[i]],step[[i]]}],{i,1,5}]

Die Interpolationsfehler an ausgesparten Datenpunkten werden kompensiert durch ein verbessertes "Gesamtbild". Die Frage ist, was das heisst, und mit welchem Mass die Abweichung gemessen wird.

Clear[interpoly];
Do[interpoly[i]=InterpolatingPolynomial[durchflussRed[i],x],{i,1,5}]

i=1;
Show[{Plot[interpoly[i],{x,0,12},PlotRange->{-200,600}],durchflussplot[i]}]

Die numerische Berechnung einer Füllstandsfunktion durch Integration und Summierung des letzten Abschnitts kann nun erneut und mit leicht besserem Resultat ausgeführt werden.

Splines → {Interpolation}

Piecewise Spline-Interpolation und Integration, Ableitungen

Spline of order 3 (cubic spline) with appropriate settings for the options InterpolationOrder and Method. Mathematica uses a generic spline method called B-Splines (Basis splines) which does not consider constraints at the boundaries (like the natural spline). B-Splines use a cleverly devised set of basis functions in order to simplify the computation of the spline formulas and in order to improve the algorithmic treatment.

Do[splinepoly[i]=Interpolation[durchfluss[i],InterpolationOrder->3],{i,1,5}]

The results are of the type InterpolationFunction. This is a functional object fully integrated into Mathematicas operations with functions. Postfix //InputForm if internal structure of such an object has to be considered.

Integrate[splinepoly[5][x],x]

splinepoly[5]//InputForm;

i=5;
Show[{Plot[{splinepoly[i][x]},{x,1.7818,11.7478},PlotRange->All],durchflussplot[i]}]

The accumulated flow datas are approximated by integration of the piecewise interpolating functions:

Integrate[splinepoly[i][x],x]

Do[
splinepolyacc[i]=(Integrate[splinepoly[i][x],x]);
splinepolydiff[i]=(D[splinepoly[i][x],{x,1}]),
{i,1,5}
]

Evaluation of the interpolating functions outside the range of data yields warnings.

i=5;
Show[Plot[{splinepolyacc[i]},{x,1.7818,11.14},PlotStyle->Hue[i/5]]]
Show[Plot[{splinepolydiff[i]},{x,1.7818,11.14},PlotStyle->Hue[i/5]]]

Plot[Evaluate@Table[splinepolyacc[i],{i,1,5}],{x,1.7818,11.7478},PlotStyle->Table[Hue[i/5],{i,1,5}],PlotRange->Automatic]

The next  GraphicsGrid shows all flow data together with accumulated (integrated) and differentiated data, respectively.

GraphicsGrid[Table[{Show[durchflussplot[i],Plot[splinepoly[i][x],{x,1.7818,11.7478},PlotStyle->{Gray,Thin}]],
Plot[splinepolydiff[i],{x,1.7818,11.7478},PlotStyle->Hue[i/5],PlotLabel->"Beschleunigung "<>ToString[i],PlotRange->Automatic,AxesLabel->{"Zeit sec", "(cm^3)/(sec^2)"}],Plot[splinepolyacc[i],{x,1.7818,11.7478},PlotStyle->Hue[i/5],PlotLabel->"Füllstand "<>ToString[i],PlotRange->Automatic,AxesLabel->{"Zeit sec", "(cm^3)"}]},{i,1,5}],ImageSize->800]

Die Tridiagonalmatrix (natural spline)

The symmetric Tridiagonalmatrix for the coefficients c[i] of the quadratic terms  (((x - x[i]))^2)   (i = 1, 2 … N-1) mit c[0]==0 and c[N]==0 is the following one.

i=1;
h[i]=Table[durchfluss[i][[k+1,1]]-durchfluss[i][[k,1]],{k,1,Length[durchfluss[i]]-1}];
Manipulate[
SparseArray[{{k_,k_}:>2( h[i][[k]]+h[i][[k+1]]),{k_,l_}/;Abs[k-l]==1:>h[i][[Max[{k+1,l+1}]]]},{n,n}]//MatrixForm,
{{n,3,"n"},1,Length[durchfluss[i]]-2,1,Appearance->"Labeled"}
]

*Ausgleichsspline und Integration (Experimental)

i=1;
x[i]=Table[durchfluss[i][[k,1]],{k,1,Length[durchfluss[i]]}];
y[i]=Table[durchfluss[i][[k,2]],{k,1,Length[durchfluss[i]]}];
n=Length[durchfluss[i]];
n=10

constraints=Flatten@Table[
{((a[k-1]+b[k-1](x-x[i][[k-1]])+c[k-1](((x-x[i][[k-1]]))^2)+d[k-1](((x-x[i][[k-1]]))^3))/.{x->x[i][[k]]})==((a[k]+b[k](x-x[i][[k]])+c[k](((x-x[i][[k]]))^2)+d[k](((x-x[i][[k]]))^3))/.{x->x[i][[k]]}),
((D[a[k-1]+b[k-1](x-x[i][[k-1]])+c[k-1](((x-x[i][[k-1]]))^2)+d[k-1](((x-x[i][[k-1]]))^3),x])/.{x->x[i][[k]]})==((D[a[k]+b[k](x-x[i][[k]])+c[k](((x-x[i][[k]]))^2)+d[k](((x-x[i][[k]]))^3),x])/.{x->x[i][[k]]}),
((D[a[k-1]+b[k-1](x-x[i][[k-1]])+c[k-1](((x-x[i][[k-1]]))^2)+d[k-1](((x-x[i][[k-1]]))^3),{x,2}])/.{x->x[i][[k]]})==((D[a[k]+b[k](x-x[i][[k]])+c[k](((x-x[i][[k]]))^2)+d[k](((x-x[i][[k]]))^3),{x,2}])/.{x->x[i][[k]]})
},{k,2,n-1}];

objective[lambda_]:=(∑_(k=1)^(n-1))(((y[i][[k]]-(a[k]+b[k](x-x[i][[k]])+c[k](((x-x[i][[k]]))^2)+d[k](((x-x[i][[k]]))^3))/.{x->x[i][[k]]}))^2)+(((y[i][[n]]-(a[n-1]+b[n-1](x-x[i][[n-1]])+c[n-1](((x-x[i][[n-1]]))^2)+d[n-1](((x-x[i][[n-1]]))^3))/.{x->x[i][[n]]}))^2)+
lambda*(∑_(k=1)^(n-1))((((((2c[k]+6d[k](((x-x[i][[k]]))^1)))^3)/(3*6d[k]))/.{x->x[i][[k+1]]})-(((((2c[k]+6d[k](((x-x[i][[k]]))^1)))^3)/(3*6d[k]))/.{x->x[i][[k]]}))

vars=Flatten@Table[{a[k],b[k],c[k],d[k]},{k,1,n-1}];

Manipulate[Module[{min,pw},
min=NMinimize[{objective[lambda],constraints},vars];
pw=Piecewise[Table[{(a[k]+b[k](x-x[i][[k]])+c[k](((x-x[i][[k]]))^2)+d[k](((x-x[i][[k]]))^3))/.Partition[Last[min],4][[k]],(x[i][[k]]≤x&&x≤x[i][[k+1]])},{k,1,n-1}]]
],
{{lambda,1,"λ"},0,100,Appearance->"Labeled"}]

Lineare Ausgleichsrechnung (Linear Least-Squares) → {Fit, FindFit}

Basis polynomials {1, ((x-(Subscript[μ, i]))/(Subscript[σ, i])), (((((x-(Subscript[μ, i]))/(Subscript[σ, i]))))^2), ⋯ , (((((x-(Subscript[μ, i]))/(Subscript[σ, i]))))^n)}. The parameters (Subscript[μ, i]) and (Subscript[σ, i]), respectively, denote the empirical mean and standard deviation of the x-data sets. Using these transformations improves numerical stability when using higher polynomial degrees.

Anwendung von Fit und Integration bzw. Ableitung

A least-squares polynomial or another best fitting curve may be found by Mathematicas command Fit or FindFit (the latter is appropriate for non-linear models). The next commands show how to compute least-squares polynomials of maximal degree n in the shifted and scaled basis polynomials {1, ((x-(Subscript[μ, i]))/(Subscript[σ, i])), (((((x-(Subscript[μ, i]))/(Subscript[σ, i]))))^2), ⋯ , (((((x-(Subscript[μ, i]))/(Subscript[σ, i]))))^n)}. The parameters (Subscript[μ, i]) and (Subscript[σ, i]), respectively, denote the empirical mean and standard deviation of the x-data sets. Using these transformations improves numerical stability when using higher polynomial degrees.

polybasis1[n_,i_]:=Table[(((((x-Mean[#[[1]]&/@durchfluss[i]])/(StandardDeviation[#[[1]]&/@durchfluss[i]]))))^j),{j,0,n}]

polybasis1[8,5]

degree=4;
Do[leastsquarepoly[i]=Fit[durchfluss[i],polybasis1[degree,i],x],{i,1,5}]

leastsquarepoly[5]

i=5;
Show[{Plot[HornerForm@leastsquarepoly[i],{x,0,12},PlotRange->{-100,700}],durchflussplot[i]}]

The accumulated flow datas are approximated by integration of least-squares polynomials. The results are normalized and stored in HornerForm.

Do[
leastsquarepolyacc[i]=HornerForm@(Integrate[leastsquarepoly[i],x]);
leastsquarepolyacc[i]=leastsquarepolyacc[i]-(leastsquarepolyacc[i]/.{x->1.7818});
leastsquarepolydiff[i]=HornerForm@(D[leastsquarepoly[i],x]),
{i,1,5}]

i=5;
Show[Plot[leastsquarepolyacc[i],{x,0,12},PlotStyle->Hue[i/5]]]

Plot[Evaluate@Table[leastsquarepolyacc[i],{i,1,5}],{x,0,12},PlotStyle->Table[Hue[i/5],{i,1,5}],PlotRange->Automatic]

The next  GraphicsGrid shows all flow data together with accumulated (integrated) and differentiated data.

GraphicsGrid[Table[{Show[durchflussplot[i],Plot[leastsquarepoly[i],{x,0,12},PlotStyle->{Gray,Thin}]],Plot[leastsquarepolydiff[i],{x,0,12},PlotStyle->Hue[i/5],PlotLabel->"Beschleunigung "<>ToString[i],PlotRange->Automatic,AxesLabel->{"Zeit sec", "(cm^3)/(sec^2)"}],
Plot[leastsquarepolyacc[i],{x,0,12},PlotStyle->Hue[i/5],PlotLabel->"Füllstand "<>ToString[i],PlotRange->Automatic,AxesLabel->{"Zeit sec", "(cm^3)"}]},{i,1,5}],ImageSize->800]

DesignMatrix und Normalgleichungen → {DesignMatrix, LinearSolve}

The DesignMatrix for the numerical  data and the polynomial basis:

{a,b}={1.7818,11.14};

polybasis[n_,i_]:=Table[(((x))^j),{j,0,n}]

polybasis[16,5]

Shifted and scaled basis polynomials {1, ((x-(Subscript[μ, i]))/(Subscript[σ, i])), (((((x-(Subscript[μ, i]))/(Subscript[σ, i]))))^2), ⋯ , (((((x-(Subscript[μ, i]))/(Subscript[σ, i]))))^n)}. The parameters (Subscript[μ, i]) and (Subscript[σ, i]), respectively, denote the empirical mean and standard deviation of the x-data sets. Using these transformations improves numerical stability when using higher polynomial degrees.

polybasis1[n_,i_]:=Table[(((((x-Mean[#[[1]]&/@durchfluss[i]])/(StandardDeviation[#[[1]]&/@durchfluss[i]]))))^j),{j,0,n}]

polybasis1[16,5]

degree=4; (* Maximalgrad der Basispolynome *)
m=5; (* Durchflussnr. *)
dm=DesignMatrix[durchfluss[m],polybasis[degree,m],x,IncludeConstantBasis->False];
Dimensions[dm]
dm//MatrixForm;

The normal matrix ...

Transpose[dm].dm//MatrixForm

y=Table[durchfluss[m][[i,2]],{i,1,Length[durchfluss[m]]}];

as=LinearSolve[Transpose[dm].dm,Transpose[dm].y]

as.polybasis[degree,m]

Show[{Plot[as.polybasis[degree,m],{x,0,12},PlotRange->{-100,700}],durchflussplot[m]}]

degree=4;
Do[leastsquarespoly[i]=Fit[durchfluss[i],polybasis[degree,i],x],{i,1,5}]

leastsquarespoly[5]//Simplify

Show[{Plot[leastsquarespoly[m],{x,a,b},PlotRange->{-100,700}],durchflussplot[m]}]

The accumulated flow datas are approximated by integration of least-squares polynomials. The results are normalized and stored in HornerForm.

Do[
leastsquarespolyacc[i]=HornerForm@(Integrate[leastsquarespoly[i],x]);
leastsquarespolyacc[i]=leastsquarespolyacc[i]-(leastsquarespolyacc[i]/.{x->1.7818});
leastsquarespolydiff[i]=HornerForm@(D[leastsquarespoly[i],x]),
{i,1,5}]

i=5;
Show[Plot[leastsquarespolyacc[i],{x,a,b},PlotStyle->Hue[i/5]]]

Plot[Evaluate@Table[leastsquarespolyacc[i],{i,1,5}],{x,a,b},PlotStyle->Table[Hue[i/5],{i,1,5}],PlotRange->Automatic]

The next GraphicsGrid shows all flow data together with accumulated (integrated) and differentiated data.

GraphicsGrid[Table[{Show[durchflussplot[i],Plot[leastsquarespoly[i],{x,a,b},PlotStyle->{Gray,Thin}]],Plot[leastsquarespolydiff[i],{x,a,b},PlotStyle->Hue[i/5],PlotLabel->"Beschleunigung "<>ToString[i],PlotRange->Automatic,AxesLabel->{"Zeit sec", "(cm^3)/(sec^2)"}],
Plot[leastsquarespolyacc[i],{x,a,b},PlotStyle->Hue[i/5],PlotLabel->"Füllstand "<>ToString[i],PlotRange->Automatic,AxesLabel->{"Zeit sec", "(cm^3)"}]},{i,1,5}],ImageSize->800]

DesignMatrix, QR-Zerlegung → {DesignMatrix, QRDecomposition}

The DesignMatrix for the numerical data and the corresponding normalized basis polynomials:

polybasis1[n_,i_]:=Table[(((((x-Mean[#[[1]]&/@durchfluss[i]])/(StandardDeviation[#[[1]]&/@durchfluss[i]]))))^j),{j,0,n}]

degree=16; (* Maximalgrad der Basispolynome *)
m=5; (* Durchflussnr. *)
dm=DesignMatrix[durchfluss[m],polybasis1[degree,m],x,IncludeConstantBasis->False];
%//MatrixForm;

The normal matrix ...

Transpose[dm].dm//MatrixForm

... and the design matrix followed by QRDecomposition :

{q,r}=QRDecomposition[dm];

In Mathematica q is row-orthogonal and r is upper right square triangle. With these matrices we can solve for the coefficients of the least-squares approximating polynomial. The coefficients refer to the basis.

y=Table[durchfluss[m][[i,2]],{i,1,Length[durchfluss[m]]}];

as=LinearSolve[r,q.y]

as.polybasis1[degree,m]

The results by QRDecomposition do not differ from those of Mathematicas Fit Command.

degree=16;
leastsquaresbasispoly=Fit[durchfluss[m],polybasis1[degree,m],x]

DesignMatrix, SVD-Zerlegung → {DesignMatrix, SingularValueDecomposition}

The DesignMatrix for the numerical data and the corresponding Newton-Basis:

polybasis[n_,i_]:=Table[(((x))^j),{j,0,n}]

polybasis1[n_,i_]:=Table[(((((x-Mean[#[[1]]&/@durchfluss[i]])/(StandardDeviation[#[[1]]&/@durchfluss[i]]))))^j),{j,0,n}]

degree=4;
polybasis1[degree,5]

(* Maximalgrad der Basispolynome *)
m=5; (* Durchflussnr. *)
dm=DesignMatrix[durchfluss[m],polybasis[degree,m],x,IncludeConstantBasis->False];
%//MatrixForm;

The normal matrix ...

Transpose[dm].dm//MatrixForm

... and the design matrix followed by SingularValueDecomposition :

{u,d,v}=SingularValueDecomposition[dm,degree+1];

The matrices u and v are orthogonal, d is diagonal. 
With these matrices we can solve for the coefficients of the least-squares approximating polynomial. The coefficients refer to the Newton-basis.

MatrixForm/@{d,v}

Dimensions/@{u,d,v}

y=Table[durchfluss[m][[i,2]],{i,1,Length[durchfluss[m]]}];

as=v.((Transpose[u].y)/Diagonal[d])

as.polybasis1[degree,m]//Expand

Show[{Plot[as.polybasis1[degree,m],{x,a,b},PlotRange->{-100,700}],durchflussplot[m]}]

degree=4;
leastsquarespoly1=Fit[durchfluss[m],polybasis1[degree,m],x]//Expand

Show[{Plot[leastsquarespoly1,{x,a,b},PlotRange->{-100,700}],durchflussplot[m]}]

Chebyshev Least-Squares → {DesignMatrix, ChebyshevT}

The DesignMatrix for the numerical Chebyshev data and the corresponding Chebyshev-Basis:

{a,b}={1.7818,11.14};

chebbasis[k_,x_]:=Table[ChebyshevT[j,-1+2((x-a)/(b-a))],{j,0,k}]

chebbasis[4,x]

degree=4; (* Maximalgrad der Basispolynome *)
m=5; (* Durchflussnr. *)
dm=DesignMatrix[durchfluss[m],chebbasis[degree,x],x,IncludeConstantBasis->False];
%//MatrixForm;

The normal matrix ...

Transpose[dm].dm//Chop//MatrixForm

y=Table[durchfluss[m][[i,2]],{i,1,Length[durchfluss[m]]}];

as=(Transpose[dm].y)/Diagonal[Transpose[dm].dm]

as.chebbasis[degree,x]

degree=4;
leastsquareschebpoly[m]=Fit[durchfluss[m],chebbasis[degree,x],x]

m=5;
Show[{Plot[HornerForm@leastsquareschebpoly[m],{x,0,12},PlotRange->{-100,700}],durchflussplot[m]}]

*Unstabile Orthogonalisierung mit modifiziertem Gram-Schmidt-Verfahren

The DesignMatrix for the numerical data and the corresponding Newton-Basis:

m=5; (* Durchflussnr. *)
newtonpoly[k_,var_]:=Product[var-x[i],{i,0,k-1}]/.Table[x[i]->durchfluss[m][[i+1,1]],{i,0,k}]

newtonpoly[16,x]

To avoid normalization a modified version of the Gram-Schmidt procedure only using Projections is defined:

myGramSchmidt[vecs_,ip___]:=Module[{ovecs=vecs},
Do[ovecs[[i]]-=Projection[ovecs[[i]],ovecs[[j]],ip],{i,2,Length[vecs]},{j,1,i-1}];
ovecs]

m=5; (* Durchflussnr. *)
degree=8;
orthonewton2=myGramSchmidt[Table[newtonpoly[k,x],{k,0,degree}],Sum[(#1/.{x->durchfluss[m][[i]][[1]]})* ( #2/.{x->durchfluss[m][[i]][[1]]}),{i,1,Length[durchfluss[m]]}]&];

orthonewton2;

dm=DesignMatrix[durchfluss[m],HornerForm/@orthonewton2,x,IncludeConstantBasis->False];

The modified Gram-Schmidt procedure is numerically not stable. The matrix of the normal equations is not diagonal.

Transpose[dm].dm//MatrixForm;

The next commands show better numeric stability for the shifted and scaled basis polynomials {1, ((x-(Subscript[μ, i]))/(Subscript[σ, i])), (((((x-(Subscript[μ, i]))/(Subscript[σ, i]))))^2), ⋯ , (((((x-(Subscript[μ, i]))/(Subscript[σ, i]))))^n)} together with Mathematicas Orthogonalize and HornerForm. The parameters (Subscript[μ, i]) and (Subscript[σ, i]), respectively, denote the empirical mean and standard deviation of the x-data sets. Using these transformations improves numerical stability when using higher polynomial degrees. Mathematicas internal command Orthogonalize is more stable. The most stable and fastest Householder method is not applicable for non-numeric data.

polybasis1[n_,i_]:=Table[(((((x-Mean[#[[1]]&/@durchfluss[i]])/(StandardDeviation[#[[1]]&/@durchfluss[i]]))))^j),{j,0,n}]

m=5; (* Durchflussnr. *)
degree=8;
ortho1=myGramSchmidt[HornerForm/@polybasis1[degree,m],Sum[(#1/.{x->durchfluss[m][[i]][[1]]})* ( #2/.{x->durchfluss[m][[i]][[1]]}),{i,1,Length[durchfluss[m]]}]&];

dm=DesignMatrix[durchfluss[m],HornerForm/@ortho1,x,IncludeConstantBasis->False];

Transpose[dm].dm//MatrixForm

m=5; (* Durchflussnr. *)
degree=8;
ortho2=Orthogonalize[HornerForm/@polybasis1[degree,m],Sum[(#1/.{x->durchfluss[m][[i]][[1]]})* ( #2/.{x->durchfluss[m][[i]][[1]]}),{i,1,Length[durchfluss[m]]}]&];

dm=DesignMatrix[durchfluss[m],HornerForm/@ortho2,x,IncludeConstantBasis->False];

Transpose[dm].dm//MatrixForm;