diff --git a/lib/integrator.py b/lib/integrator.py index a20a639..eebfd00 100755 --- a/lib/integrator.py +++ b/lib/integrator.py @@ -24,7 +24,7 @@ def dp_dt(m_array, q_array): dp_array[np.isnan(dp_array)] = 0. return dp_array -def frogleap(duration, step, dyn_syst, display=False): +def frogleap(duration, step, dyn_syst, recover_param=False, display=False): """ Leapfrog integrator for first order partial differential equations. iteration : half-step drift -> full-step kick -> half-step drift @@ -34,6 +34,7 @@ def frogleap(duration, step, dyn_syst, display=False): q_array = dyn_syst.get_positions() p_array = dyn_syst.get_momenta() + if display: d = DynamicUpdate() d.min_x, d.max_x = -1.5*np.abs(q_array).max(), +1.5*np.abs(q_array).max() @@ -57,6 +58,7 @@ def frogleap(duration, step, dyn_syst, display=False): for i, body in enumerate(dyn_syst.bodylist): body.q = q_array[i] body.p = p_array[i] - body.v = body.p/body.m - - return dyn_syst + if body.m != 0.: + body.v = body.p/body.m + if recover_param: + return dyn_syst diff --git a/main.py b/main.py index e5e2df3..7f141e3 100755 --- a/main.py +++ b/main.py @@ -24,7 +24,8 @@ def main(): bodylist.append(Body(m[i], q[i], v[i])) dyn_syst = System(bodylist) - new_dyn_syst = frogleap(10, 0.01, dyn_syst, display=True) + new_dyn_syst = frogleap(10, 0.1, dyn_syst, display=False) + print(np.all(new_dyn_syst.get_positions() == dyn_syst.get_positions())) return 0 if __name__ == '__main__':