How can I Implement Runge Kutta 4 to plot the trajectory of a given point in a vector field?
5 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
Hello there,
I'm trying to approximate the trajectory of a particle with initial coordinates a,b and an initial velocity vector <P,Q> on a vector field of the form f(x,y) = <u(x,y), v(x,y)>.
Is there any form to calculate this using RK4? The output should be a plot of the vector field and the trajectory of the particle.
I have tried everything but I still got no solution to this problem, maybe you guys could tell me a different approach to solve this using Matlab.
Thanks for reading at this point and attending my request, have a nice day.
2 Commenti
Risposta accettata
Chunru
il 27 Lug 2021
As you have the intial point (0,0) and the intial speed (0,0), so the solution never progress to a different point. If you change the speed formla or initial point, then you can see something different.
[a,b] = meshgrid(0:0.2:2,0:0.2:2);
u = cos(a).*b;
v = sin(a).*b;
izq=0;
der=2;
h=0.2;
t=(izq:h:der);
%x0 = [0; 0];
x0 = [0.2; 0.2];
[t, y] = ode23(@odefun,t,x0);
figure
%plot(t,y);
plot(y(:,1), y(:,2))
hold on
quiver(a,b,u,v);
hold off
function y = odefun(t, x)
%u = cos(a).*b;
%v = sin(a).*b;
y(1,1) = cos(x(1)) .* x(2); % u(x, y)
y(2,1) = sin(x(1)) .* x(2); % v(x, y)
end
2 Commenti
Chunru
il 27 Lug 2021
Read https://blogs.mathworks.com/cleve/2014/05/26/ordinary-differential-equation-solvers-ode23-and-ode45/
For an ODE , you need to specify the function f(t, x). For 2D ODE, you specify f(t, x, y). In [t, y] = ode23(@odefun,t,x0) , @(odefun) define f(t, x, y) as a function handle; t is the time points and x0 is the initial value. Then ode23 is doing something similar to Runge Kutta you have implemented.
Più risposte (2)
Chunru
il 27 Lug 2021
Matlab has ode23 and ode45.
https://blogs.mathworks.com/cleve/2014/05/26/ordinary-differential-equation-solvers-ode23-and-ode45/
doc ode23 or ode45 for solving the differential equations.
0 Commenti
Chunru
il 27 Lug 2021
tspan = 0:.1:10;
x0 = [0; 0];
[t,x] = ode23(@odefun,tspan,x0);
subplot(211); plot(t, x); xlabel('t'); legend('x', 'y')
subplot(212); plot(x(:,1), x(:,2)); xlabel('x'); ylabel('y')
function y = odefun(t, x)
y(1,1) = (x(1).^2 + x(2).^2) * cos(x(1)) + .1; % u(x, y)
y(2,1) = (x(1).^2 + x(2).^2) * sin(x(2)) + .2; % v(x, y)
end
0 Commenti
Vedere anche
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!