필터 지우기
필터 지우기

My ellipse wont close after running the code .

조회 수: 2 (최근 30일)
Panda05
Panda05 2023년 12월 2일
댓글: Panda05 2023년 12월 3일
Hy guys . In this task , i have to implement periodic boundary conditons to my ellipse . However , the boundary conditions are not being implemented and i dont know what mistake i'm making if you comment them, the whole code still runs and you can see that the ellipse wont close at the end because the boundary conditions are not working at all. Any suggestions and solutions on how to fix this are highly appreciated.
Below is the code and it has to be saved in 3 files to run it since i use 2 functions and the 3rd file i call the 2 functions.
%---------------------------------Example (i use to call both functions)--------------------------------------
% Example with an ellipse
theta = linspace(0, 2 * pi, 8);
t_x = cos(theta);
t_y = sin(theta);
% Interpolate x-coordinate
[b_x, c_x, d_x] = periodic_cubic_spline_interpolation(theta, t_x);
% Interpolate y-coordinate
[b_y, c_y, d_y] = periodic_cubic_spline_interpolation(theta, t_y);
% Evaluate the spline at points
t_eval = linspace(0, 2 * pi, 500);
x_eval = arrayfun(@(t) evaluate_cubic_spline(theta, t_x, b_x, c_x, d_x, t), t_eval);
y_eval = arrayfun(@(t) evaluate_cubic_spline(theta, t_y, b_y, c_y, d_y, t), t_eval);
% Plot the results
figure;
plot(t_x, t_y, 'ko', 'LineWidth', 1.5);
hold on;
grid on;
plot(x_eval, y_eval, 'r-', 'LineWidth', 2);
title('Cubic Spline Interpolation of Ellipse');
legend('True Ellipse', 'Cubic Spline');
xlabel('X-axis');
ylabel('Y-axis');
%---------------------------------------function 1--------------------------------------------------------
function [b, c, d] = periodic_cubic_spline_interpolation(t, x)
n = length(t);
h = diff(t);
alpha = zeros(n, 1);
beta = zeros(n, 1);
for i = 2:n - 1
alpha(i) = 3 * (x(i + 1) - x(i)) / h(i) - 3 * (x(i) - x(i - 1)) / h(i - 1);
beta(i) = 2 * (h(i - 1) + h(i));
end
% Periodic boundary conditions------%(these are not being implemented!!!)%--------------------------------
alpha(1) = 3 * (x(2) - x(1)) / h(1) - 3 * (x(1) - x(n)) / h(n - 1);
alpha(n) = 3 * (x(1) - x(n)) / h(n - 1) - 3 * (x(n) - x(n - 1)) / h(n - 2);
beta(1) = 2 * (h(n - 1) + h(1));
beta(n) = 2 * (h(n - 1) + h(n - 2));
l = zeros(n, 1);
mu = zeros(n, 1);
z = zeros(n, 1);
l(1) = 1;
mu(1) = 0;
z(1) = 0;
for i = 2:n - 1
l(i) = 2 * (t(i + 1) - t(i - 1)) - h(i - 1) * mu(i - 1);
mu(i) = h(i) / l(i);
z(i) = (alpha(i) - h(i - 1) * z(i - 1)) / l(i);
end
l(n) = 1;
z(n) = 0;
c = zeros(n, 1);
b = zeros(n, 1);
d = zeros(n, 1);
for j = n - 1:-1:1
c(j) = z(j) - mu(j) * c(j + 1);
b(j) = (x(j + 1) - x(j)) / h(j) - h(j) * (c(j + 1) + 2 * c(j)) / 3;
d(j) = (c(j + 1) - c(j)) / (3 * h(j));
end
end
%------------------------------------function 2---------------------------------------------------------
function spline_eval = evaluate_cubic_spline(t, x, b, c, d, t_eval)
n = length(t);
j = max(min(find(t >= t_eval, 1) - 1, n - 2), 1);
tj = t(j);
spline_eval = x(j) + b(j) * (t_eval - tj) + c(j) * (t_eval - tj) ^ 2 + d(j) * (t_eval - tj) ^ 3;
end

답변 (1개)

Matt J
Matt J 2023년 12월 2일
편집: Matt J 2023년 12월 2일
Can't you just use the native spline() command?
theta = linspace(0, 2 * pi, 8);
t_x = 2*cos(theta);
t_y = sin(theta);
dt_x=0; %end slope x
dt_y=(t_y(2)-t_y(1))./(theta(2)-theta(1)); %end slope y
theta_eval = linspace(0, 2 * pi, 500);
t_eval=periodicSpline(theta,[t_x;t_y], [dt_x; dt_y], theta_eval);
plot(t_x, t_y, 'ko', 'LineWidth', 1.5);
hold on;
grid on;
plot(t_eval(1,:), t_eval(2,:), 'r-', 'LineWidth', 2); axis equal; axis padded
title('Cubic Spline Interpolation of Ellipse');
legend('True Ellipse', 'Cubic Spline');
xlabel('X-axis');
ylabel('Y-axis');
hold off
function out = periodicSpline(tdata,xydata,xyslope, tquery)
xydata(:,end)=xydata(:,1);
xydata=[xyslope, xydata,xyslope];
out=spline(tdata,xydata,tquery);
end
  댓글 수: 8
Torsten
Torsten 2023년 12월 3일
편집: Torsten 2023년 12월 3일
Make the ansatz
p_j(x) = a_j + b_j*(x-x_j) + c_j*(x-x_j)^2 + d_j*(x-x_j)^3 (j=1,...,N-1)
where N is the number of x-data points.
So in total you have 4*(N-1) unknowns a_j, b_j, c_j and d_j that are to be determined from the 4*(N-1) conditions
p_j(x_j) = y_j (j=1,...,N-1)
p_(N-1)(x_N) = p_1(x_1)
p_j(x_j) = p_(j-1)(x_j) (j=2,...N-1)
p_(N-1)'(x_N) = p_1'(x_1)
p_j'(x_j) = p_(j-1)'(x_j) (j=2,...,N-1)
p_(N-1)''(x_N) = p_1''(x_1)
p_j''(x_j) = p_(j-1)''(x_j) (j=2,...,N-1)
This gives a linear system of equations in the unknowns a_j, b_j, c_j and d_j. Solve it.
Panda05
Panda05 2023년 12월 3일
Thank you so much everyone ! I tried to implement that and it worked . May God bless you

댓글을 달려면 로그인하십시오.

카테고리

Help CenterFile Exchange에서 Spline Postprocessing에 대해 자세히 알아보기

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by