Ive written the following code in MatLab but I need help creating a graph that accurately depicts whats occuring. I am using potential energy and spring constants to calculate the effects of celestial bodies orbiting. please let me know if you notice any mistakes in my code. I would like the graph to be a accurate depiction of the Orbit of the Earth moon and sun based upon the formulas in my code. $\%$ Define masses

m$1=1\mathrm{.}9891$e$30$; $\%$ Mass of particle $1($Sun$)$

m$2=5\mathrm{.}97219$e$24$; $\%$ Mass of particle $2($Earth$)$

m$3=7\mathrm{.}34767$e$22$; $\%$ Mass of particle $3($Moon$)$

$\%$ Define spring constants

k$12=0\mathrm{.}1$; $\%$ Spring constant between particle $1$ and particle $2$

k$23=0\mathrm{.}5$; $\%$ Spring constant between particle $2$ and particle $3$

$\%$ Define initial positions

x$1=0$; y$1=0$; $\%$ Initial position of particle $1($Sun$)$

x$2=149\mathrm{.}6$e$9$; y$2=0$; $\%$ Initial position of particle $2($Earth$)$

x$3=$ x$2+384\mathrm{.}4$e$6$; y$3=0$; $\%$ Initial position of particle $3($Moon$)$

$\%$ Define initial velocities

vx$1=0$; vy$1=0$; $\%$ Initial velocity of particle $1($Sun$)$

vx$2=0$; vy$2=29\mathrm{.}783$e$3$; $\%$ Initial velocity of particle $2($Earth$)$

vx$3=0$; vy$3=29\mathrm{.}783$e$3+1\mathrm{.}022$e$3$; $\%$ Initial velocity of particle $3($Moon$)$

$\%$ Define time step and number of iterations

dt $=3600$; $\%$ Time step

n $=365*24$; $\%$ Number of iterations $(1$ year$)$

$\%$ Initialize arrays to store positions

x$1\_$array $=$ zeros$(1,$ n$)$; y$1\_$array $=$ zeros$(1,$ n$)$;

x$2\_$array $=$ zeros$(1,$ n$)$; y$2\_$array $=$ zeros$(1,$ n$)$;

x$3\_$array $=$ zeros$(1,$ n$)$; y$3\_$array $=$ zeros$(1,$ n$)$;

$\%$ Store initial positions

x$1\_$array$(1)=$ x$1$; y$1\_$array$(1)=$ y$1$;

x$2\_$array$(1)=$ x$2$; y$2\_$array$(1)=$ y$2$;

x$3\_$array$(1)=$ x$3$; y$3\_$array$(1)=$ y$3$;

$\%$ Perform the iterations and calculate positions

for i $=2$:n

$\%$ Calculate distances between particles

r$12=$ sqrt$(($x$2-$ x$1)^2+($y$2-$ y$1)^2)$;

r$23=$ sqrt$(($x$3-$ x$2)^2+($y$3-$ y$2)^2)$;

$\%$ Calculate forces using the spring potential

F$12=-$k$12*($r$12-1)$; $\%$ Assuming equilibrium distance r$0=1$ for particle $1-2$

F$23=-$k$23*($r$23-1)$; $\%$ Assuming equilibrium distance r$0=1$ for particle $2-3$

$\%$ Calculate accelerations

ax$1=$ F$12*($x$2-$ x$1)/$ r$12/$ m$1$;

ay$1=$ F$12*($y$2-$ y$1)/$ r$12/$ m$1$;

ax$2=($F$12*($x$1-$ x$2)/$ r$12+$ F$23*($x$3-$ x$2)/$ r$23)/$ m$2$;

ay$2=($F$12*($y$1-$ y$2)/$ r$12+$ F$23*($y$3-$ y$2)/$ r$23)/$ m$2$;

ax$3=$ F$23*($x$2-$ x$3)/$ r$23/$ m$3$;

ay$3=$ F$23*($y$2-$ y$3)/$ r$23/$ m$3$;

$\%$ Update velocities

vx$1=$ vx$1+$ ax$1*$ dt; vy$1=$ vy$1+$ ay$1*$ dt;

vx$2=$ vx$2+$ ax$2*$ dt; vy$2=$ vy$2+$ ay$2*$ dt;

vx$3=$ vx$3+$ ax$3*$ dt; vy$3=$ vy$3+$ ay$3*$ dt;

$\%$ Update positions

x$1=$ x$1+$ vx$1*$ dt; y$1=$ y$1+$ vy$1*$ dt;

x$2=$ x$2+$ vx$2*$ dt; y$2=$ y$2+$ vy$2*$ dt;

x$3=$ x$3+$ vx$3*$ dt; y$3=$ y$3+$ vy$3*$ dt;

$\%$ Store positions in arrays

x$1\_$array$($i$)=$ x$1$; y$1\_$array$($i$)=$ y$1$;

x$2\_$array$($i$)=$ x$2$; y$2\_$array$($i$)=$ y$2$;

x$3\_$array$($i$)=$ x$3$; y$3\_$array$($i$)=$ y$3$;

end