How can i make my data more efficient in time when running script?
1 view (last 30 days)
Show older comments
I'm calculating the orbit of Earth and other planets around the sun using force of gravity. The script i have runs quick enough using a single planet, but it then begins to take longer as i add more planets (obviously).
I'm hoping if someone can give me some help on how to use the zeros function or another way to make the time taken to process the calculations more efficiently. I'm using a for loop to calculate and update each position through each iteration.
here is an example of the data:
days=365*48; % number of steps to make
stepsize=1/48; %(dt) or in Days
sunm =1.9891e30; % Sun Mass (kg)
G=1.4879e-34; % AU^3/(kg*Day^2)
%%Earth Initial Conditions
Ex=-9.320752400360321E-01; % Earth pos at x AU
Ey=-3.683263876537671E-01; % Earth pos at y AU
Ez=1.235087446387007E-05; % Earth pos at z AU
E_Vx=6.036562219130709E-03;% Earth x velocity AU/day
E_Vy=-1.606778488778610E-02;% Earth y velocity AU/day
E_Vz=6.657490105653882E-07; % Earth z velocity AU/day
%%Earth orbit
for count=1:days
% Save the current position of earth into a matrix
finalposition(count,1)=Ex;
finalposition(count,2)=Ey;
finalposition(count,3)=Ez;
%Calculate current Acceleration
R=norm([Ex,Ey,Ez]);
% assuming the sun is always at [0,0,0]
Atot=G*sunm/R^2;
unitvector=[Ex,Ey,Ez]/R;
Ax=-Atot*unitvector(1); %acceleration in the x
Ay=-Atot*unitvector(2); %acceleration in the y
Az=-Atot*unitvector(3); %acceleration in the z
%Update the position of Earth
Ex=Ex+E_Vx*stepsize+0.5*Ax*stepsize^2;
Ey=Ey+E_Vy*stepsize+0.5*Ax*stepsize^2;
Ez=Ez+E_Vz*stepsize+0.5*Ax*stepsize^2;
%Update the velocity of Earth
E_Vx=E_Vx+Ax*stepsize;
E_Vy=E_Vy+Ay*stepsize;
E_Vz=E_Vz+Az*stepsize;
end
hold on
plot3(finalposition(:,1),finalposition(:,2),finalposition(:,3),'b')
Thank you very much for your time, Alex :)
3 Comments
Titus Edelhofer
on 14 Apr 2015
Hi,
I would strongly suggest to use a "real" ODE solver instead of the update you are using. Computing the orbit requires quite some good precision to give meaningful results. Take a look at e.g. ode45 ...
Titus
Accepted Answer
the cyclist
on 14 Apr 2015
Edited: the cyclist
on 14 Apr 2015
Most importantly, not for speed, but for correctness, is that it looks like your formulas for "Update the position of the Earth" incorrectly use Ax for the y and z components.
Probably the biggest thing to help the speed will be to preallocate finalposition. Put this line ahead of your for loop:
finalposition = zeros(days,3);
This assigns memory for that variable ahead, rather than growing that array.
For speed, you could also just define one vector A, and do your updates like
A=-Atot*unitvector
Ditto for the other components.
More Answers (0)
See Also
Categories
Find more on Gravitation, Cosmology & Astrophysics in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!