<<PhysicalConstants`

kappa=GravitationalConstant[[1]];

casTed=AbsoluteTime[Date[]]

Δt=120;
casMinus=casTed-Δt;
casPlus=casTed+Δt;

rZeme0=AstronomicalData["Earth",{"Position",DateList[casTed]}];

vZeme0=(1/(2*Δt))(AstronomicalData["Earth",{"Position",DateList[casPlus]}]-AstronomicalData["Earth",{"Position",DateList[casMinus]}]);

mZeme=AstronomicalData["Earth","Mass"];
mSlunce=AstronomicalData["Sun","Mass"];
silaSlunceNaZemi=kappa*((mZeme*mSlunce)/((Norm[-rZeme[t]+rSlunce[t]])^3))*(-rZeme[t]+rSlunce[t])/.rSlunce[t]->0

rov=mZeme* rZeme''[t]==silaSlunceNaZemi
pocPodm={rZeme[0]==rZeme0,rZeme'[0]==vZeme0}
rovs=Flatten[{rov,pocPodm}]

tmax=3600*24*365.25;
res1=NDSolve[rovs,rZeme,{t,0,tmax}][[1]]

fce=rZeme/.res1;

rovs

cara2=Graphics[{Black,Thickness[0.006],Line[souradnice2D/@Range[0,tmax,0.2tmax]]}];

souradnice2D[t_]:=Take[fce[t],2];
bod[t_]:=Graphics[{Red,PointSize[0.1],Point[souradnice2D[t]]}];
cara2=Graphics[{Black,Thickness[0.006],Line[souradnice2D/@Range[0,tmax,0.01tmax]]}];
rn=5*(10^11);
Manipulate[
Show[bod[t],cara2,PlotRange->{{-rn,rn},{-rn,rn}}],
{t,0,tmax},ControlPlacement->Top
];

silaSlunceNaZemi=kappa*((mZeme*mSlunce)/((Norm[-rZeme[t]+rSlunce[t]])^3))*(-rZeme[t]+rSlunce[t]);
rovZeme={mZeme* rZeme''[t]==silaSlunceNaZemi};
rovSlunce={mSlunce* rSlunce''[t]==-silaSlunceNaZemi};
pocPodm={rZeme[0]==rZeme0,rZeme'[0]==vZeme0,rSlunce[0]=={0,0,0},rSlunce'[0]=={0,0,0}};
rovs=Union[rovSlunce,rovZeme,pocPodm];
tmax=4*3600*24*365.25;
res1=NDSolve[rovs,{rZeme,rSlunce},{t,0,tmax}][[1]];

TreeForm[ rZeme''[t]]

Union[Cases[rovs,a_Symbol[t]:>a,{0,Infinity}]]

nezname=Union[Cases[rovs,a_Symbol[t]:>a,{0,Infinity}]]

souradnice2DSlunce[t_]:=Take[rSlunce[t]/.res1,2];
bodSlunce[t_]:=Graphics[{Red,PointSize[0.1],Point[souradnice2DSlunce[t]]}];
rn=7*(10^6);
cara2=Graphics[{Black,Thickness[0.006],Line[souradnice2DSlunce/@Range[0,tmax,0.01tmax]]}];
Manipulate[
Show[bodSlunce[t],cara2,PlotRange->{{-rn,rn},{-rn,1.5rn}}],
{t,0,tmax},ControlPlacement->Top
]