In [6]:
import numpy as np
from resonance.nonlinear_systems import SingleDoFNonLinearSystem


To apply arbitrary forcing to a single degree of freedom linear or nonlinear system, you can do so with SingleDoFNonLinearSystem (SingleDoFLinearSystem does not support arbitrary forcing...yet).

Add constants, a generalized coordinate, and a generalized speed to the system.

In [7]:
sys = SingleDoFNonLinearSystem()

In [8]:
sys.constants['m'] = 100 # kg
sys.constants['c'] = 1.1*1.2*0.5/2
sys.constants['k'] = 10
sys.constants['Fo'] = 1000 # N
sys.constants['Ft'] = 100 # N/s
sys.constants['to'] = 3.0  # s

sys.coordinates['x'] = 0.0
sys.speeds['v'] = 0.0


Create a function that evaluates the first order form of the non-linear equations of motion. In this case:

$$\dot{x} = v \\ m\dot{v} + c \textrm{sgn}(v)v^2 + k \textrm{sgn}(x)x^2 = F(t)$$

Make the arbitrary forcing term, $F$, an input to this function.

In [9]:
def eval_eom(x, v, m, c, k, F):
xdot = v
vdot = (F - np.sign(v)*c*v**2 - np.sign(x)*k*x**2) / m
return xdot, vdot


Note that you cannot add this to the system because F has not been defined.

In [10]:
sys.diff_eq_func = eval_eom

---------------------------------------------------------------------------
KeyError                                  Traceback (most recent call last)
<ipython-input-10-bf0e21615aa6> in <module>()
----> 1 sys.diff_eq_func = eval_eom

/opt/conda/lib/python3.6/site-packages/resonance/nonlinear_systems.py in diff_eq_func(self, func)
83         # NOTE : This will throw an error if the function's args are not in the
84         # system.
---> 85         [self._get_par_vals(k) for k in getargspec(func).args]
86         self._diff_eq_func = func
87         self._check_diff_eq()

/opt/conda/lib/python3.6/site-packages/resonance/nonlinear_systems.py in <listcomp>(.0)
83         # NOTE : This will throw an error if the function's args are not in the
84         # system.
---> 85         [self._get_par_vals(k) for k in getargspec(func).args]
86         self._diff_eq_func = func
87         self._check_diff_eq()

/opt/conda/lib/python3.6/site-packages/resonance/system.py in _get_par_vals(self, par_name)
410             return self.measurements[par_name]
411         else:
--> 412             raise KeyError(msg.format(par_name))
413
414     def _check_meas_func(self, func):

KeyError: 'F is not in constants, coordinates, speeds, or measurements.'

To rememdy this, create a function that returns the input value given the appropriate constants and time.

In [11]:
def eval_step_input(Fo, to, time):
if time < to:
return 0.0
else:
return Fo

In [12]:
import matplotlib.pyplot as plt
%matplotlib widget

In [13]:
ts = np.linspace(0, 10)
plt.plot(ts, eval_step_input(5.0, 3.0, ts))

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
1 ts = np.linspace(0, 10)
----> 2 plt.plot(ts, eval_step_input(5.0, 3.0, ts))

<ipython-input-11-ec3c74db3f90> in eval_step_input(Fo, to, time)
1 def eval_step_input(Fo, to, time):
----> 2     if time < to:
3         return 0.0
4     else:
5         return Fo

ValueError: The truth value of an array with more than one element is ambiguous. Use a.any() or a.all()
In [14]:
ts < 3.0

Out[14]:
array([ True,  True,  True,  True,  True,  True,  True,  True,  True,
True,  True,  True,  True,  True,  True, False, False, False,
False, False, False, False, False, False, False, False, False,
False, False, False, False, False, False, False, False, False,
False, False, False, False, False, False, False, False, False,
False, False, False, False, False])
In [15]:
def eval_step_input(Fo, to, time):
F = np.empty_like(time)
for i, ti in enumerate(time):
if ti < to:
F[i] = 0.0
else:
F[i] = Fo
return F

In [16]:
plt.plot(ts, eval_step_input(5.0, 3.0, ts))

Out[16]:
[<matplotlib.lines.Line2D at 0x7fedc991ad68>]
In [17]:
eval_step_input(5.0, 3.0, ts)

Out[17]:
array([0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 5., 5.,
5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5.,
5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5.])
In [18]:
eval_step_input(5.0, 3.0, 7.0)

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-18-a99dc016e877> in <module>()
----> 1 eval_step_input(5.0, 3.0, 7.0)

<ipython-input-15-9792dfd92762> in eval_step_input(Fo, to, time)
1 def eval_step_input(Fo, to, time):
2     F = np.empty_like(time)
----> 3     for i, ti in enumerate(time):
4         if ti < to:
5             F[i] = 0.0

TypeError: 'float' object is not iterable
In [19]:
def eval_step_input(Fo, to, time):
if np.isscalar(time):
if time < to:
return 0.0
else:
return Fo
else:
F = np.empty_like(time)
for i, ti in enumerate(time):
if ti < to:
F[i] = 0.0
else:
F[i] = Fo
return F

In [20]:
eval_step_input(5.0, 3.0, 7.0)

Out[20]:
5.0
In [21]:
eval_step_input(5.0, 3.0, ts)

Out[21]:
array([0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 5., 5.,
5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5.,
5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5.])
In [22]:
True * 5.0

Out[22]:
5.0
In [23]:
False * 5.0

Out[23]:
0.0
In [24]:
(ts >= 3.0)*5.0

Out[24]:
array([0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 5., 5.,
5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5.,
5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5.])
In [25]:
(5.0 >= 3.0)*5.0

Out[25]:
5.0
In [26]:
def eval_step_input(Fo, to, time):
return (time >=to)*Fo

In [27]:
eval_step_input(5.0, 3.0, ts)

Out[27]:
array([0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 5., 5.,
5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5.,
5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5., 5.])
In [28]:
eval_step_input(5.0, 3.0, 7.0)

Out[28]:
5.0
In [29]:
sys.add_measurement('F', eval_step_input)

In [30]:
sys.diff_eq_func = eval_eom

In [31]:
traj = sys.free_response(20.0)

In [32]:
traj.plot(subplots=True)

Out[32]:
array([<matplotlib.axes._subplots.AxesSubplot object at 0x7fedc9891518>,
<matplotlib.axes._subplots.AxesSubplot object at 0x7fedc97e86a0>,
<matplotlib.axes._subplots.AxesSubplot object at 0x7fedc97f60f0>,
<matplotlib.axes._subplots.AxesSubplot object at 0x7fedc97bc518>],
dtype=object)
In [33]:
def eval_ramp_input(Ft, to, time):
return (time >= to)*(Ft*time - Ft*to)

In [34]:
del sys.measurements['F']

In [35]:
sys.add_measurement('F', eval_ramp_input)

In [36]:
sys.measurements

Out[36]:
{'F': -0.0}
In [37]:
traj = sys.free_response(20.0)

In [38]:
traj.plot(subplots=True)

Out[38]:
array([<matplotlib.axes._subplots.AxesSubplot object at 0x7fedc96a6828>,
<matplotlib.axes._subplots.AxesSubplot object at 0x7fedc9660fd0>,
<matplotlib.axes._subplots.AxesSubplot object at 0x7fedc969b470>,
<matplotlib.axes._subplots.AxesSubplot object at 0x7fedc96520f0>],
dtype=object)