Menu
▾
▴
Re: [Matplotlib-users] Lorenz - solution
|
From: Prahas D. N. <pra...@gm...> - 2015-03-10 16:41:56
|
@ Adam -- thanks! @ Everyone: given the Lorenz code shared, is there a way to generate a log file of the x,y,z points generated? Thanks in advance. --Prahas On Tue, Mar 10, 2015 at 8:59 AM, Adam Hughes <hug...@gm...> wrote: > That's pretty swag! > > On Tue, Mar 10, 2015 at 11:49 AM, Benjamin Root <ben...@ou...> wrote: >> >> +1000!! >> >> Great job! Would you mind if I clean it up a bit and add it to the >> mplot3d/animation gallery? Full credit, of course. >> >> >> On Tue, Mar 10, 2015 at 11:30 AM, Prahas David Nafissian >> <pra...@gm...> wrote: >>> >>> Friends, >>> >>> I thought you'd like to see the solution. >>> >>> Many thanks to Jake Vanderplas for his code and teachings: >>> >>> >>> https://jakevdp.github.io/blog/2013/02/16/animating-the-lorentz-system-in-3d/ >>> >>> If you start a new IP Notebook session, run as your first entry: >>> >>> %pylab >>> >>> and then copy and paste the text below and run it, you should be good to >>> go >>> (on a Mac, at least). >>> >>> There are several parameters I've changed from his original, and I've >>> commented as I've changed. The original code is at the link above. >>> >>> There is one error in his code -- I've documented it below. >>> >>> Again, thanks to the community, Jake, and Ben Root. >>> >>> --Prahas >>> >>> ****************** >>> >>> import numpy as np >>> from scipy import integrate >>> >>> from matplotlib import pyplot as plt >>> from mpl_toolkits.mplot3d import Axes3D >>> from matplotlib.colors import cnames >>> from matplotlib import animation >>> >>> # orig value of N_traj was 20 -- very cool this way. >>> >>> N_trajectories = 1 >>> >>> def lorentz_deriv((x, y, z), t0, sigma=10., beta=8./3, rho=28.0): >>> """Compute the time-derivative of a Lorentz system.""" >>> return [sigma * (y - x), x * (rho - z) - y, x * y - beta * z] >>> >>> # Choose random starting points, uniformly distributed from -15 to 15 >>> >>> np.random.seed(1) >>> >>> # changing from -15,30 to 10,5 below starts the drawing in the middle, >>> # rather than getting the long line from below.... >>> # if using N_Traj > 1, return to orig values. >>> >>> # x0 = -15 + 30 * np.random.random((N_trajectories, 3)) >>> >>> x0 = 10 + 5 * np.random.random((N_trajectories, 3)) >>> >>> >>> # Solve for the trajectories >>> >>> # orig values: 0,4,1000 >>> # 3rd value -- lower it, it gets choppier. >>> # 2nd value -- increase it -- more points, but speedier. >>> >>> # change middle num from 4 to 15 -- this adds points!!!!!!!! >>> >>> t = np.linspace(0, 40, 3000) >>> x_t = np.asarray([integrate.odeint(lorentz_deriv, x0i, t) >>> for x0i in x0]) >>> >>> # Set up figure & 3D axis for animation >>> fig = plt.figure() >>> ax = fig.add_axes([0, 0, 1, 1], projection='3d') >>> >>> # changing off to on below adds axises. slows it down but you >>> # can fix that with interval value in the animation call >>> >>> ax.axis('on') >>> >>> # choose a different color for each trajectory >>> colors = plt.cm.jet(np.linspace(0, 1, N_trajectories)) >>> >>> # set up lines and points -- this is a correction from >>> # the orig jake code. the next four lines... >>> >>> lines = [ax.plot([], [], [], '-', c=c)[0] >>> for c in colors] >>> pts = [ax.plot([], [], [], 'o', c=c)[0] >>> for c in colors] >>> >>> # prepare the axes limits >>> ax.set_xlim((-25, 25)) >>> ax.set_ylim((-35, 35)) >>> ax.set_zlim((5, 55)) >>> >>> # set point-of-view: specified by (altitude degrees, azimuth degrees) >>> ax.view_init(30, 0) >>> >>> # initialization function: plot the background of each frame >>> def init(): >>> for line, pt in zip(lines, pts): >>> line.set_data([], []) >>> line.set_3d_properties([]) >>> >>> pt.set_data([], []) >>> pt.set_3d_properties([]) >>> return lines + pts >>> >>> # animation function. This will be called sequentially with the frame >>> number >>> def animate(i): >>> # we'll step two time-steps per frame. This leads to nice results. >>> >>> i = (2 * i) % x_t.shape[1] >>> >>> for line, pt, xi in zip(lines, pts, x_t): >>> x, y, z = xi[:i].T >>> line.set_data(x, y) >>> line.set_3d_properties(z) >>> >>> pt.set_data(x[-1:], y[-1:]) >>> pt.set_3d_properties(z[-1:]) >>> >>> # changed 0.3 to 0.05 below -- this slows the rotation of the >>> view. >>> # changed 30 to 20 below >>> # changing 20 to (20 + (.1 * i)) rotates on the Z axis. trippy. >>> >>> ax.view_init(10, 0.1 * i) >>> # ax.view_init(10, 100) >>> fig.canvas.draw() >>> return lines + pts >>> >>> # instantiate the animator. I've deleted the blit switch (for Mac) >>> # enlarging frames=500 works now -- it failed before because I didn't >>> give it >>> # enough data -- by changing the t=np.linspace line above I generate >>> more points. >>> # interval larger slows it down >>> # changed inteval from 30 to 200, frames from 500 to 3000 >>> >>> anim = animation.FuncAnimation(fig, animate, init_func=init, >>> frames=3000, interval=200) >>> >>> # Save as mp4. This requires mplayer or ffmpeg to be installed. >>> COMPLEX!!!!! >>> # Instead, use a screen record program: Quicktime on the Mac; MS >>> Expression Encoder on PC. >>> # anim.save('PDNlorentz_attractor.mp4', fps=15, extra_args=['-vcodec', >>> 'libx264']) >>> >>> plt.show() >>> >>> ******************************** >>> >>> >>> ------------------------------------------------------------------------------ >>> Dive into the World of Parallel Programming The Go Parallel Website, >>> sponsored >>> by Intel and developed in partnership with Slashdot Media, is your hub >>> for all >>> things parallel software development, from weekly thought leadership >>> blogs to >>> news, videos, case studies, tutorials and more. Take a look and join the >>> conversation now. http://goparallel.sourceforge.net/ >>> _______________________________________________ >>> Matplotlib-users mailing list >>> Mat...@li... >>> https://lists.sourceforge.net/lists/listinfo/matplotlib-users >> >> >> >> >> ------------------------------------------------------------------------------ >> Dive into the World of Parallel Programming The Go Parallel Website, >> sponsored >> by Intel and developed in partnership with Slashdot Media, is your hub for >> all >> things parallel software development, from weekly thought leadership blogs >> to >> news, videos, case studies, tutorials and more. Take a look and join the >> conversation now. http://goparallel.sourceforge.net/ >> _______________________________________________ >> Matplotlib-users mailing list >> Mat...@li... >> https://lists.sourceforge.net/lists/listinfo/matplotlib-users >> > |
