接下来我们要编写以下的语句在循环中调用自定义函数rk4.m，并为其提供相应的参数
其中'gravrk'表示的是公式（4）中等号右边的函数形式，在这个问题中应为：
将其编写成如下的自定义函数gravrk(s,t,GM)，并存为gravrk.m
到此，我们就有了完整的计算程序，接下来需要检验一下它的准确性。
1）轨道的稳定性及机械能守恒
宇宙中天体的运行轨道形状和周期是非常稳定的，正是因为这样，我们的祖先往往以此来制定他们的时间标准。这里，我们运用以上的程序计算几个周期的运行轨道，可以由轨道的重合程度来判断计算结果的准确性。同时，引力是一种保守力，在保守力场中运动的物体机械能是守恒的（不考虑相对论效应）。运动天体的机械能包括两部分，即动能
和势能
,因此我们还可以通过对机械能的研究来判断计算结果的准确性。
