%matplotlib inline
from sympy import *
init_printing()
from sympy import *
x, t, z, nu = symbols('x t z nu')
diff(sin(x)*exp(x), x)
from sympy import symbols
from sympy.plotting import plot,plot3d
from sympy.abc import *
x = symbols('x')
plot(x,(x,-5,5))
from sympy import symbols
from sympy.plotting import plot
p1 = plot(x*x,xlabel='teste',line_color=0.0)
p2 = plot(x)
p1.extend(p2)
p1
x+1
plot3d((x**2 + y**2, (x, -5, 5), (y, -5, 5)),(x*y, (x, -3, 3), (y, -3, 3)))
plot3d(sin(x*10)*cos(y*5) - x*y,(x,-10,10),(y,-10,10))
import numpy as np
import scipy.integrate as spi
import matplotlib.pyplot as plt
m = 1. # particle's mass
k = 1. # drag coefficient
g = 9.81 # gravity acceleration
# We want to evaluate the system on 30 linearly
# spaced times between t=0 and t=3.
t = np.linspace(0., 3., 30)
# We simulate the system for different values of k.
for k in np.linspace(0., 1., 5):
# We simulate the system and evaluate $v$ on
# the given times.
v = spi.odeint(f, v0, t, args=(k,))
# We plot the particle's trajectory.
plt.plot(v[:,0], v[:,1], 'o-',
mew=1, ms=8, mec='w',
label='k={0:.1f}'.format(k))
plt.legend()
plt.xlim(0, 12)
plt.ylim(0, 6)
# The initial position is (0, 0).
v0 = np.zeros(4)
# The initial speed vector is oriented
# to the top right.
v0[2] = 4.
v0[3] = 10.
def f(v, t0, k):
# v has four components: v=[u, u'].
u, udot = v[:2], v[2:]
# We compute the second derivative u'' of u.
udotdot = -k/m * udot
udotdot[1] -= g
# We return v'=[u', u''].
return np.r_[udot, udotdot]
m = 1. # particle's mass
b = 1. # drag coefficient
k= 1. #spring constant
# The initial position is (0, 0).
v0 = np.zeros(2)
# The initial speed vector is oriented
# to the top right.
v0[0] = 1.
v0[1] = 0.
def f(v, t0, *args):
# v has four components: v=[u, u'].
b,k=args
u, udot = v[0], v[1]
# We compute the second derivative u'' of u.
udotdot = -b/m*udot-k/m*u
# We return v'=[u', u''].
return np.array([udot,udotdot])
t = np.linspace(0., 10., 300)
v = spi.odeint(f, v0, t, args=(4,2))
plt.plot(t,v[:,0], '-',
mew=1, ms=8, mec='w')
plt.ylabel('$x$')
plt.xlabel('$t$')
plt.xlabel('$x$')
plt.ylabel('$\dot{x}$')
plt.plot(v[:,0],v[:,1], '-',
mew=1, ms=8, mec='w')
from random import random
from math import cos,sin
t = np.linspace(0., 30., 300)
for theta in np.linspace(0, 6.28, 100):
v0=np.array([float(3*cos(theta)),float(3*sin(theta))])
# We simulate the system and evaluate $v$ on
# the given times.
v = spi.odeint(f, v0, t, args=(2,1))
# We plot the particle's trajectory.
plt.plot(v[:,0], v[:,1], '-',mew=1, ms=8, mec='w')
plt.legend()
plt.xlabel('$x$')
plt.ylabel('$\dot{x}$')
plt.title('$b=2$')
plt.xlim(-2, 2)
plt.ylim(-2, 2)
t = np.linspace(0., 30., 300)
nmax=6
fig, axarr = plt.subplots(nmax, nmax,figsize=(16,12))
for i in range(nmax):
for j in range(nmax):
for theta in np.linspace(0, 6.28, 100):
v0=np.array([float(3*cos(theta)),float(3*sin(theta))])
# We simulate the system and evaluate $v$ on
# the given times.
v = spi.odeint(f, v0, t, args=(i+1,nmax-j))
# We plot the particle's trajectory.
axarr[i,j].plot(v[:,0], v[:,1], '-',mew=1, ms=8, mec='w')
if (4*(nmax-j)<=(i+1)**2):
color='yellow'
else:
color='orange'
axarr[nmax-i-1,nmax-j-1].set_facecolor(color)
axarr[nmax-i-1,nmax-j-1].set_xlim(-2,2)
axarr[nmax-i-1,nmax-j-1].set_ylim(-2,2)
axarr[nmax-i-1,nmax-j-1].set_title('$b$={:2.1f} $k$={:2.1f}'
.format((i+1)*1.0,(nmax-j)*1.0))
plt.setp([axarr[i, j].get_xticklabels() for i in range(nmax) for j in range(nmax)], visible=False)
plt.setp([axarr[i, j].get_yticklabels() for i in range(nmax) for j in range(nmax)], visible=False)
axarr[nmax-1,2].set_xlabel(' Atrito',size=24)
axarr[3,0].set_ylabel(' Constante da Mola',size=24)
plt.suptitle('Osciladores Amortecidos',size=24)
plt.show()
def ff(v, t0, *args):
# v has four components: v=[u, u'].
b,k,A,omega=args
u, udot = v[0], v[1]
# We compute the second derivative u'' of u.
udotdot = -b/m*udot-k/m*u+A*cos(omega*t0)/m
# We return v'=[u', u''].
return np.array([udot,udotdot])
t = np.linspace(0., 30., 300)
v = spi.odeint(ff, v0, t, args=(4,2,1,1))
plt.plot(t,v[:,0], '-',
mew=1, ms=8, mec='w')
plt.ylabel('$x$')
plt.xlabel('$t$')
plt.xlabel('$x$')
plt.ylabel('$\dot{x}$')
plt.plot(v[:,0],v[:,1], '-',
mew=1, ms=8, mec='w')
t = np.linspace(0., 300., 300)
A=2
omega=2
b=2
k=1
for theta in np.linspace(0, 6.28, 100):
v0=np.array([float(3*cos(theta)),float(3*sin(theta))])
# We simulate the system and evaluate $v$ on
# the given times.
v = spi.odeint(ff, v0, t, args=(b,k,A,omega))
# We plot the particle's trajectory.
plt.plot(v[:,0], v[:,1], '-',mew=1, ms=8, mec='w')
plt.legend()
plt.xlabel('$x$',size=16)
plt.ylabel('$\dot{x}$',size=16)
plt.title('$b$={} $k$={} $A$={} $\omega$={}'.format(b,k,A,omega))
plt.xlim(-2, 2)
plt.ylim(-2, 2)
tmax=3000
steps=30000
t = np.linspace(0., tmax, steps)
def maxinv(omega):
v = spi.odeint(ff, v0, t, args=(0.1,4,1,omega))
return [omega,max(v[-10000:,0]),min(-v[-10000:,0])]
Amax=np.array([maxinv(omega) for omega in np.linspace(0.1, 5, num=100)])
plt.plot(Amax[:,0],Amax[:,1])
plt.plot(Amax[:,0],-Amax[:,2])
plt.show
v = spi.odeint(ff, v0, t, args=(2,5,1,0.1))
plt.plot(t[-1000:],v[-1000:,0], '-',
mew=1, ms=8, mec='w')
plt.ylabel('$x$')
plt.xlabel('$t$')
print(max(v[-10000:,0]))
print(min(v[-10000:,0]))