How do I make a decaying oscilating function in python?

1.1k views Asked by At

I have a code in python to represent the energy decay in a damped oscilator, it reads like this:

def E(wt, Q):
    return (np.e**(-x/Q))*(1-(1/2*Q)*np.sin(2*x))
x = np.linspace(0,20,1000)
y0 = E(x,2)
y1 = E(x,4)
y2 = E(x,8)
y3 = E(x,16)
plt.plot(x, y0, 'p', label=r'$Q=2$')
plt.plot(x, y1, 'r', label=r'$Q=4$')
plt.plot(x, y2, 'g', label=r'$Q=8$')
plt.plot(x, y3, 'b', label=r'$Q=16$')
plt.xlabel(r'$wt$')
plt.ylabel(r'$E$')
plt.title (r"$E(t)  -vs.- wt$")
plt.show()

yet it should look like this: https://www.dropbox.com/s/o2mmmi8v6kdnn2v/good_graph.png?dl=0 what am I doing wrong? I have the right function

2

There are 2 answers

0
Alexander McFarlane On BEST ANSWER

Fixed Equation

def E(wt, Q):
    return np.exp(-x/float(Q)) * ( 1. - (1./2./float(Q))*np.sin(2.* x) )

Your original equation

def E(wt, Q):
    return (np.e**(-x/Q))*(1-(1/2*Q)*np.sin(2*x))

Errors

  1. Unused Variables You never use wt
  2. BODMAS You don't have a decay parameter set up correctly, so it will oscillate too much. (1/2*Q) when you mean (1/2/Q)
  3. Integer Division (Only for Python < 3) You divide by integers which will floor your division. You need to convert the integer values to floats e.g. (1/2./float(Q))
  4. Color Parameter You want to make a purple line by passing p in plt.plot(x, y0, 'p', label=r'$Q=2$') but p creates the strange point plot behaviour. To fix this pass in the color names explicitly e.g. plt.plot(x, y0, color='purple', label=r'$Q=2$')

Edited Layout - see below

Off Topic

To make your title nicer:

from matplotlib import rc
rc('font', **{'family': 'serif', 'serif': ['Computer Modern']})
rc('text', usetex=True)
plt.title (r"$E(t)$ vs. $w_t$")
plt.xlabel(r'$w_t$')
plt.ylabel(r'$E(t)$')

Full Code

from matplotlib import rc
import bumpy as np
import matplotlib.pyplot as plot

# set up fonts
rc('font', **{'family': 'serif', 'serif': ['Computer Modern']})
rc('text', usetex=True)

# set up labels
plt.title (r"$E(t)$ vs. $w_t$")
plt.xlabel(r'$w_t$')
plt.ylabel(r'$E(t)$')

# plot and store plots in yList, zip color labels into loop
yList = []
for i,color in zip(xrange(1,5),['purple', 'r', 'g', 'b']):
    y = E(x, 2**i)
    yList.append(y)
    plt.plot(x, y, color = color, label= r'$Q=%s$' % i)
pt.show()
0
jcoppens On

It seems the function isn't quite right. I suspect the original formula was right but you didn't encode it right in Python.

import numpy as np
import pylab as plt
def E(x, Q):
    return np.exp(-x/Q) * (1-(1.0/2.0/Q) * np.sin(2*x))
def main():
    x = np.linspace(0,20,1000)
    y0 = E(x,2)
    y1 = E(x,4)
    y2 = E(x,8)
    y3 = E(x,16)
    plt.plot(x, y0, 'p', label=r'$Q=2$')
    plt.plot(x, y1, 'r', label=r'$Q=4$')
    plt.plot(x, y2, 'g', label=r'$Q=8$')
    plt.plot(x, y3, 'b', label=r'$Q=16$')
    plt.xlabel(r'$wt$')
    plt.ylabel(r'$E$')
    plt.title (r"$E(t)  -vs.- wt$")
    plt.show()    
    return 0
if __name__ == '__main__':
    main()*

enter image description here