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)
<ipython-input-13-02dad508aefb> in <module>()
      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)