Plot Lagrange parabolas for a Simpson's composite rule
显示 更早的评论
This is my code:
clearvars; close all; clc;
f =@ (x) 2 + sin((pi/20)*x) + sin((pi/3)*x);
x0 = input("Inferior limit in x to aproximate: ") ;
xm = input("Exterior limit in x to approximate: ");
n = input("Number of parabolas for the Simpson's composite rule: ");
m = 2*n;
h = (xm-x0)/(m);
X = x0:h:xm;
Y = f(X);
Zint = zeros(1,n);
Zext = zeros(1,n+1);
evenidx =@ (v) v(2:2:end);
oddidx =@ (v) v(1:2:end);
Xint = evenidx(X);
Xext = oddidx(X);
Yint = f(Xint);
Yext = f(Xext);
fplot(f, [x0,xm], "b") % Function f(x)
hold on
plot([Xint; Xint], [Zint; Yint], "g") % This marks the interior lines in green
plot([Xext; Xext], [Zext; Yext], "r") % This marks the exterior lines in red (the limits of the parabolas)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for i = 1:2:m
P =@ (x) ((x-X(i+1))*(x-X(i+2)))/((X(i)-X(i+1))*(X(i)-x(i+2)))*Y(i) + ...
((x-X(i))*(x-X(i+2)))/(X(i+1)-X(i))*(X(i+1)-X(i+2))*Y(i+1) + ... % Within here lies my problem
((x-X(i))*(x-X(i+1)))/(X(i+2)-X(i))*(X(i+2)-X(i+1))*Y(i+2);
flot(P, [X(i), X(i+2)], "r")
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
In this code I provide the function and the intervals I want in my Simpson's rule and I expect to plot the Simpson's rule graphically. My question comes from the last point. I was trying to find out an automated way to plot the parabolas between one point X(i) and the point X(i+2) while touching the point X(i+1). For that I tried to draw an fplot of P being the Lagrange polynomial formula inside a loop. As I can see, this cannot be done in Matlab, at least in this way.
There is any way to draw the parabolas betweeen the red lines? Thank you so much!
PD: The function f(x) is totally arbitrary and can be changed as one likes
采纳的回答
clearvars; close all; clc;
f =@ (x) 2 + sin((pi/20)*x) + sin((pi/3)*x);
x0 = 0;%input("Inferior limit in x to aproximate: ") ;
xm = 1;%input("Exterior limit in x to approximate: ");
n = 5;%input("Number of parabolas for the Simpson's composite rule: ");
m = 2*n;
h = (xm-x0)/(m);
X = x0:h:xm;
Y = f(X);
Zint = zeros(1,n);
Zext = zeros(1,n+1);
evenidx =@ (v) v(2:2:end);
oddidx =@ (v) v(1:2:end);
Xint = evenidx(X);
Xext = oddidx(X);
Yint = f(Xint);
Yext = f(Xext);
fplot(f, [x0,xm], "b") % Function f(x)
hold on
plot([Xint; Xint], [Zint; Yint], "g") % This marks the interior lines in green
plot([Xext; Xext], [Zext; Yext], "r") % This marks the exterior lines in red (the limits of the parabolas)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for i = 1:2:m
P =@ (x) ((x-X(i+1)).*(x-X(i+2)))/((X(i)-X(i+1))*(X(i)-X(i+2)))*Y(i) + ...
((x-X(i)).*(x-X(i+2)))/((X(i+1)-X(i))*(X(i+1)-X(i+2)))*Y(i+1) + ... % Within here lies my problem
((x-X(i)).*(x-X(i+1)))/((X(i+2)-X(i))*(X(i+2)-X(i+1)))*Y(i+2);
fplot(P, [X(i)-2, X(i+2)+2], "r")
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
hold off
更多回答（0 个）
类别
在 帮助中心 和 File Exchange 中查找有关 Eigenvalue Problems 的更多信息
Translated by 
