small changes in the examples

examples/rossler_system.py

@@ -1,27 +1,55 @@
 import numpy as np
 import multiflap as mf
 import matplotlib.pyplot as plt
-
-x = [10., 10., 3.6]
+from mpl_toolkits.mplot3d import Axes3D
+from scipy.integrate import odeint
+x = [.1, 5., 3.6]
 
 time_array = np.linspace(0, 180, 90000)
 mymodel = mf.Rossler(a=0.2, b=0.2, c=5.7)
 
-ms_obj =  mf.MultipleShootingPeriod(x, M=2, period_guess= 5., t_steps=50000, model=mymodel)
+ms_obj =  mf.MultipleShootingPeriod(x, M=30, period_guess= 5.,
+                                    t_steps=5000, model=mymodel, integrator='odeint')
 
-mysol = mf.SolverPeriod(ms_obj = ms_obj).lma(5.)
+mysol = mf.SolverPeriod(ms_obj = ms_obj).lma()
 
 jac = mysol[4]
 
 eigenvalues, eigenvectors = np.linalg.eig(jac)
 
+# plot the trajectory related to the initial condition
+[_, initial_traj] = ms_obj.get_mappedpoint(x, 0, 100)
+initial_points = ms_obj.get_initial_guess()
 
 sol_array = mysol[3].space
 sol_time = mysol[3].time
 period = sol_time[-1]
 
-plt.plot( sol_time, sol_array[:,0], label = "D1")
-plt.plot( sol_time, sol_array[:,1], label = "D2")
-plt.plot( sol_time, sol_array[:,2], label = "R")
+initial_traj = initial_traj.x
+
+fig1 = plt.figure(1)
+ax = fig1.gca(projection='3d')
+ax.set_xlabel('$x$')
+ax.set_ylabel('$y$')
+ax.set_zlabel('$z$')
+ax.plot(sol_array[:, 0],
+        sol_array[:, 1],
+        sol_array[:, 2],color = 'b')
+ax.scatter(initial_points[:, 0],
+        initial_points[:, 1],
+        initial_points[:, 2],color = 'red', alpha=0.5)
+plt.legend()
+plt.show()
+fig2 = plt.figure(2)
+ax = fig2.gca(projection='3d')
+ax.set_xlabel('$x$')
+ax.set_ylabel('$y$')
+ax.set_zlabel('$z$')
+ax.plot(sol_array[:, 0],
+        sol_array[:, 1],
+        sol_array[:, 2],color = 'b', label='limit cycle')
+ax.plot(initial_traj[:, 0],
+        initial_traj[:, 1],
+        initial_traj[:, 2],color = 'red', alpha=0.5, label='initial point trajectory')
 plt.legend()
 plt.show()

multiflap/ms_package/multiple_shooting.py

@@ -5,7 +5,7 @@
 
 class MultipleShooting:
 
-    def __init__(self, x0, M = 2, period = 0.25, 
+    def __init__(self, x0, M = 2, period = 0.25,
                  t_steps = 70, model = None, option_jacobian = 'analytical'):
         self.point_number = M
         self.dim = model.dimension  # number of states of the ODE system
@@ -337,7 +337,6 @@ def get_ms_scheme(self, x0):
                (i*self.dim):(i*self.dim)+self.dim] = jacobian
 
 
-        #trajectory = np.asanyarray(complete_solution)
         # Last block of the scheme
         MS[-self.dim:, 0:self.dim] = -np.eye(self.dim)
         MS[-self.dim:, -self.dim:] = np.eye(self.dim)