> restart:
> with(plots):with(LinearAlgebra):
Warning, the name changecoords has been redefined
Aproksymacja
Aproksymacja funkcji dyskretnej
> n:=9;
n := 9
> X:=array(0..n,[0,1,2,3,4,5,6,7,8,9]):
> Y:=array(0..n,[-1.2,-8,-22.1,-37.6,-49.5,-58.7,-54.0,-39.3,-3.5,
54.1]):
> plot(zip((x,y)->[x,y],X,Y),style=point,symbol=cross,symbolsize=3
0);
> pkt:=%:
>
Generowanie i rozwiazanie ukladu równan na wspolczynniki wielomianu aproksymujacego a[j]
n m
S = yk - xkj ]2
"[ "a
j
k =0 j=0
> m:=3; # Digits:=20; # dla m=9
m := 3
> S:=add((Y[k]-add(a[j]*X[k]^j,j=0..m))^2,k=0..n);
2 2 2
S := ( -1.2 - a0) + ( -8 - a0 - a1 - a2 - a3 ) + ( -22.1 - a0 - 2 a1 - 4 a2 - 8 a3)
2 2
+ ( -37.6 - a0 - 3 a1 - 9 a2 - 27 a3) + (-49.5 - a0 - 4 a1 - 16 a2 - 64 a3 )
2 2
+ ( -58.7 - a0 - 5 a1 - 25 a2 - 125 a3) + (-54.0 - a0 - 6 a1 - 36 a2 - 216 a3 )
2 2
+ ( -39.3 - a0 - 7 a1 - 49 a2 - 343 a3) + (-3.5 - a0 - 8 a1 - 64 a2 - 512 a3 )
2
+ ( 54.1 - a0 - 9 a1 - 81 a2 - 729 a3)
> for i from 0 to m do
 r||i:=diff(S,a[i])=0;
end do;
r0 := 20 a0 + 90 a1 + 570 a2 + 4050 a3 + 439.6 = 0
r1 := 90 a0 + 570 a1 + 4050 a2 + 30666 a3 + 1593.4 = 0
r2 := 570 a0 + 4050 a1 + 30666 a2 + 241650 a3 + 4811.8 = 0
r3 := 4050 a0 + 30666 a1 + 241650 a2 + 1956810 a3 - 1595.0 = 0
> roz:=solve({r||(0..m)},{seq(a[i],i=0..m)});
roz := { a2 = -5.098892774, a1 = -2.753399378, a0 = -1.175104895, a3 = 0.6760683761 }
> assign(roz):
> g1:=add(a[j]*x^j,j=0..m); # wielomian aproksymujÄ…cy
g1 := -1.175104895 - 2.753399378 x - 5.098892774 x2 + 0.6760683761 x3
> p:=plot(g1,x=0..n):
> display([pkt,p]);
> S;
4.588004950
S
µ =
b - a
> epsilon:=evalf(sqrt(S/(X[n]-X[0])));
µ := 0.7139876244
> a:='a':
Zapis macierzowy ukladu równan
nn
i+ j
Mi+1, j+1 = , bi+1 = yk xki
"x "
k
k =0 k =0
>
> M:=Matrix(m+1):b:=Vector(m+1):
> for i from 0 to m do
b[i+1]:=add(Y[k]*X[k]^i,k=0..n):
 for j from 0 to m do
M[i+1,j+1]:=add(X[k]^(i+j),k=0..n):
end do:
end do:
> M,b;
îÅ‚ 10 45 285 2025Å‚Å‚ îÅ‚ -219.8 Å‚Å‚
ïÅ‚ śł ïÅ‚ śł
ïÅ‚ śł ïÅ‚ śł
45 285 2025 15333śł ïÅ‚ -796.7
ïÅ‚ śł
ïÅ‚ śł, ïÅ‚-2405.9śł
śł
285 2025 15333 120825śł ïÅ‚
ïÅ‚
ïÅ‚ śł ïÅ‚ śł
ðÅ‚2025 15333 120825 978405ûÅ‚ ðÅ‚ 797.5 ûÅ‚
> a:=LinearSolve(M,b);
îÅ‚-1.17510489510479842Å‚Å‚
ïÅ‚ śł
ïÅ‚-2.75339937839944104śł
śł
ïÅ‚
a := ïÅ‚
śł
ïÅ‚-5.09889277389276428śł
ïÅ‚ śł
ðÅ‚0.676068376068375664ûÅ‚
> g2:=add(a[j+1]*x^j,j=0..m); # wielomian aproksymujÄ…cy wyznaczony
z układu równan A.a=b
g2 := -1.175104895 - 2.75339937839944104 x - 5.09889277389276428 x2
+ 0.676068376068375664 x3
> 'g1'=g1;
g1 = -1.175104895 - 2.753399378 x - 5.098892774 x2 + 0.6760683761 x3
> a:='a':b:='b':
Komenda w Maple'u
> with(CurveFitting);
Warning, the name LeastSquares has been rebound
[BSpline, BSplineCurve, Interactive, LeastSquares, PolynomialInterpolation,
RationalInterpolation, Spline, ThieleInterpolation ]
> type(X,list);
false
> X:=convert(X,list); Y:=convert(Y,list);
X := [0, 1, 2, 3, 4, 5, 6, 7, 8, 9 ]
Y := [-1.2, -8, -22.1, -37.6, -49.5, -58.7, -54.0, -39.3, -3.5, 54.1 ]
> type(X,list);
true
> gM:=LeastSquares(X,Y,x,curve=add(a[i]*x^i,i=0..m));
gM := -1.175104895 - 2.75339937839955874 x - 5.09889277389273054 x2
+ 0.676068376068373222 x3
> 'g2'=g2;
g2 = -1.175104895 - 2.75339937839944104 x - 5.09889277389276428 x2
+ 0.676068376068375664 x3
>
Aproksymacja funkcji okresowej
> f:=x->piecewise(x<-1,-2-x,x<1,x,2-x);
f := x piecewise( x < -1, -2 - x, x < 1, x, 2 - x )
> plot(f,-2..2,scaling=constrained);
> fokr:=piecewise(x<-3,x+4,x<-1,-2-x,x<1,x,x<3,2-x,-4+x); # kilka
okresow na wykresie
Å„Å‚ x + 4 x < -3
ôÅ‚
ôÅ‚ -2 - x x < -1
ôÅ‚
ôÅ‚
fokr := òÅ‚ x x < 1
ôÅ‚
ôÅ‚
2 - x x < 3
ôÅ‚
ôÅ‚
ół -4 + x otherwise
> plot(fokr,x=-5..5,scaling=constrained);
T
2
22Ä„ix
ai = f (x) cos dx, i = 0,1,...
+"
TT
-T
2
T
2
22Ä„ix
bi = f (x)sin dx, i = 1, 2,...
+"
TT
-T
2
> T:=4;
T := 4
> m:=10; # 5, 10, 15, 25, 50
m := 10
> for i from 0 to m do
a[i]:=2/T*evalf(Int(f(x)*cos(2*Pi/T*i*x),x=-T/2..T/2)):
b[i]:=2/T*evalf(Int(f(x)*sin(2*Pi/T*i*x),x=-T/2..T/2)):
end do:
m
1 2Ä„i 2Ä„i
g(x) = a0 + cos x + bi sin x)
"(a
2 i TT
i=1
> g:=1/2*a[0]+add(a[i]*cos(2*Pi/T*i*x)+b[i]*sin(2*Pi/T*i*x),i=1..m
);
ëÅ‚ Ä„ x öÅ‚ ëÅ‚ 3 Ä„ x öÅ‚ ëÅ‚ 5 Ä„ x öÅ‚
g := 0.8105694690 sinìÅ‚ ÷Å‚
ìÅ‚ ÷Å‚ - 0.09006327436 sinìÅ‚ ÷Å‚ ìÅ‚ ÷Å‚
ìÅ‚ ÷Å‚ + 0.03242277876 sinìÅ‚ ÷Å‚
íÅ‚ 2 Å‚Å‚ íÅ‚ 2 Å‚Å‚ íÅ‚ 2 Å‚Å‚
ëÅ‚ 7 Ä„ x öÅ‚ ëÅ‚ 9 Ä„ x öÅ‚
- 0.01654223406 sinìÅ‚ ÷Å‚ + 0.01000703048 sinìÅ‚ ÷Å‚ - 0.1000000000 10-13 sin( 5 Ä„ x )
ìÅ‚ ÷Å‚ ìÅ‚ ÷Å‚
íÅ‚ 2 Å‚Å‚ íÅ‚ 2 Å‚Å‚
> plot([fokr,g],x=-5..5,color=[red,blue],scaling=constrained,numpo
ints=200); # porównanie
>
>