models Module

class Meijaard2007Model(parameter_set)

Bases: _Model

Carvallo-Whipple model presented in [Meijaard2007]. It is both linear and the minimal model in terms of states and coordinates that fully describe the vehicle’s dynamics: self-stability and non-minimum phase behavior.

Parameters:
parameter_setParameterSet

The paramter_set.to_parameterization('meijaard2007') must return a dictionary that maps floats to the parameter keys containing:

  • IBxx : x moment of inertia of the frame/rider [kg*m**2]

  • IBxz : xz product of inertia of the frame/rider [kg*m**2]

  • IBzz : z moment of inertia of the frame/rider [kg*m**2]

  • IFxx : x moment of inertia of the front wheel [kg*m**2]

  • IFyy : y moment of inertia of the front wheel [kg*m**2]

  • IHxx : x moment of inertia of the handlebar/fork [kg*m**2]

  • IHxz : xz product of inertia of the handlebar/fork [kg*m**2]

  • IHzz : z moment of inertia of the handlebar/fork [kg*m**2]

  • IRxx : x moment of inertia of the rear wheel [kg*m**2]

  • IRyy : y moment of inertia of the rear wheel [kg*m**2]

  • c : trail [m]

  • g : acceleration due to gravity [m/s**2]

  • lam : steer axis tilt [rad]

  • mB : frame/rider mass [kg]

  • mF : front wheel mass [kg]

  • mH : handlebar/fork assembly mass [kg]

  • mR : rear wheel mass [kg]

  • rF : front wheel radius [m]

  • rR : rear wheel radius [m]

  • v : speed [m/s]

  • w : wheelbase [m]

  • xB : x distance to the frame/rider center of mass [m]

  • xH : x distance to the frame/rider center of mass [m]

  • zB : z distance to the frame/rider center of mass [m]

  • zH : z distance to the frame/rider center of mass [m]

Attributes:
input_varslist of strings

Ordered list of ASCII strings that name the model’s input variables.

state_varslist of strings

Ordered list of ASCII strings that name the model’s state variables.

input_vars_latexlist of raw strings

Ordered list of LaTeX strings that name the model’s input variables.

state_vars_latexlist of raw strings

Ordered list of LaTeX strings that name the model’s state variables.

References

[Meijaard2007]

Meijaard J.P, Papadopoulos Jim M, Ruina Andy and Schwab A.L, 2007, Linearized dynamics equations for the balance and steer of a bicycle: a benchmark and review, Proc. R. Soc. A., 463:1955–1982 http://doi.org/10.1098/rspa.2007.1857

calc_eigen(left=False, **parameter_overrides)

Returns the right (or left) eigenvalues and eigenvectors of the linear model.

Parameters:
leftboolean, optional

If true, the left eigenvectors will be returned, i.e. A.T*v=lam*v.

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
evalsndarray, shape(4,) or shape (n,4)

Eigenvalues.

evecsndarray, shape(4,4) or shape (n,4,4)

Eigenvectors, each columns are eigenvectors and are associated with same index of the eigenvalues.

Examples

>>> from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
>>> from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
>>> from bicycleparameters.models import Meijaard2007Model
>>> p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
>>> m = Meijaard2007Model(p)
>>> evals, evecs = m.calc_eigen()
>>> evals
array([-6.74423162+0.j        , -2.9146438 +0.j        ,
        3.39244999+0.61085077j,  3.39244999-0.61085077j])
>>> evecs
array([[ 0.00197344+0.j        , -0.2953538 +0.j        ,
         0.04320146-0.0753826j ,  0.04320146+0.0753826j ],
       [ 0.14665803+0.j        ,  0.13447333+0.j        ,
        -0.26053575+0.04691255j, -0.26053575-0.04691255j],
       [-0.01330934+0.j        ,  0.86085111+0.j        ,
         0.1926063 -0.22934205j,  0.1926063 +0.22934205j],
       [-0.98909574+0.j        , -0.39194186+0.j        ,
        -0.91251108+0.j        , -0.91251108-0.j        ]])
form_reduced_canonical_matrices(**parameter_overrides)

Returns the canonical speed and gravity independent matrices for the Whipple-Carvallo bicycle model linearized about the nominal upright configuration.

Parameters:
**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
Mndarray, shape(2,2) or shape(n,2,2)

Mass matrix.

C1ndarray, shape(2,2) or shape(n,2,2)

Velocity independent damping matrix.

K0ndarray, shape(2,2) or shape(n,2,2)

Gravity independent part of the stiffness matrix.

K2ndarray, shape(2,2) or shape(n,2,2)

Velocity squared independent part of the stiffness matrix.

Notes

The canonical matrices complete the following equation:

M*q'' + v*C1*q' + [g*K0 + v**2*K2]*q = f

where:

  • q = [phi, delta]

  • f = [Tphi, Tdelta]

phi

Bicycle roll angle.

delta

Steer angle.

Tphi

Roll torque.

Tdelta

Steer torque.

v

Bicylce longitudinal speed.

g

Acceleration due to gravity.

Examples

>>> from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
>>> from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
>>> from bicycleparameters.models import Meijaard2007Model
>>> p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
>>> m = Meijaard2007Model(p)
>>> M, C1, K0, K2 = m.form_reduced_canonical_matrices()
>>> M
array([[102.78013216,   1.53582801],
       [  1.53582801,   0.24890226]])
>>> C1
array([[ 0.       , 26.3947333],
       [-0.4503006,  1.037066 ]])
>>> K0
array([[-89.32195981,  -1.74159477],
       [ -1.74159477,  -0.67769624]])
>>> K2
array([[ 0.        , 74.12543   ],
       [ 0.        ,  1.57021553]])
>>> M, _, _, _ = m.form_reduced_canonical_matrices(mB=150.0)
>>> M
array([[176.52178763,   2.69074048],
       [  2.69074048,   0.26699004]])
form_state_space_matrices(**parameter_overrides)

Returns the A and B matrices for the Whipple-Carvallo model linearized about the upright constant velocity configuration.

Parameters:
**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
Andarray, shape(4,4) or shape(n,4,4)

The state matrix.

Bndarray, shape(4,2) or shape(n,4,2)

The input matrix.

Notes

A and B describe the Whipple model in state space form:

x' = A * x + B * u

where the states are:

x = |roll angle | = |phi     |
    |steer angle|   |delta   |
    |roll rate  |   |phidot  |
    |steer rate |   |deltadot|

and the inputs are:

u = |roll torque | = |Tphi  |
    |steer torque|   |Tdelta|

Examples

>>> from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
>>> from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
>>> from bicycleparameters.models import Meijaard2007Model
>>> p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
>>> m = Meijaard2007Model(p)
>>> A, B = m.form_state_space_matrices()
>>> A
array([[ 0.        ,  0.        ,  1.        ,  0.        ],
       [ 0.        ,  0.        ,  0.        ,  1.        ],
       [ 8.26150335, -0.9471634 , -0.02977958, -0.21430735],
       [17.66475151, 26.24590352,  1.99289841, -2.84419587]])
>>> B
array([[ 0.        ,  0.        ],
       [ 0.        ,  0.        ],
       [ 0.01071772, -0.06613267],
       [-0.06613267,  4.42570676]])
input_units = ['N-m', 'N-m']
input_vars = ['Tphi', 'Tdelta']
input_vars_latex = ['T_\\phi', 'T_\\delta']
plot_eigenvalue_parts(ax=None, colors=None, show_stable_regions=True, hide_zeros=False, sort_modes=True, **parameter_overrides)

Returns a matplotlib axis of the real and imaginary parts of the eigenvalues plotted against the provided parameter.

Parameters:
axAxes

Matplotlib axes.

colorssequence, len(4)

Matplotlib colors for the 4 modes.

show_stable_regionsboolean, optional

If true, a grey shaded background will indicate stable regions.

sort_modesboolean, optional

Sorting the modes does not always work so you can disable it.

hide_zerosboolean or float, optional

If true, real or imaginary parts that are smaller than 1e-12 will not be plotted. Providing a float will set the tolerance.

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Examples

import numpy as np
from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
from bicycleparameters.models import Meijaard2007Model
p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
m = Meijaard2007Model(p)
m.plot_eigenvalue_parts(v=np.linspace(0.0, 10.0, num=101))

(Source code, png, hires.png, pdf)

../_images/models-1.png
plot_eigenvectors(hide_zeros=False, **parameter_overrides)

Plots the components of the eigenvectors in the real and imaginary plane.

Parameters:
**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
axesndarray, shape(n, 4)

Polar plot axes for each eigenvector (columns). The rows correspond to a varied parameter.

Examples

import numpy as np
from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
from bicycleparameters.models import Meijaard2007Model
p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
m = Meijaard2007Model(p)
m.plot_eigenvectors(v=[1.0, 3.0, 5.0])

(Source code, png, hires.png, pdf)

../_images/models-2.png
plot_mode_simulations(times, hide_zeros=False, **parameter_overrides)

Returns matplotlib subplot axes with a simulation of each mode.

Parameters:
timesarray_like, shape(n,)

Monotonic increasing time values to simulate over.

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
axesndarray, shape(4,2)

Subplot axes with the modes on the rows and the angles in the first column and the angular rates in the second column.

Examples

import numpy as np
from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
from bicycleparameters.models import Meijaard2007Model
p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
m = Meijaard2007Model(p)
times = np.linspace(0.0, 5.0, num=51)
m.plot_mode_simulations(times, v=6.0)

(Source code, png, hires.png, pdf)

../_images/models-3.png
plot_simulation(times, initial_conditions, input_func=None, **parameter_overrides)

Returns the state and input trajectories at each time value.

Parameters:
timesarray_like, shape(n,)

Monotonic increasing time values to simulate over.

initial_conditionsarray_like, shape(4,)

Initial values of the states.

input_funcfunction

Takes form u = f(t, x) where u is array_like, shape(2,).

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
axesndarray, shape(3,)

Three subplots that plot the input trajectories, state angle trajectories, and state angular rates.

Examples

import numpy as np
from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
from bicycleparameters.models import Meijaard2007Model
p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
m = Meijaard2007Model(p)
times = np.linspace(0.0, 5.0, num=51)
x0 = np.deg2rad([10.0, 5.0, 0.0, 0.0])
m.plot_simulation(times, x0, v=6.0)

(Source code, png, hires.png, pdf)

../_images/models-4.png
simulate(times, initial_conditions, input_func=None, **parameter_overrides)

Returns the state and input trajectories at each time value.

Parameters:
timesarray_like, shape(n,)

Monotonic increasing time values to simulate over.

initial_conditionsarray_like, shape(4,)

Initial values of the states.

input_funcfunction

Takes form u = f(t, x) where u is array_like, shape(2,).

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
statesndarray, shape(n, 4)

State trajectories over n time values.

inputsndatrray, shape(n, 2)

Input trajectories over n time values.

simulate_modes(times, **parameter_overrides)

Returns simulation results showing the behavior of each eigenmode.

Parameters:
timesarray_like, shape(n,)

Monotonic increasing time values to simulate over.

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
resultsndarray, shape(4, n, 4)

State trajectories for each mode with the shape corresponding to (mode, time, state).

state_units = ['rad', 'rad', 'rad/s', 'rad/s']
state_vars = ['phi', 'delta', 'phidot', 'deltadot']
state_vars_latex = ['\\phi', '\\delta', '\\dot{\\phi}', '\\dot{\\delta}']
class Meijaard2007WithFeedbackModel(parameter_set)

Bases: Meijaard2007Model

Linear Carvallo-Whipple bicycle model that includes full state feedback to drive all states to zero. With two inputs (roll torque and steer torque) and four states (roll angle, steer angle, roll rate, steer rate) there are eight control gain parameters in addition to the parameters defined in [Meijaard2007].

The states are:

x = |roll angle         | = |phi     |
    |steer angle        |   |delta   |
    |roll angular rate  |   |phidot  |
    |steer angular rate |   |deltadot|

The inputs are:

u = |roll torque | = |Tphi  |
    |steer torque|   |Tdelta|

Applying full state feedback gives this controller:

u = -K*x = -|kTphi_phi, kTphi_del, kTphi_phid, kTphi_deld|*|phi     |
            |kTdel_phi, kTdel_del, kTdel_phid, kTdel_deld| |delta   |
                                                           |phidot  |
                                                           |deltadot|

This represents the new model:

x' = (A - B*K)*x + B*u

so steer and roll torque can be applied in parallel to the feedback control.

calc_eigen(left=False, **parameter_overrides)

Returns the right (or left) eigenvalues and eigenvectors of the linear model.

Parameters:
leftboolean, optional

If true, the left eigenvectors will be returned, i.e. A.T*v=lam*v.

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
evalsndarray, shape(4,) or shape (n,4)

Eigenvalues.

evecsndarray, shape(4,4) or shape (n,4,4)

Eigenvectors, each columns are eigenvectors and are associated with same index of the eigenvalues.

Examples

>>> from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
>>> from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
>>> from bicycleparameters.models import Meijaard2007Model
>>> p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
>>> m = Meijaard2007Model(p)
>>> evals, evecs = m.calc_eigen()
>>> evals
array([-6.74423162+0.j        , -2.9146438 +0.j        ,
        3.39244999+0.61085077j,  3.39244999-0.61085077j])
>>> evecs
array([[ 0.00197344+0.j        , -0.2953538 +0.j        ,
         0.04320146-0.0753826j ,  0.04320146+0.0753826j ],
       [ 0.14665803+0.j        ,  0.13447333+0.j        ,
        -0.26053575+0.04691255j, -0.26053575-0.04691255j],
       [-0.01330934+0.j        ,  0.86085111+0.j        ,
         0.1926063 -0.22934205j,  0.1926063 +0.22934205j],
       [-0.98909574+0.j        , -0.39194186+0.j        ,
        -0.91251108+0.j        , -0.91251108-0.j        ]])
form_reduced_canonical_matrices(**parameter_overrides)

Returns the canonical speed and gravity independent matrices for the Whipple-Carvallo bicycle model linearized about the nominal upright configuration.

Parameters:
**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
Mndarray, shape(2,2) or shape(n,2,2)

Mass matrix.

C1ndarray, shape(2,2) or shape(n,2,2)

Velocity independent damping matrix.

K0ndarray, shape(2,2) or shape(n,2,2)

Gravity independent part of the stiffness matrix.

K2ndarray, shape(2,2) or shape(n,2,2)

Velocity squared independent part of the stiffness matrix.

Notes

The canonical matrices complete the following equation:

M*q'' + v*C1*q' + [g*K0 + v**2*K2]*q = f

where:

  • q = [phi, delta]

  • f = [Tphi, Tdelta]

phi

Bicycle roll angle.

delta

Steer angle.

Tphi

Roll torque.

Tdelta

Steer torque.

v

Bicylce longitudinal speed.

g

Acceleration due to gravity.

Examples

>>> from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
>>> from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
>>> from bicycleparameters.models import Meijaard2007Model
>>> p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
>>> m = Meijaard2007Model(p)
>>> M, C1, K0, K2 = m.form_reduced_canonical_matrices()
>>> M
array([[102.78013216,   1.53582801],
       [  1.53582801,   0.24890226]])
>>> C1
array([[ 0.       , 26.3947333],
       [-0.4503006,  1.037066 ]])
>>> K0
array([[-89.32195981,  -1.74159477],
       [ -1.74159477,  -0.67769624]])
>>> K2
array([[ 0.        , 74.12543   ],
       [ 0.        ,  1.57021553]])
>>> M, _, _, _ = m.form_reduced_canonical_matrices(mB=150.0)
>>> M
array([[176.52178763,   2.69074048],
       [  2.69074048,   0.26699004]])
form_state_space_matrices(**parameter_overrides)

Returns the A and B matrices for the Carvallo-Whipple model linearized about the upright constant velocity configuration with a full state feedback steer controller to drive the states to zero.

Returns:
Andarray, shape(4,4) or shape(n,4,4)

The state matrix.

Bndarray, shape(4,2) or shape(n,4,2)

The input matrix.

Notes

A, B, and K describe the model in state space form:

x' = (A - B*K)*x + B*u

where:

x = |phi     | = |roll angle         |
    |delta   |   |steer angle        |
    |phidot  |   |roll angular rate  |
    |deldot  |   |steer angular rate |

K = | kTphi_phi kTphi_del kTphi_phid kTphi_deld |
    | kTdel_phi kTdel_del kTdel_phid kTdel_deld |

u = |Tphi  | = |roll torque |
    |Tdelta|   |steer torque|
input_units = ['N-m', 'N-m']
input_vars = ['Tphi', 'Tdelta']
input_vars_latex = ['T_\\phi', 'T_\\delta']
plot_eigenvalue_parts(ax=None, colors=None, show_stable_regions=True, hide_zeros=False, sort_modes=True, **parameter_overrides)

Returns a matplotlib axis of the real and imaginary parts of the eigenvalues plotted against the provided parameter.

Parameters:
axAxes

Matplotlib axes.

colorssequence, len(4)

Matplotlib colors for the 4 modes.

show_stable_regionsboolean, optional

If true, a grey shaded background will indicate stable regions.

sort_modesboolean, optional

Sorting the modes does not always work so you can disable it.

hide_zerosboolean or float, optional

If true, real or imaginary parts that are smaller than 1e-12 will not be plotted. Providing a float will set the tolerance.

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Examples

import numpy as np
from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
from bicycleparameters.models import Meijaard2007Model
p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
m = Meijaard2007Model(p)
m.plot_eigenvalue_parts(v=np.linspace(0.0, 10.0, num=101))

(Source code, png, hires.png, pdf)

../_images/models-5.png
plot_eigenvectors(hide_zeros=False, **parameter_overrides)

Plots the components of the eigenvectors in the real and imaginary plane.

Parameters:
**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
axesndarray, shape(n, 4)

Polar plot axes for each eigenvector (columns). The rows correspond to a varied parameter.

Examples

import numpy as np
from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
from bicycleparameters.models import Meijaard2007Model
p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
m = Meijaard2007Model(p)
m.plot_eigenvectors(v=[1.0, 3.0, 5.0])

(Source code, png, hires.png, pdf)

../_images/models-6.png
plot_gains(axes=None, **parameter_overrides)

Plots the gains versus a single varying parameter. The parameter_overrides should contain one parameter that is an array, other than the eight gains. That parameter will be used for the x axis. The gains can be either arrays of the same length or scalars.

Parameters:
axesarray_like, shape(2, 4)

Matplotlib axes set to plot to.

parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
axesndarray, shape(2, 4)

Array of matplotlib axes.

Examples

import numpy as np
from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
from bicycleparameters.models import Meijaard2007WithFeedbackModel
p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
m = Meijaard2007WithFeedbackModel(p)
m.plot_gains(v=np.linspace(0.0, 10.0, num=101),
             kTdel_phid=-10.0*np.linspace(0.0, 5.0, num=101))

(Source code, png, hires.png, pdf)

../_images/models-7.png
plot_mode_simulations(times, hide_zeros=False, **parameter_overrides)

Returns matplotlib subplot axes with a simulation of each mode.

Parameters:
timesarray_like, shape(n,)

Monotonic increasing time values to simulate over.

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
axesndarray, shape(4,2)

Subplot axes with the modes on the rows and the angles in the first column and the angular rates in the second column.

Examples

import numpy as np
from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
from bicycleparameters.models import Meijaard2007Model
p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
m = Meijaard2007Model(p)
times = np.linspace(0.0, 5.0, num=51)
m.plot_mode_simulations(times, v=6.0)

(Source code, png, hires.png, pdf)

../_images/models-8.png
plot_simulation(times, initial_conditions, input_func=None, **parameter_overrides)

Returns the state and input trajectories at each time value.

Parameters:
timesarray_like, shape(n,)

Monotonic increasing time values to simulate over.

initial_conditionsarray_like, shape(4,)

Initial values of the states.

input_funcfunction

Takes form u = f(t, x) where u is array_like, shape(2,).

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
axesndarray, shape(3,)

Three subplots that plot the input trajectories, state angle trajectories, and state angular rates.

Examples

import numpy as np
from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
from bicycleparameters.models import Meijaard2007Model
p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
m = Meijaard2007Model(p)
times = np.linspace(0.0, 5.0, num=51)
x0 = np.deg2rad([10.0, 5.0, 0.0, 0.0])
m.plot_simulation(times, x0, v=6.0)

(Source code, png, hires.png, pdf)

../_images/models-9.png
simulate(times, initial_conditions, input_func=None, **parameter_overrides)

Returns the state and input trajectories at each time value.

Parameters:
timesarray_like, shape(n,)

Monotonic increasing time values to simulate over.

initial_conditionsarray_like, shape(4,)

Initial values of the states.

input_funcfunction

Takes form u = f(t, x) where u is array_like, shape(2,).

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
statesndarray, shape(n, 4)

State trajectories over n time values.

inputsndatrray, shape(n, 2)

Input trajectories over n time values.

simulate_modes(times, **parameter_overrides)

Returns simulation results showing the behavior of each eigenmode.

Parameters:
timesarray_like, shape(n,)

Monotonic increasing time values to simulate over.

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
resultsndarray, shape(4, n, 4)

State trajectories for each mode with the shape corresponding to (mode, time, state).

state_units = ['rad', 'rad', 'rad/s', 'rad/s']
state_vars = ['phi', 'delta', 'phidot', 'deltadot']
state_vars_latex = ['\\phi', '\\delta', '\\dot{\\phi}', '\\dot{\\delta}']
class MooreRiderLean2012Model(parameter_set)

Bases: _Model

Carvallo-Whipple model with an extra rigid body that represents the rider’s upper body which can lean relative to the bicycle’s rear frame.

This specific model is described in [Moore2012] and the equations of motion linearized about the nominal upright configuration at a constant longitudinal speed were extracted from the source code from the dissertation and integrated here.

12 States:

  • \(q_1\), 0: First Cartesian coordinate in the ground plane to rear wheel contact [m]

  • \(q_2\), 1: Second Cartesian coordinate in the ground plane to rear wheel contact [m]

  • \(q_3\), 2: rear frame yaw angle [rad]

  • \(q_4\), 3: rear frame roll angle [rad]

  • \(q_6\), 4: rear wheel angle relative to rear frame [rad]

  • \(q_7\), 5: steer angle (front frame relative to rear frame) [rad]

  • \(q_8\), 6: front wheel angle relative to front frame [rad]

  • \(q_9\), 7: rider lean angle relative to rear frame [rad]

  • \(u_4\), 8: rear frame roll rate [rad/s]

  • \(u_6\): rear wheel angular rate relative to rear frame [rad/s]

  • \(u_7\), 9: steer angular rate relative to the rear frame [rad/s]

  • \(u_9\), 10: rider lean angular rate relative to the rear frame [rad/s]

3 Inputs:

  • \(T_4\), 0: roll torque (between rear frame and ground) [Nm]

  • \(T_7\), 1: steer torque (between rear and front frames [Nm]

  • \(T_9\), 2: rider lean torque (between rider upper body and rear frame) [Nm]

calc_eigen(left=False, **parameter_overrides)

Returns the right (or left) eigenvalues and eigenvectors of the linear model.

Parameters:
leftboolean, optional

If true, the left eigenvectors will be returned, i.e. A.T*v=lam*v.

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
evalsndarray, shape(4,) or shape (n,4)

Eigenvalues.

evecsndarray, shape(4,4) or shape (n,4,4)

Eigenvectors, each columns are eigenvectors and are associated with same index of the eigenvalues.

Examples

>>> from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
>>> from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
>>> from bicycleparameters.models import Meijaard2007Model
>>> p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
>>> m = Meijaard2007Model(p)
>>> evals, evecs = m.calc_eigen()
>>> evals
array([-6.74423162+0.j        , -2.9146438 +0.j        ,
        3.39244999+0.61085077j,  3.39244999-0.61085077j])
>>> evecs
array([[ 0.00197344+0.j        , -0.2953538 +0.j        ,
         0.04320146-0.0753826j ,  0.04320146+0.0753826j ],
       [ 0.14665803+0.j        ,  0.13447333+0.j        ,
        -0.26053575+0.04691255j, -0.26053575-0.04691255j],
       [-0.01330934+0.j        ,  0.86085111+0.j        ,
         0.1926063 -0.22934205j,  0.1926063 +0.22934205j],
       [-0.98909574+0.j        , -0.39194186+0.j        ,
        -0.91251108+0.j        , -0.91251108-0.j        ]])
form_state_space_matrices(**parameter_overrides)

Returns the A and B matrices for the model linearized about the upright constant longitudinal speed configuration.

Parameters:
**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
Andarray, shape(12,12) or shape(n,12,12)

The state matrix.

Bndarray, shape(12,3) or shape(n,12,3)

The input matrix.

Notes

A and B describe the model in state space form:

x' = A * x + B * u

Examples

>>> from bicycleparameters.parameter_dicts import mooreriderlean2012_browser_jason
>>> from bicycleparameters.parameter_sets import MooreRiderLean2012ParameterSet
>>> from bicycleparameters.models import MooreRiderLean2012Model
>>> p = MooreRiderLean2012ParameterSet(mooreriderlean2012_browser_jason)
>>> m = MooreRiderLean2012Model(p)
>>> A, B = m.form_state_space_matrices(v=4.0)
>>> A
array([[ 0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
         0.00000000e+00, -8.26759379e-33,  0.00000000e+00,
         0.00000000e+00,  0.00000000e+00, -0.00000000e+00,
        -3.40958859e-01, -0.00000000e+00,  0.00000000e+00],
       [ 0.00000000e+00,  0.00000000e+00,  4.00000000e+00,
         0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
         0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
        -0.00000000e+00, -0.00000000e+00,  0.00000000e+00],
       [ 0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
         0.00000000e+00,  0.00000000e+00,  3.28701304e+00,
         0.00000000e+00,  0.00000000e+00, -0.00000000e+00,
         0.00000000e+00,  5.63565404e-02,  0.00000000e+00],
       [ 0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
         0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
         0.00000000e+00,  0.00000000e+00,  1.00000000e+00,
         0.00000000e+00,  0.00000000e+00,  0.00000000e+00],
       [ 0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
         0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
         0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
         1.00000000e+00,  0.00000000e+00,  0.00000000e+00],
       [ 0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
         0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
         0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
         0.00000000e+00,  1.00000000e+00,  0.00000000e+00],
       [ 0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
        -0.00000000e+00,  3.30703751e-32,  0.00000000e+00,
         0.00000000e+00,  0.00000000e+00, -0.00000000e+00,
         9.92516037e-01, -0.00000000e+00,  0.00000000e+00],
       [ 0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
         0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
         0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
         0.00000000e+00,  0.00000000e+00,  1.00000000e+00],
       [ 0.00000000e+00,  0.00000000e+00, -0.00000000e+00,
         1.11095692e+01,  0.00000000e+00, -1.52341011e+01,
        -0.00000000e+00, -8.46683843e+00, -2.04638864e-01,
         0.00000000e+00, -1.34590803e+00,  0.00000000e+00],
       [ 0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
         0.00000000e+00, -0.00000000e+00,  0.00000000e+00,
         0.00000000e+00,  0.00000000e+00, -0.00000000e+00,
        -0.00000000e+00, -0.00000000e+00, -0.00000000e+00],
       [ 0.00000000e+00,  0.00000000e+00, -0.00000000e+00,
         1.04939762e+01,  0.00000000e+00,  5.61738098e+00,
        -0.00000000e+00,  1.68341653e+01,  8.26892076e+00,
         0.00000000e+00, -1.01980345e+01,  0.00000000e+00],
       [ 0.00000000e+00,  0.00000000e+00,  0.00000000e+00,
        -1.47851073e+01, -0.00000000e+00,  2.03755703e+01,
         0.00000000e+00,  5.11495540e+01,  4.37754540e-01,
        -0.00000000e+00,  2.54094086e+00, -0.00000000e+00]])
>>> B
array([[ 0.        ,  0.        ,  0.        ],
       [ 0.        ,  0.        ,  0.        ],
       [ 0.        ,  0.        ,  0.        ],
       [ 0.        ,  0.        ,  0.        ],
       [ 0.        ,  0.        ,  0.        ],
       [ 0.        ,  0.        ,  0.        ],
       [ 0.        ,  0.        ,  0.        ],
       [ 0.        ,  0.        ,  0.        ],
       [ 0.02883808, -0.11361237, -0.09393174],
       [ 0.        , -0.        , -0.        ],
       [-0.11361237,  4.59077823,  0.24303462],
       [-0.09393174,  0.24303462,  0.48717313]])
input_units = ['N-m', 'N-m', 'N-m']
input_vars = ['T4', 'T7', 'T9']
input_vars_latex = ['T_4', 'T_7', 'T_9']
plot_eigenvalue_parts(ax=None, colors=None, show_stable_regions=True, hide_zeros=False, sort_modes=True, **parameter_overrides)

Returns a matplotlib axis of the real and imaginary parts of the eigenvalues plotted against the provided parameter.

Parameters:
axAxes

Matplotlib axes.

colorssequence, len(4)

Matplotlib colors for the 4 modes.

show_stable_regionsboolean, optional

If true, a grey shaded background will indicate stable regions.

sort_modesboolean, optional

Sorting the modes does not always work so you can disable it.

hide_zerosboolean or float, optional

If true, real or imaginary parts that are smaller than 1e-12 will not be plotted. Providing a float will set the tolerance.

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Examples

import numpy as np
from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
from bicycleparameters.models import Meijaard2007Model
p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
m = Meijaard2007Model(p)
m.plot_eigenvalue_parts(v=np.linspace(0.0, 10.0, num=101))

(Source code, png, hires.png, pdf)

../_images/models-10.png
plot_eigenvectors(hide_zeros=False, **parameter_overrides)

Plots the components of the eigenvectors in the real and imaginary plane.

Parameters:
**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
axesndarray, shape(n, 4)

Polar plot axes for each eigenvector (columns). The rows correspond to a varied parameter.

Examples

import numpy as np
from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
from bicycleparameters.models import Meijaard2007Model
p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
m = Meijaard2007Model(p)
m.plot_eigenvectors(v=[1.0, 3.0, 5.0])

(Source code, png, hires.png, pdf)

../_images/models-11.png
plot_mode_simulations(times, hide_zeros=False, **parameter_overrides)

Returns matplotlib subplot axes with a simulation of each mode.

Parameters:
timesarray_like, shape(n,)

Monotonic increasing time values to simulate over.

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
axesndarray, shape(4,2)

Subplot axes with the modes on the rows and the angles in the first column and the angular rates in the second column.

Examples

import numpy as np
from bicycleparameters.parameter_dicts import meijaard2007_browser_jason
from bicycleparameters.parameter_sets import Meijaard2007ParameterSet
from bicycleparameters.models import Meijaard2007Model
p = Meijaard2007ParameterSet(meijaard2007_browser_jason, True)
m = Meijaard2007Model(p)
times = np.linspace(0.0, 5.0, num=51)
m.plot_mode_simulations(times, v=6.0)

(Source code, png, hires.png, pdf)

../_images/models-12.png
plot_simulation(times, initial_conditions, input_func=None, **parameter_overrides)

Returns the state and input trajectories at each provided time value.

Parameters:
timesarray_like, shape(m,)

Monotonic increasing time values to simulate over.

initial_conditionsarray_like, shape(12,)

Initial values of the states: [q1, q2, q3, q4, q6, q7, q8, q9, u4, u6, u7, u9].

input_funcfunction

Takes form r = f(t, x) where r is array_like, shape(4,). with r = [T4, T7, T9].

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
axesndarray, shape(4,)

Four subplots that plot the input trajectories, distance trajectories, state angle trajectories, and state angular rates.

Examples

import numpy as np

from bicycleparameters.parameter_dicts import mooreriderlean2012_browser_jason
from bicycleparameters.parameter_sets import MooreRiderLean2012ParameterSet
from bicycleparameters.models import MooreRiderLean2012Model
p = MooreRiderLean2012ParameterSet(mooreriderlean2012_browser_jason)
m = MooreRiderLean2012Model(p)
times = np.linspace(0.0, 5.0, num=51)
x0 = np.zeros(12)
x0[5] = np.deg2rad(5.0)  # steer angle
m.plot_simulation(times, x0, v=20.0)

(Source code, png, hires.png, pdf)

../_images/models-13.png
simulate(times, initial_conditions, input_func=None, **parameter_overrides)

Returns the state and input trajectories at each time value.

Parameters:
timesarray_like, shape(n,)

Monotonic increasing time values to simulate over.

initial_conditionsarray_like, shape(4,)

Initial values of the states.

input_funcfunction

Takes form u = f(t, x) where u is array_like, shape(2,).

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
statesndarray, shape(n, 4)

State trajectories over n time values.

inputsndatrray, shape(n, 2)

Input trajectories over n time values.

simulate_modes(times, **parameter_overrides)

Returns simulation results showing the behavior of each eigenmode.

Parameters:
timesarray_like, shape(n,)

Monotonic increasing time values to simulate over.

**parameter_overridesdictionary

Parameter keys that map to floats or array_like of floats shape(n,). All keys that map to array_like must be of the same length.

Returns:
resultsndarray, shape(4, n, 4)

State trajectories for each mode with the shape corresponding to (mode, time, state).

state_units = ['m', 'm', 'rad', 'rad', 'rad', 'rad', 'rad', 'rad', 'rad/s', 'rad/s', 'rad/s']
state_vars = ['q1', 'q2', 'q3', 'q4', 'q6', 'q7', 'q8', 'q9', 'u4', 'u6', 'u7', 'u9']
state_vars_latex = ['q_1', 'q_2', 'q_3', 'q_4', 'q_6', 'q_7', 'q_8', 'q_9', 'u_4', 'u_6', 'u_7', 'u_9']