Trajectorys of 4 objects on 3D plane not symmetrical
조회 수: 2 (최근 30일)
이전 댓글 표시
Hi all,
Sorry i'm a beginner at Matlab and am having issues with getting 4 spots on a 3D plane with the intitial conditions (0,0,0), (1,0,3), (1,2,0), (0,2,3) to symmetrically meet at (0.5,1,1.5).
My 2D simulation appeared like this:
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/999315/image.png)
But from my 3D simulation they appeared like this:
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/999320/image.png)
Im unsure what line of code I need to change to fix this from my code.
My 3D code is given below:
clc, clear, clear all
ti=0;
tf=1000;
tspan= linspace(ti,tf,1000000);
f0 = [0; 0; 0; 1; 2; 0; 1; 0; 3; 0; 2; 3];
[t, f] = ode45(@g, tspan, f0);
%Trajectories are done in 3 coordinates x,y,z
i_1=(f(:,1)); %Trajectory of bug 1
j_1=(f(:,5));
k_1=(f(:,9));
i_2=(f(:,2)); %Trajectory of bug 2
j_2=(f(:,6));
k_2=(f(:,10));
i_3=(f(:,3)); %Trajectory of bug 3
j_3=(f(:,7));
k_3=(f(:,11));
i_4=(f(:,4)); %Trajectory of bug 4
j_4=(f(:,8));
k_4=(f(:,12));
figure
plot3(i_1, j_1, k_1)
hold on
plot3(i_2, j_2, k_2)
plot3(i_3, j_3, k_3)
plot3(i_4, j_4, k_4)
hold off
title('test')
xlabel('x')
ylabel('y')
zlabel('z')
function dxdt = g(t,f)
a = sqrt((f(2)-f(1))^2 + (f(6)-f(5))^2 + (f(10)-f(9))^2);
b = sqrt((f(3)-f(2))^2 + (f(7)-f(6))^2 + (f(11)-f(10))^2);
c = sqrt((f(4)-f(3))^2 + (f(8)-f(7))^2 + (f(12)-f(11))^2);
d = sqrt((f(1)-f(4))^2 + (f(5)-f(8))^2 + (f(9) -f(12))^2);
dxdt = [ ...
(f(2) - f(1)) / a;
(f(6) - f(5)) / a;
(f(10) - f(9)) / a;
(f(3) - f(2)) / b;
(f(7) - f(6)) / b;
(f(11) - f(10)) / b;
(f(4) - f(3)) / c;
(f(8) - f(7)) / c;
(f(12) - f(11)) / c;
(f(1) - f(4)) / d;
(f(5) - f(8)) / d;
(f(9) - f(12)) / d];
end
%My 2D code is here:
tinitial = 0;
tfinal = 1000;
tspan = linspace(tinitial, tfinal, 1000000);
f0 = [0; 0; 1; 1; 0; 2; 2; 0];
[t, f] = ode45(@g, tspan, f0);
j = (f(:, 1));
k = (f(:, 5));
a = (f(:, 2));
b = (f(:, 6));
d = (f(:, 3));
e = (f(:, 7));
h = (f(:, 4));
i = (f(:, 8));
figure
hold on
plot3(j, k, tspan)
plot3(a, b, tspan)
plot3(d, e, tspan)
plot3(h, i, tspan)
title('test')
xlabel('x')
ylabel('y')
zlabel('t')
function dfdt = g(t,f)
dfdt = [(f(2) - f(1))/(sqrt((f(2)-f(1))^2 + (f(6)-f(5))^2));
(f(3)-f(2))/(sqrt((f(3)-f(2))^2 + (f(7)-f(6))^2));
(f(4) - f(3))/(sqrt((f(4)-f(3))^2 + (f(8)-f(7))^2));
(f(1) - f(4))/(sqrt((f(1)-f(4))^2 + (f(5)-f(8))^2));
(f(6) - f(5))/(sqrt((f(2)-f(1))^2 + (f(6)-f(5))^2));
(f(7) - f(6))/(sqrt((f(3)-f(2))^2 + (f(7)-f(6))^2));
(f(8) - f(7))/(sqrt((f(4)-f(3))^2 + (f(8)-f(7))^2));
(f(5) - f(8))/(sqrt((f(1)-f(4))^2 + (f(5)-f(8))^2))];
end
Cheers
댓글 수: 2
Torsten
2022년 5월 16일
I you change the plot command in the 2d-case as you did in the 3d-case, are the curves still as expected ?
figure
plot3(j,k,tspan)
hold on
plot3(a,b,tspan)
plot3(d,e,tspan)
plot3(h,i,tspan)
hold off
title('test')
xlabel('x')
ylabel('y')
zlabel('t')
채택된 답변
Bjorn Gustavsson
2022년 5월 16일
편집: Bjorn Gustavsson
2022년 5월 16일
Your main problem (or rather: step to take) is that you do the calculations element-by-element in the different arrays. That is not the way to go through life. You will have a much easier time if you start working with the arrays using vector-algebra. If you do that correctly there is no real difference between the 2-D and 3-D case. For your problem I guess that at each time-step you want to move the points towards the next point, some fraction of the distance, or something like that. For this you should use something like (for moving point 2 towards point 1):
dr21 = r1 - r2;
l_21 = norm(dr21);
lstep = l_21/30;
r2_next = r2 + (lstep/l_21)*dr21;
Then you have to put this into your loop and do it for all 4 points.
This way of working with vector-algebra removes (almost) all the thinking you have to do for a problem like this. 5-10 minutes with pen and paper sketching arrows and scribling trivial additions and subtractions will cut down the time required by factors of 10. Not just a single factor of 10, factors...
HTH
댓글 수: 6
Bjorn Gustavsson
2022년 5월 18일
Ah, sorry about that. I just realised that I confused/combined your question with someone elses where the task was to compare the ODE-solution with a simple fixed-size stepping solution. For that "someonelses problem" I thought that a step-length of 1/30 might be sensible if starting at coordinates of single-digit magnitudes. Just cut out that constant - in your setting I guess it would simply be something like a scalling of time.
Bjorn Gustavsson
2022년 5월 18일
Yes
r2 = f0([1 5 9]+1);
r2 == f0([1 5 9]+1);
Matlab's indexing capabilities are ridiculously handy! It lets you do things like:
A = [1 2;3 4];
A = zeros(size(A));
A([1 3 4]) = rand(1,3);
C = zeros(size(A));
C(:) = randn(4,1);
which can be increadibly useful - not that these examples are, but these capacities have been very handy for problems I've had to solve.
Glad that the suggestion helped. If you have time once you've solved the problem we'd be happy to see the result - everyone likes a nice figure.
추가 답변 (0개)
참고 항목
카테고리
Help Center 및 File Exchange에서 Pulsed Waveforms에 대해 자세히 알아보기
제품
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!