Modules

MMS Module

dyn_sys

Started on Tue Feb 15 17:25:59 2022

@author: Vincent MAHE

Analyse systems of coupled nonlinear equations using the Method of Multiple Scales (MMS). This sub-module defines the dynamical system.

class oscilate.MMS.dyn_sys.Dynamical_system(t, x, Eq, omegas, F=0, fF=None)

Bases: object

The dynamical system studied. See Dynamical systems of interest for a detailed description of the dynamical system.

Parameters:

• t (sympy.Symbol) – time \(t\).

• x (sympy.Function or list of sympy.Function) – Unknown(s) of the problem.

• Eq (sympy.Expr or list of sympy.Expr) – System’s equations without forcing, which can be defined separately (see parameters F and fF). Eq is the unforced system of equations describing the system’s dynamics.

• omegas (sympy.Symbol or list of sympy.Symbol) – The natural frequency of each oscillator.

• F (sympy.Symbol or 0, optional) – Forcing amplitude \(F\). Default is 0.

• fF (sympy.Expr or list of sympy.Expr, optional) – For each oscillator, specify the coefficient multiplying the forcing terms in the equation. It can be used to define parametric forcing. Typically, if the forcing is \(x F \cos(\omega t)\), then fF = x. Default is a list of 1, so the forcing is direct.

complex_coordinates()

Introduce the complex coordinates, denoted \(z_i\) for oscillator \(i\), and defined as

\[\begin{split}\begin{cases} x_i(t) & = z_i(t) + \bar{z}_i(t), \\ \dot{x}_i(t) & = \textrm{j} \omega_{i} (z_i(t) - \bar{z}_i(t)), \end{cases}\end{split}\]

where \(\bar{\bullet}\) denotes the transpose.

complex_form()

Rewrite the dynamical system in complex form, resulting in a system of complex 1st order coupled ODEs from the initial real 2nd order ODEs. See Dynamical systems of interest for details.

class oscilate.MMS.dyn_sys.Forcing(F, fF)

Bases: object

Define the forcing on the system as

• A forcing amplitude F,

• Forcing coefficients fF, used to introduce parametric forcing or simply weight the harmonic forcing.

For the \(i^\textrm{th}\) oscillator, denoting fF[i] as \(f_{F,i}(\boldsymbol{x}(t), \dot{\boldsymbol{x}}(t), \ddot{\boldsymbol{x}}(t))\), the forcing term on that oscillator is \(f_{F,i} F \cos(\omega t)\).

mms

Started on Tue Feb 15 17:25:59 2022

@author: Vincent MAHE

Analyse systems of coupled nonlinear equations using the Method of Multiple Scales (MMS). This sub-module defines the multiple scales system from the dynamical one, and allows to apply the MMS.

oscilate.MMS.mms.Chain_rule_d2fdt2(f, tS, eps)

Apply the chain rule to express second order time derivatives in terms of the time scales’ derivatives.

Parameters:

• f (sympy.Function) – Time scales-dependent function \(f(\boldsymbol{t})\).

• tS (list) – Time scales \(\boldsymbol{t}^\intercal = [t_0, \cdots, t_{N_e}]\).

• eps (sympy.Symbol) – Small parameter \(\epsilon\).

Returns:

d2fdt2 – \(\mathrm{d}^2 f/ \mathrm{d}t^2\) expressed in terms of the time scales.

Return type:

sympy.Function

Notes

Consider a time scales-dependent function \(f(t_0, t_1, ...)\), where \(t_0\) is the fast time and \(t_1, ...\) are the slow times. The Chain Rule is applied to give the expression of \(\mathrm{d}^2 f/ \mathrm{d}t^2\) in terms of the time scales.

oscilate.MMS.mms.Chain_rule_dfdt(f, tS, eps)

Apply the chain rule to express first order time derivatives in terms of the time scales’ derivatives.

Parameters:

• f (sympy.Function) – Time scales-dependent function \(f(\boldsymbol{t})\).

• tS (list) – Time scales \(\boldsymbol{t}^\intercal = [t_0, \cdots, t_{N_e}]\).

• eps (sympy.Symbol) – Small parameter \(\epsilon\).

Returns:

dfdt – \(\mathrm{d} f/ \mathrm{d}t\) expressed in terms of the time scales.

Return type:

sympy.Function

Notes

Consider a time scales-dependent function \(f(t_0, t_1, ...)\), where \(t_0\) is the fast time and \(t_1, ...\) are the slow times. The Chain Rule is applied to give the expression of \(\mathrm{d} f/ \mathrm{d}t\) in terms of the time scales.

class oscilate.MMS.mms.Forcing_MMS(F, f_order, fF, forcing_term)

Bases: object

Define the forcing on the system as

• A forcing amplitude F

• A scaling order f_order for the forcing

• Forcing coefficients fF

• Forcing terms (direct or parametric) forcing_term

class oscilate.MMS.mms.Multiple_scales_system(dynamical_system, eps, Ne, omega_ref, sub_scaling, ratio_omegaMMS=1, eps_pow_0=0, ratio_omega_osc=None, detunings=0)

Bases: object

The multiple scales system. See Application of the MMS for a detailed description of the dynamical system.

Parameters:

• dynamical_system (Dynamical_system) – The dynamical system.

• eps (sympy.Symbol) – Small perturbation parameter \(\epsilon\).

• Ne (int) – Truncation order of the asymptotic series and order of the slowest time scale.

• omega_ref (sympy.Symbol) – Reference frequency \(\omega_{\textrm{ref}}\) of the MMS. Not necessarily the frequency around which the MMS is going to be applied, see ratio_omegaMMS.

• sub_scaling (list of tuples) – Substitutions to do to scale the equations. Links small parameters to their scaled counterpart through \(\epsilon\).

• ratio_omegaMMS (int or sympy.Rational, optional) –

  Specify the frequency omegaMMS around which the MMS is going to be applied in terms of \(\omega_{\textrm{ref}}\). Denoting ratio_omegaMMS as \(r_{\textrm{MMS}}\), this means that

  \[\omega_{\textrm{MMS}} = r_{\textrm{MMS}} \omega_{\textrm{ref}}.\]

  Use ratio_omegaMMS=Rational(p,q) for

  \[q \omega_{\textrm{MMS}} = p \omega_{\textrm{ref}}\]

  to get better-looking results than the float \(p/q\). Default is 1.

• eps_pow_0 (int, optional) –

  Order of the leading order term in the asymptotic series of each oscillators’ response. For the \(i^{\textrm{th}}\) oscillator and denoting eps_pow_0 as \(\lambda_0\), this means that

  \[x_i = \epsilon^{\lambda_0} x_{i,0} + \epsilon^{\lambda_0+1} x_{i,1} + \cdots.\]

  Default is 0.

• ratio_omega_osc (list of int or sympy.Rational or None, optional) –

  Specify the natural frequencies of the oscillators \(\omega_i\) in terms of the reference frequency \(\omega_{\textrm{ref}}\). Denoting ratio_omega_osc[i] as \(r_i\), this means that

  \[\omega_i \approx r_i \omega_{\textrm{ref}}.\]

  Use ratio_omega_osc[i]=Rational(p,q) for

  \[q \omega_{i} \approx p \omega_{\textrm{ref}}\]

  to get better-looking results than the float \(p/q\). Default is None for each oscillator, so the \(\omega_i\) are arbitrary and there are no internal resonances. Detuning can be introduced through the detunings keyword argument.

• detunings (list of sympy.Symbol or int, optional) –

  The detuning of each oscillator. Denoting detunings[i] as \(\delta_i\), this means that

  \[\omega_i = r_i \omega_{\textrm{ref}} + \delta_i.\]

  Default is 0 for each oscillator.

absolute_phase()

Introduce an absolute phase \(\psi\) defined as (for oscillator \(i\))

\[\psi = \frac{r_i}{r_{\textrm{MMS}}} \omega t - \beta_i,\]

such that the leading order, homogeneous solution for oscillator \(i\) takes the form

\[x^{\textrm{h}}_{i,0}(t) = a_i(t) \cos(\psi).\]

Then, deduce the solution on \(\psi\) from that on \(\beta\). The introduction of that absolute phase is useful to write the time response in a more compact way, and to derive the instantaneous frequency of oscillation.

apply_MMS_shared()

Apply the MMS. This method contains the operations that are common to the oscillator and complex form. These include a change of phase variables to make the system autonomous, a separation of solvability conditions into real and imaginary parts, and their reconstitution, resulting in the modulation equations. See apply_MMS() or apply_MMS() for details.

autonomous_phases()

Define phase coordinates that render an autonomous system.

Notes

Define new phase coordinates \(\beta_i\) to transform nonautonomous equations into autonomous ones. The \(\beta_i\) are defined as

\[\beta_i = \frac{r_i}{r_{\textrm{MMS}}} \sigma t_1 - \phi_i,\]

where we recall that \(\omega = r_{\textrm{MMS}} \omega_{\textrm{ref}} + \epsilon \sigma\) and \(\omega_{i,0} = r_i \omega_{\textrm{ref}}\). See details on this choice in Application of the MMS.

find_harmonics()

Determine the harmonics contained in the MMS solutions.

Returns:

harmonics – list of the harmonics appearing in the MMS solutions.

Return type:

list

modulation_equations()

Derive the modulation equations of the polar coordinates system.

Notes

Derive the modulation equations of the polar coordinates system (defined in polar_coordinates() and autonomous_phases()) from the secular conditions. For oscillator \(i\) and at order \(j\), these are defined as

\[\begin{split}\begin{cases} \textrm{D}_{j} a_i & = f_{a_i}^{(j)}(\boldsymbol{a}, \boldsymbol{\beta}), \\ a_i \textrm{D}_{j}\beta_i & = f_{\beta_i}^{(j)}(\boldsymbol{a}, \boldsymbol{\beta}), \end{cases}\end{split}\]

where \(\boldsymbol{a}\) and \(\boldsymbol{\beta}\) are vectors containing the polar amplitudes and phases.

The aim here is to compute all the \(f_{a_i}^{(j)}\) and \(f_{\beta_i}^{(j)}\). This is done by:

1. Introducing polar coordinates in the secular terms

2. Splitting the real and imaginary parts of the (complex) secular terms

3. Using the autonomous phase coordinates

4. Collecting the terms governing the slow amplitude and phase dynamics.

oscillators_frequencies()

Gives the expression of every oscillator frequency in terms of the reference frequency, possibly with a detuning.

Notes

For the \(i^\textrm{th}\) oscillator, this leads to

\[\omega_i = r_i \omega_{\textrm{ref}} + \delta_i .\]

An associated first-order natural frequency \(\omega_{i,0}\) is defined by neglecting the detuning \(\delta_i\), which is at least of order \(\epsilon\), resulting in

\[\omega_{i,0} = r_i \omega_{\textrm{ref}}.\]

polar_coordinates()

Define the polar coordinates.

Notes

Define polar coordinates such that, for oscillator \(i\), the complex amplitude of the homogeneous leading order solution is defined as

\[A_i(\boldsymbol{t}_s) = \frac{1}{2} a_i(\boldsymbol{t}_s) e^{\textrm{j} \phi_i(\boldsymbol{t}_s)},\]

where \(a_i(\boldsymbol{t}_s)\) and \(\phi_i(\boldsymbol{t}_s)\) are the solution’s amplitude and phase, respectively.

reconstitution()

Use the reconstitution method to combine the modulation equations at each order. This reconstitution is based on the chain rule relation

\[\dfrac{\textrm{d}(\bullet)}{\textrm{d}t} = \sum_{i=0}^{N_e} \epsilon^{i} \mathrm{D}_i (\bullet) + \mathcal{O}(\epsilon^{N_e+1}).\]

For oscillator \(i\), this results in the reconstituted modulation equations

\[\begin{split}\begin{cases} \dfrac{\textrm{d} a_i}{\textrm{d} t} & = f_{a_i}(\boldsymbol{a}, \boldsymbol{\beta}), \\ a_i \dfrac{\textrm{d} \beta_i}{\textrm{d} t} & = f_{\beta_i}(\boldsymbol{a}, \boldsymbol{\beta}), \end{cases}\end{split}\]

Note that some MMS approaches do not apply this reconstitution step.

sol_x_harmonics()

Extract the harmonics caracteristics from the x solution.

sol_x_polar(orders_polar=0)

Write the solutions using the polar coordinates and \(\cos\) and \(\sin\) functions.

Parameters:

orders_polar (str or int or list of int, optional) – The orders at which the solutions will be rewritten in polar form. If "all", then all solution orders will be rewritten. If int, then only a single order will be rewritten. If list of int, then the listed orders will be rewritten. Default is 0, so only the leading order solution will be rewritten.

solve_instantaneous_frequency()

Compute the instantaneous frequency from the absolute phase through \(\omega_\text{NL} = \dot{\psi}\).

solve_slow_time(IC)

Compute the slow time evolution of the amplitude and phase of an oscillator.

Parameters:

IC (dict(), optional) – Initial conditions on \(a\) and \(\beta\), given as a dict with keys a and beta. The values associated to these keys are dict of the form \(\{\dot{a}(0), a_{\text{IC}} \}\) and \(\{\dot{\beta}(0), \beta_{\text{IC}} \}\) Default is an empty dict, subsequently filled with some symbols for the \(a_{\text{IC}}\) and \(\beta_{\text{IC}}\).

solve_transient(solve_dof=None, IC={})

Solve the transient response of an oscillator.

Parameters:

• solve_dof (None or int, optional) – The oscillator to solve for. If None, no oscillator is solved for. Default is None.

• IC (dict(), optional) – See solve_slow_time().

Notes

Find the transient solution for a given oscillator with the other oscillators’ amplitude set to 0. Setting the other oscillators’ amplitude to zero is done using the method substitution_solve_dof().

solve_x_transient()

Compute the displacement \(x\) associated to a transient trajectory.

substitution_solve_dof(solve_dof)

Set every oscillator amplitude to 0 except the one to solve for.

Notes

If one wants to solve for \(a_i\), then the system is evaluated for \(a_j=0, \; \forall j \neq i\).

time_scales()

Define the time scales.

Notes

The time scales are defined as (see Application of the MMS)

\[t_0 = t, \quad t_1 = \epsilon t, \quad \cdots, \quad t_{N_e} = \epsilon^{N_e} t.\]

Substitutions from the physical time \(t\) to the time scales \(t_i, \; i=0, ..., N_e\) are also prepared. Note that \(t_0\) is refered-to as the fast time as it captures the oscillations. Other time scales are refered-to as slow time as they capture the modulations.

class oscilate.MMS.mms.Sol_MMS

Bases: object

Solutions obtained when applying the MMS.

class oscilate.MMS.mms.Sol_harmonics(x, omega, t)

Bases: object

Harmonic decomposition of the solution.

Parameters:

x (Expr) – The x(t) solution of an oscillator

class oscilate.MMS.mms.Sol_transient

Bases: object

Solutions obtained when applying the MMS and solving for the transient response.

class oscilate.MMS.mms.Substitutions_MMS(sub_t, sub_scaling, sub_omega, sub_sigma)

Bases: object

Substitutions used in the MMS.

oscilate.MMS.mms.cartesian_to_polar(y, sub_polar, sub_phase=None)

Rewrites an expression or a matrix y from cartesian to polar coordinates.

Parameters:

• y (sympy.Expr or sympy.Matrix) – A sympy expression or matrix written in cartesian coordinates.

• sub_polar (list) – A list of substitutions to perform to go from cartesian to polar coordinates.

• sub_phase (list, optional) – Additional substitutions to try and get rid of phases, so that only the amplitude remains in the expression. Default is None.

Returns:

yp – The initial expression or matrix written in polar coordinates.

Return type:

sympy.Expr or sympy.Matrix

oscilate.MMS.mms.get_harmonics(x, omega, t)

Extract the harmonic numbers from x(t).

Parameters:

• x (Expr) – A periodic signal.

• omega (Symbol) – Oscillation frequency.

• t (Symbol) – Time.

Returns:

harmonics – The harmonics contained in x.

Return type:

list

oscilate.MMS.mms.rescale(expr, mms)

Rescales a scaled expression.

Parameters:

• expr (sympy.Expr) – An unscaled expression, i.e. an expression appearing at some order of \(\epsilon\).

• mms (Multiple_scales_system) – The mms system, containing substitutions to scale an expression.

Returns:

expr_scaled – The scaled expression.

Return type:

sympy.Expr

oscilate.MMS.mms.scale_parameters(param, scaling, eps)

Scale parameters with the scaling parameter \(\epsilon\).

Parameters:

• param (list of sympy.Symbol and/or sympy.Function) – Unscaled parameters.

• scaling (list of int or float) – The scaling for each parameter.

• eps (sympy.Symbol) – Small parameter \(\epsilon\).

Returns:

• param_scaled (list of sympy.Symbol and/or sympy.Function) – Scaled parameters.

• sub_scaling (list of 2 lists of tuple) – Substitutions from scaled to unscaled parameters and vice-versa.

  • \(1^{\text{st}}\) list: The substitutions to do to introduce the scaled parameters in an expression.

  • \(2^{\text{nd}}\) list: The substitutions to do to reintroduce the unscaled parameters in a scaled expression.

Notes

For a given parameter \(p\) and a scaling order \(\lambda\), the associated scaled parameter \(\tilde{p}\) is

\[p = \epsilon^{\lambda} \tilde{p} .\]

mms_oscillator

Started on Tue Feb 15 17:25:59 2022

@author: Vincent MAHE

Analyse systems of coupled nonlinear equations using the Method of Multiple Scales (MMS). This sub-module defines the multiple scales system from the dynamical one with an oscillator form, i.e. 2nd order real differential equations. It also carries out the MMS derivations, still in oscillator form.

class oscilate.MMS.mms_oscillator.Coord_MMS_oscillator(mms)

Bases: object

The coordinates used in the MMS.

class oscilate.MMS.mms_oscillator.Multiple_scales_oscillator(dynamical_system, eps, Ne, omega_ref, sub_scaling, ratio_omegaMMS=1, eps_pow_0=0, ratio_omega_osc=None, detunings=0)

Bases: Multiple_scales_system

The multiple scales system in oscillator form, i.e. 2nd order differential equations. See the parent class Multiple_scales_system for a description of the input parameters.

apply_MMS(orders_polar=0)

Apply the MMS.

Parameters:

orders_polar (str, int or list of int, optional) – The orders at which the solutions will be rewritten in polar form. See sol_x_polar().

Notes

The application of the MMS is operated by successively calling the following methods.

1. system_t0(): An equivalent system is written in terms of the fast time scale \(t_0\). This introduces the temporary unknowns \(\tilde{x}_{i,j}(t_0)\) and allows the use of dsolve().

2. sol_order_0(): Leading order solutions \(x_{i,0}(\boldsymbol{t})\) are defined.

3. secular_analysis(): The leading order solutions are introduced in the equations and the secular terms at each order are identified. Cancelling those secular terms is a condition for bounded solutions. It leads to a system of complex modulation equations governing the slow time evolution of the complex amplitude of the homogeneous leading order solutions. Each equation takes the form

   \[\textrm{D}_{j} A_i(\boldsymbol{t}_s) = f_{A_i}^{(j)}(\boldsymbol{A}, t_1).\]

   After cancelling the secular terms the higher order equations are solved successively to express the higher order solutions \(x_{i,j}(\boldsymbol{t}),\; j>0\) in terms of the leading order ones.

4. apply_MMS_shared(): The functions that are shared among various MMS forms (oscillator and complex). This calls

   1. autonomous_phases(): The phase coordinates are changed from \(\phi_i(\boldsymbol{t}_s)\) to \(\beta_i(\boldsymbol{t}_s)\) to eliminate the slow time \(t_1\) from the secular terms. This will be used afterwards to obtain an autonomous system.

   2. modulation_equations(): The secular conditions are split into real and imaginary parts, polar coordinates are used and the autonomous phases are introduced, resulting in an autonomous system of modulation equations on polar coordinates. Equations come by two, one representing the amplitude modulation while the other represents the phase’s, such that

      \[\begin{split}\begin{cases} \textrm{D}_{j} a_i(\boldsymbol{t}_s) & = f_{a_i}^{(j)}(\boldsymbol{a}, \boldsymbol{\beta}), \\ a_i \textrm{D}_{j} \beta_i(\boldsymbol{t}_s) & = f_{\beta_i}^{(j)}(\boldsymbol{a}, \boldsymbol{\beta}). \end{cases}\end{split}\]

      This is the key result of the MMS.

   3. reconstitution(): The modulations on each time scale are now combined to reintroduce the physical time, resulting in a system of the form

      \[\begin{split}\begin{cases} \dfrac{\textrm{d}}{dt} a_i(t) & = f_{a_i}(\boldsymbol{a}, \boldsymbol{\beta}), \\ a_i \dfrac{\textrm{d}}{dt} \beta_i(t) & = f_{\beta_i}(\boldsymbol{a}, \boldsymbol{\beta}). \end{cases}\end{split}\]

      This is known as the reconstitution method.

5. sol_x_polar(): The leading and higher order solutions are rewritten in terms of polar coordinates using \(\cos\) and \(\sin\) functions.

asymptotic_series(dynamical_system, eps_pow_0=0)

Define the asymptotic series on the oscillators’ coordinates x.

Notes

The series expansion for oscillator \(i\) (and for a leading order term \(\epsilon^0 = 1\)) takes the form (see Application of the MMS)

\[x_i(t) = x_{i,0}(t) + \epsilon x_{i,1}(t) + \epsilon^2 x_{i,2}(t) + \cdots + \epsilon^{N_e} x_{i,N_e}(t) + \mathcal{O}(\epsilon^{N_e+1}).\]

On top of introducing the terms of the asymptotic series, this function prepares substitutions from

1. dof \(x_i(t)\) to temporary \(t\)-dependent asymptotic terms \(x_{i,j}(t)\), such that

   \[x_i(t) = x_{i,0}(t) + \epsilon x_{i,1}(t) + \epsilon^2 x_{i,2}(t) + \cdots + \epsilon^{N_e} x_{i,N_e}(t),\]

2. Temporary \(x_{i,j}(t)\) to the time scales-dependent terms \(x_{i,j}(\boldsymbol{t})\), such that

   \[x_i(\boldsymbol{t}) = x_{i,0}(\boldsymbol{t}) + \epsilon x_{i,1}(\boldsymbol{t}) + \epsilon^2 x_{i,2}(\boldsymbol{t}) + \cdots + \epsilon^{N_e} x_{i,N_e}(\boldsymbol{t}).\]

compute_EqO(dynamical_system)

Compute the system of equations for each oscillator at each order of \(\epsilon\). This system is described in Application of the MMS.

The output EqO is a list of lists:

• The \(1^{\text{st}}\) level lists are associated to the equations for each oscillator,

• The \(2^{\text{nd}}\) level lists are associated to the orders of \(\epsilon\) from the lowest to the highest order.

Parameters:

dynamical_system (Dynamical_system) – The dynamical system.

forcing_MMS(dynamical_system)

Rewrite the forcing terms.

Parameters:

dynamical_system (Dynamical_system) – The dynamical system.

Notes

The initial forcing terms are

\[f_{F,i}(\boldsymbol{x}(t), \dot{\boldsymbol{x}}(t), \ddot{\boldsymbol{x}}(t)) F \cos(\omega t), \quad i=1,...,N.\]

Rewritting them involves

1. Replacing the \(x_i(t)\) by their series expansions written in terms of time scales,

2. Scaling the forcing and the parameters in \(f_{F,i}\) if any,

3. Truncating terms whose order is larger than the largest order retained in the MMS,

4. Rewrite \(\cos(\omega t)\) as

   \[\cos(\omega t) = \frac{1}{2} e^{\mathrm{j}(\omega_{\textrm{MMS}} + \epsilon \sigma)t} + cc = \frac{1}{2} e^{\mathrm{j}(\omega_{\textrm{MMS}}t_0 + \sigma t_1)} + cc.\]

secular_analysis()

Identify the secular terms in the MMS equations. This allows to:

1. Compute the modulation equations of the complex amplitudes \(A_i(\boldsymbol{t}_s)\), coming from the elimination of the secular terms,

2. Derive nonsecular MMS equations, i.e. MMS equations with the secular terms cancelled,

3. Use the nonsecular equations to express the higher order solutions \(x_{i,j}(\boldsymbol{t}),\; j>0\) in terms of the \(\boldsymbol{A}(\boldsymbol{t}_s)\).

sol_higher_order(EqO_t0, xO_t0, io, ix)

Compute the higher order solutions \(x_{i,j}(\boldsymbol{t}_s),\; j>0\).

Parameters:

• EqO_t0 (list of list of sympy.Expr) – The MMS equations at each order and for each oscillator written with \(t_0\) as the only independent variable.

• xO_t0 (list of list of sympy.Function) – Oscillators’ response at each order written in terms of \(t_0\) only, \(\tilde{x}_{i,j}(t_0)\).

• io (int) – The order of \(\epsilon\).

• ix (int) – The dof number.

sol_order_0()

Compute the leading-order solutions for each oscillator.

Notes

For oscillator \(i\), the homogeneous solution takes the general form

\[x_{i,0}^{\textrm{h}}(\boldsymbol{t}) = A_i(\boldsymbol{t}_s) e^{\textrm{j} \omega_{i,0} t_0} + cc,\]

where \(A_i(\boldsymbol{t}_s)\) is an unknown complex amplitude.

If the oscillator is subject to hard forcing (i.e. forcing appears at leading order), then the particular solution

\[x_{i,0}^{\textrm{p}}(t_0, t_1) = B_i e^{\textrm{j} \omega t} + cc = B_i e^{\textrm{j} (\omega_{\textrm{MMS}} t_0 + \sigma t_1)} + cc\]

is also taken into account. \(B_i\) is a time-independent function of the forcing parameters.

system_t0()

Rewrite the system with the fast time scale as the only independent variable.

Notes

This is a trick to use dsolve(), which only accepts functions of 1 variable, to solve higher order equations. The higher order equations are rewritten in terms of temporary coordinates \(\tilde{x}_{i,j}(t_0)\) in place of \(x_{i,j}(\boldsymbol{t})\), with \(i,j\) denoting the oscillator number and \(\epsilon\) order, respectively. This is equivalent to temporary considering that \(\boldsymbol{A}(\boldsymbol{t}_s)\) does not depend on the slow times, which is of no consequence as there are no slow time derivatives appearing in the higher order equations at this stage. Indeed, they were either substituted using the complex modulation equations, or they disappeared when eliminating the secular terms.

mms_complex

Created on Mon Apr 20 19:41:59 2026

@author: Vincent MAHE

Analyse systems of coupled nonlinear equations using the Method of Multiple Scales (MMS). This sub-module defines the multiple scales system from the dynamical one with a complex form, i.e. 1st order complex differential equations. It also carries out the MMS derivations, still in complex form, and ends up transforming back to the oscillator form.

class oscilate.MMS.mms_complex.Coord_MMS_complex(mms)

Bases: object

The coordinates used in the MMS.

class oscilate.MMS.mms_complex.Multiple_scales_complex(dynamical_system, eps, Ne, omega_ref, sub_scaling, ratio_omegaMMS=1, eps_pow_0=0, ratio_omega_osc=None, detunings=0)

Bases: Multiple_scales_system

The multiple scales system in complex form, i.e. 1st order complex differential equations. See the parent class Multiple_scales_system for a description of the input parameters.

apply_MMS(orders_polar=0)

Apply the MMS in complex form.

Parameters:

orders_polar (str, int or list of int, optional) – The orders at which the solutions will be rewritten in polar form. See sol_x_polar().

Notes

The application of the MMS is operated by successively calling the following methods.

1. system_t0(): An equivalent system is written in terms of the fast time scale \(t_0\). This introduces the temporary unknowns \(\tilde{z}_{i,j}(t_0)\) and allows the use of dsolve().

2. sol_order_0(): Leading order solutions \(z_{i,0}(\boldsymbol{t})\) are defined.

3. secular_analysis(): The leading order solutions are introduced in the equations and the secular terms at each order are identified. Cancelling those secular terms is a condition for bounded solutions. It leads to a system of complex modulation equations governing the slow time evolution of the complex amplitude of the homogeneous leading order solutions. Each equation takes the form

   \[\textrm{D}_{j} A_i(\boldsymbol{t}_s) = f_{A_i}^{(j)}(\boldsymbol{A}, t_1).\]

   After cancelling the secular terms the higher order equations are solved successively to express the higher order solutions \(z_{i,j}(\boldsymbol{t}),\; j>0\) in terms of the leading order ones.

4. apply_MMS_shared(): The functions that are shared among various MMS forms (oscillator and complex). This calls

   1. autonomous_phases(): The phase coordinates are changed from \(\phi_i(\boldsymbol{t}_s)\) to \(\beta_i(\boldsymbol{t}_s)\) to eliminate the slow time \(t_1\) from the secular terms. This will be used afterwards to obtain an autonomous system.

   2. modulation_equations(): The secular conditions are split into real and imaginary parts, polar coordinates are used and the autonomous phases are introduced, resulting in an autonomous system of modulation equations on polar coordinates. Equations come by two, one representing the amplitude modulation while the other represents the phase’s, such that

      \[\begin{split}\begin{cases} \textrm{D}_{j} a_i(\boldsymbol{t}_s) & = f_{a_i}^{(j)}(\boldsymbol{a}, \boldsymbol{\beta}), \\ a_i \textrm{D}_{j} \beta_i(\boldsymbol{t}_s) & = f_{\beta_i}^{(j)}(\boldsymbol{a}, \boldsymbol{\beta}). \end{cases}\end{split}\]

      This is the key result of the MMS.

   3. reconstitution(): The modulations on each time scale are now combined to reintroduce the physical time, resulting in a system of the form

      \[\begin{split}\begin{cases} \dfrac{\textrm{d}}{dt} a_i(t) & = f_{a_i}(\boldsymbol{a}, \boldsymbol{\beta}), \\ a_i \dfrac{\textrm{d}}{dt} \beta_i(t) & = f_{\beta_i}(\boldsymbol{a}, \boldsymbol{\beta}). \end{cases}\end{split}\]

      This is known as the reconstitution method.

5. sol_x_polar(): The leading and higher order solutions are rewritten in terms of polar coordinates using \(\cos\) and \(\sin\) functions.

asymptotic_series(dynamical_system, eps_pow_0=0)

Define the asymptotic series on the complex oscillators’ coordinates z.

Notes

The series expansion for the complex coordinate of oscillator \(i\) (and for a leading order term \(\epsilon^0 = 1\)) takes the form (see Application of the MMS)

\[z_i(t) = z_{i,0}(t) + \epsilon z_{i,1}(t) + \epsilon^2 z_{i,2}(t) + \cdots + \epsilon^{N_e} z_{i,N_e}(t) + \mathcal{O}(\epsilon^{N_e+1}).\]

On top of introducing the terms of the asymptotic series, this function prepares substitutions from

1. dof \(z_i(t)\) to temporary \(t\)-dependent asymptotic terms \(z_{i,j}(t)\), such that

   \[z_i(t) = z_{i,0}(t) + \epsilon z_{i,1}(t) + \epsilon^2 z_{i,2}(t) + \cdots + \epsilon^{N_e} z_{i,N_e}(t),\]

2. Temporary \(z_{i,j}(t)\) to the time scales-dependent terms \(z_{i,j}(\boldsymbol{t})\), such that

   \[z_i(\boldsymbol{t}) = z_{i,0}(\boldsymbol{t}) + \epsilon z_{i,1}(\boldsymbol{t}) + \epsilon^2 z_{i,2}(\boldsymbol{t}) + \cdots + \epsilon^{N_e} z_{i,N_e}(\boldsymbol{t}).\]

compute_EqO(dynamical_system)

Compute the system of complex equations for each oscillator at each order of \(\epsilon\). This system is described in Application of the MMS.

The output EqO is a list of lists:

• The \(1^{\text{st}}\) level lists are associated to the equations in complex form for each oscillator,

• The \(2^{\text{nd}}\) level lists are associated to the orders of \(\epsilon\) from the lowest to the highest order.

Parameters:

dynamical_system (Dynamical_system) – The dynamical system.

forcing_MMS(dynamical_system)

Rewrite the forcing terms.

Parameters:

dynamical_system (Dynamical_system) – The dynamical system.

Notes

The initial forcing terms are

\[f_{F,i}(\boldsymbol{x}(t), \dot{\boldsymbol{x}}(t), \ddot{\boldsymbol{x}}(t)) F \cos(\omega t), \quad i=1,...,N,\]

and after the complex transformation

The initial forcing terms are

\[f_{F_z,i}(\boldsymbol{z}(t), \dot{\boldsymbol{z}}(t)) F \cos(\omega t), \quad i=1,...,N,\]

where

\[f_{F_z,i}(\boldsymbol{z}(t), \dot{\boldsymbol{z}}(t)) = - \frac{\mathrm{j}}{2 \omega_i} f_{F,i}(\boldsymbol{x}(\boldsymbol{z}(t)), \dot{\boldsymbol{\boldsymbol{z}(t)}}(t), \ddot{\boldsymbol{\boldsymbol{z}(t), \dot{\boldsymbol{z}}(t)}}(t))\]

Rewritting them involves

1. Replacing the \(z_i(t)\) by their series expansions written in terms of time scales,

2. Scaling the forcing and the parameters in \(f_{F_z,i}\) if any,

3. Truncating terms whose order is larger than the largest order retained in the MMS,

4. Rewrite \(\cos(\omega t)\) as

   \[\cos(\omega t) = \frac{1}{2} e^{\mathrm{j}(\omega_{\textrm{MMS}} + \epsilon \sigma)t} + cc = \frac{1}{2} e^{\mathrm{j}(\omega_{\textrm{MMS}}t_0 + \sigma t_1)} + cc.\]

secular_analysis()

Identify the secular terms in the MMS equations. This allows to:

1. Compute the modulation equations of the complex amplitudes \(A_i(\boldsymbol{t}_s)\), coming from the elimination of the secular terms,

2. Derive nonsecular MMS equations, i.e. MMS equations with the secular terms cancelled,

3. Use the nonsecular equations to express the higher order complex solutions \(z_{i,j}(\boldsymbol{t}),\; j>0\) in terms of themath:boldsymbol{A}(boldsymbol{t}_s).

sol_higher_order(EqO_t0, zO_t0, io, ix)

Compute the higher order complex solutions \(z_{i,j}(\boldsymbol{t}_s),\; j>0\).

Parameters:

• EqO_t0 (list of list of sympy.Expr) – The MMS equations at each order and for each oscillator written with \(t_0\) as the only independent variable.

• zO_t0 (list of list of sympy.Function) – Oscillators’ response at each order written in terms of \(t_0\) only, \(\tilde{z}_{i,j}(t_0)\).

• io (int) – The order of \(\epsilon\).

• ix (int) – The dof number.

sol_order_0()

Compute the leading-order solutions for each oscillator.

Notes

For oscillator \(i\), the homogeneous solution takes the general form

\[z_{i,0}^{\textrm{h}}(\boldsymbol{t}) = A_i(\boldsymbol{t}_s) e^{\textrm{j} \omega_{i,0} t_0},\]

where \(A_i(\boldsymbol{t}_s)\) is an unknown complex amplitude.

If the oscillator is subject to hard forcing (i.e. forcing appears at leading order), then the particular solution

\[z_{i,0}^{\textrm{p}}(t_0, t_1) = B_{i,+} e^{\textrm{j} \omega t} + B_{i,-} e^{-\textrm{j} \omega t} = B_{i,+} e^{\textrm{j} (\omega_{\textrm{MMS}} t_0 + \sigma t_1)} + B_{i,-} e^{-\textrm{j} (\omega_{\textrm{MMS}} t_0 + \sigma t_1)}\]

is also taken into account. \(B_{i, \pm}\) are time-independent functions of the forcing parameters.

sol_xO()

Compute the x solutions at each order from the z ones.

sol_z_complex()

Compute the z solution from the zO solutions at each order. Keep it in terms of the complex coordinates A.

sol_z_polar()

Compute the z solution from the zO solutions at each order and introduce the polar coordinates

system_t0()

Rewrite the system with the fast time scale as the only independent variable.

Notes

This is a trick to use dsolve(), which only accepts functions of 1 variable, to solve higher order equations. The higher order equations are rewritten in terms of temporary coordinates \(\tilde{z}_{i,j}(t_0)\) in place of \(z_{i,j}(\boldsymbol{t})\), with \(i,j\) denoting the oscillator number and \(\epsilon\) order, respectively. This is equivalent to temporary considering that \(\boldsymbol{A}(\boldsymbol{t}_s)\) does not depend on the slow times, which is of no consequence as there are no slow time derivatives appearing in the higher order equations at this stage. Indeed, they were either substituted using the complex modulation equations, or they disappeared when eliminating the secular terms.

steady_state

Started on Tue Feb 15 17:25:59 2022

@author: Vincent MAHE

Analyse systems of coupled nonlinear equations using the Method of Multiple Scales (MMS). This sub-module defines the steady state system from the multiple scales one, and contains the functions for steady state analysis.

class oscilate.MMS.steady_state.Coord_SS

Bases: object

The coordinates used in the steady state analysis.

class oscilate.MMS.steady_state.Forcing_SS(mms)

Bases: object

Define the forcing on the system.

class oscilate.MMS.steady_state.Sol_LC(ss, mms)

Bases: object

Solutions obtained when evaluating the limit cycle (LC) of the system.

class oscilate.MMS.steady_state.Sol_SS(ss, mms)

Bases: object

Solutions obtained when evaluating at steady state.

class oscilate.MMS.steady_state.Sol_bbc(ss, mms)

Bases: object

Solutions obtained when evaluating the backbone curve of the response.

class oscilate.MMS.steady_state.Sol_forced

Bases: object

Solutions obtained when evaluating the forced response of the system.

class oscilate.MMS.steady_state.Stability

Bases: object

Stability analysis parameters and outputs.

class oscilate.MMS.steady_state.Steady_state(mms)

Bases: object

Steady state analysis of the multiple scales system. See Evaluation at steady state for a detailed description of the dynamical system.

Parameters:

mms (Multiple_scales_system) – The multiple scales system.

Jacobian_cartesian()

Compute the Jacobian of the modulation equations systems expressed in cartesian coordinates (see stability_analysis_forced()).

Returns:

J – Jacobian of the cartesian system.

Return type:

sympy.Matrix

Jacobian_polar()

Compute the Jacobian of the modulation equations systems expressed in polar coordinates (see stability_analysis_forced()).

Returns:

J – Jacobian of the polar system.

Return type:

sympy.Matrix

additional_cartesian_substitutions(fp, fq)

Reformulate the already-existing substitutions from polar to cartesian to try and substitute leftover polar terms.

Parameters:

• fp (list of sympy.Expr) – Evolution functions for the cartesian coordinates \(p_i\). There are polar coordinates remaining.

• fq (list of sympy.Expr) – Evolution functions for the cartesian coordinates \(q_i\). There are polar coordinates remaining.

Returns:

• fp (list of sympy.Expr) – Evolution functions for the cartesian coordinates \(p_i\). Additional substitutions were performed to get rid of polar coordinates.

• fq (list of sympy.Expr) – Evolution functions for the cartesian coordinates \(q_i\). Additional substitutions were performed to get rid of polar coordinates.

bifurcation_curves(detJ, sol_type='forced', var_a=False, var_sig=True, solver=<function solve_poly2>)

Compute bifurcation curves, defined by the simple bifurcation points of the slow time system obtained for any forcing frequency and amplitude.

Parameters:

• detJ (sympy.Expr) – The determinant of the matrix.

• sol_type (str, optional) – The solution type, among {“forced”, “bbc”, “LC”} Default is “forced”.

• var_a (bool, optional) – Consider the \(i^{\textrm{th}}\) oscillator’s amplitude \(a_i\) as the variable and find the bifurcation curve as an expression for \(a_i\). detJ is rarely a quadratic polynomial in \(a_i\), so this can rarely be computed easily. Default is False.

• var_sig (bool, optional) – Consider the detuning \(\sigma\) as the variable and find the bifurcation curve as an expression for \(\sigma\). detJ is often a quadratic polynomial in \(\sigma\), so this can often be computed. Default is True.

• solver (function, optional) – The solver to use to compute the bifurcation curves. Available are solver called as solve(expr, x), which solve expr=0 for x. solve() can be used but is sometimes slow. Default is solve_poly2().

Returns:

• bif_a (list) – The bifurcation curves for \(a_i^2\).

• bif_sig (list) – The bifurcation curves for \(\sigma\).

Notes

Simple bifurcations are associated to real eigenvalues of the Jacobian matrix crossing the imaginary axis, hence passing through 0. As the determinant is the product of the eigenvalues, the zeros of the determinant give the bifurcation points of the system.

cartesian_coordinates()

Define cartesian coordinates from the polar ones.

Notes

The homogeneous leading order solution for oscillator \(i\) expressed in polar coordinates takes the form

\[\begin{split}\begin{split} x_{i,0}^{\textrm{h}}(t) & = a_i \cos\left(\frac{r_i}{r_{\textrm{MMS}}}\omega t - \beta_i \right), \\ & = a_i \cos(\beta_i) \cos\left(\frac{r_i}{r_{\textrm{MMS}}}\omega t\right) + a_i \sin(\beta_i) \sin\left(\frac{r_i}{r_{\textrm{MMS}}}\omega t\right) \end{split}\end{split}\]

The polar coordinates are defined as

\[\begin{split}\begin{cases} p_i = a_i \cos (\beta_i), \\ q_i = a_i \sin (\beta_i), \end{cases}\end{split}\]

such that the leading order solution can be written as

\[x_{i,0}^{\textrm{h}}(t) = p_i \cos\left(\frac{r_i}{r_{\textrm{MMS}}}\omega t\right) + q_i \sin\left(\frac{r_i}{r_{\textrm{MMS}}}\omega t\right).\]

check_cartesian_substitutions(a, beta, fp, fq)

Check if substitutions from polar to cartesian coordinates are complete.

Parameters:

• a (list of sympy.Symbol) – Amplitudes of the leading order solutions.

• beta (list of sympy.Symbol) – Phases of the leading order solutions.

• fp (list of sympy.Expr) – Evolution functions for the cartesian coordinates \(p_i\).

• fq (list of sympy.Expr) – Evolution functions for the cartesian coordinates \(q_i\).

Returns:

substitution_OK – True if substitutions are complete. False otherwise.

Return type:

bool

eigenvalues(A, detA=None, trA=None)

Computes the eigenvalues of a matrix \(\textrm{A}\).

Parameters:

• A (sympy.Matrix) – The matrix whose eigenvalues are to be computed.

• detA (sympy.Expr) – Determinant of A. Default is None.

• trA (sympy.Expr) – Trace of A. Default is None.

Returns:

eigvals – The eigenvalues of \(\textrm{A}\).

Return type:

list

modulation_equations_SS(mms)

Evaluate the modulation equations at steady state (polar system).

modulation_equations_cartesian()

Compute the modulation equations of the cartesian coordinates system.

Notes

Write the modulation equations using the cartesian coordinates (defined in cartesian_coordinates()). For oscillator \(i\), this results in

\[\begin{split}\begin{cases} \dfrac{\textrm{d} p_i}{\textrm{d} t} & = f_{p_i}(\boldsymbol{p}, \boldsymbol{q}), \\ \dfrac{\textrm{d} q_i}{\textrm{d} t} & = f_{q_i}(\boldsymbol{p}, \boldsymbol{q}), \end{cases}\end{split}\]

where \(\boldsymbol{p}(t)\) and \(\boldsymbol{q}(t)\) are vectors containing the cartesian coordinates.

polar_coordinates_SS(mms)

Introduce time-independent amplitudes and phases (polar coordinates).

solve_F()

Solve the forced response in terms of the forcing amplitude \(F\). Returns \(F(a_i)\)

solve_LC(solve_dof=None, betai=0)

Find the limit cycle (LC) of a given oscillator with the other oscillators’ amplitude set to 0.

Parameters:

• solve_dof (None or int, oprtional) – The oscillator number to solve for. If None, no oscillator is solved for. Default is None.

• betai (int or Expr, optional) – The initial phase. Default is 0.

solve_LC_amplitude()

Compute the amplitude of the homogeneous, leading order solution on the limit cycle.

solve_LC_frequency()

Compute the oscillation frequency on the limit cycle.

solve_LC_x()

Compute the displacement \(x\) on the limit cycle.

solve_LC_xmax()

Compute the maximum displacement \(x_{\textrm{max}}\) associated to the limit cycle. The choice of the time \(t_{\textrm{max}}\) associated to \(x_{\textrm{max}}\) is based on the fact that the leading order solution for oscillator \(i\), given by

\[x^{\textrm{h}}_{i,0}(t) = a_i \cos\left( \frac{r_i}{r_{\textrm{MMS}}} \omega t - \beta_i\right),\]

dominates, such that

\[t_{\textrm{max}} = \beta_{i, \text{LC}} \frac{r_{\textrm{MMS}}}{r_i \omega},\]

where \(\beta_{i, \text{LC}}\) denotes the phase value associated to the Limit cycle.

solve_a()

Solve the forced response in terms of the oscillator’s amplitude. For readability, the output actually returned is \(a_i^2(\sigma, F)\).

solve_bbc(c=[], solve_dof=None)

Find the backbone curve (bbc) of a given oscillator with the other oscillators’ amplitude set to 0.

Parameters:

• c (list, optional) – Damping terms. They will be set to 0 to compute the backbone curve. Note that these are the scaled damping terms. Default is [].

• solve_dof (None or int, oprtional) – The oscillator number to solve for. If None, no oscillator is solved for. Default is None.

Notes

The backbone curve describes the frequency of free oscillations as a function of the oscillator’s amplitude. In the presence of small damping (as in the case in the MMS) the frequency of free oscillations is close from the resonance frequency. The backbone curve can therefore be interpreted as the backbone of the forced response.

The backbone curve of oscillator \(i\) typically takes the form

\[\omega_{\textrm{bbc}}^{(i)} = k\omega_{i} + f_{\textrm{bbc}}^{(i)} (a_i),\]

where \(k=1,\; k<1,\; k>1\) are associated to direct, superharmonic and subharmonic responses, respectively.

solve_bbc_x()

Compute the displacement \(x\) on the backbone curve.

solve_bbc_xmax()

Compute the maximum displacement \(x_{\textrm{max}}\), associated to the backbone curve. The choice of the time \(t_{\textrm{max}}\) associated to \(x_{\textrm{max}}\) is based on the fact that the leading order solution for oscillator \(i\), given by

\[x^{\textrm{h}}_{i,0}(t) = a_i \cos\left( \frac{r_i}{r_{\textrm{MMS}}} \omega t - \beta_i\right),\]

dominates, such that

\[t_{\textrm{max}} = \beta_{i, \text{bbc}} \frac{r_{\textrm{MMS}}}{r_i \omega},\]

where \(\beta_{i, \text{bbc}}\) denotes the phase value associated to the backbone curve. Note that for a direct response,

\[\beta_{i, \text{bbc}} = \frac{\pi}{2},\]

such that the forcing and displacement are in phase quadrature on the backbone curve. This phase value at resonance is multiplied by a rational in the case of secondary resonance.

solve_forced(solve_dof=None)

Solve the forced response of an oscillator.

Parameters:

solve_dof (None or int, optional) – The oscillator to solve for. If None, no oscillator is solved for. Default is None.

Notes

Find the steady state solution for a given oscillator with the other oscillators’ amplitude set to 0.

To do so, one must choose an oscillator to chose for, say oscillator \(i\). Then, the following methods are called:

1. substitution_solve_dof(): Set the other oscillators’ amplitude to 0, i.e. \(a_j = 0 \; \forall j \neq i\).

2. solve_phase(): express the oscillator’s phase \(\beta_i\) as a function of its amplitude \(a_i\).

3. solve_sigma(): find the expression of \(\sigma(a_i)\).

4. solve_a(): find the expression of \(a_i (\sigma, F)\).

5. solve_F(): find the expression of \(F(a_i)\).

solve_phase()

Find solutions for the oscillator’s phase \(\beta_i\). The solutions actually returned are \(\sin(\beta_i)\) and \(\cos(\beta_i)\).

solve_sigma()

Solve the forced response in terms of the detuning \(\sigma\). Returns \(\sigma(a_i)\). It is recalled that \(\omega = \omega_{\textrm{MMS}} + \epsilon \sigma\).

stability_analysis_forced(coord='cartesian', rewrite_polar=True, eigenvalues=False, bifurcation_curves=False, trace_curves=False, analyse_blocks=False, kwargs_bif={})

Evaluate the stability of a steady state, forced solution. See Stability analysis for a detailed description of the dynamical system.

Parameters:

• coord (str, optional) – Either "cartesian" or "polar". Specifies the coordinates to use for the stability analysis. "cartesian" is recommended as it prevents divisions by 0, which occur when at least one of the oscillator has a 0 ampliutude. Default is "cartesian".

• rewrite_polar (bool, optional) – Rewrite the Jacobian’s determinant and trace in polar coordinates (if computed using cartesian ones). This can be time consuming and the current back substitutions from cartesian to polar coordinates are not always sufficient. Default is True.

• eigenvalues (bool, optional) – Compute the eigenvalues of the Jacobian. Default is False.

• bifurcation_curves (bool, optional) – Compute the bifurcation curves. Default is False.

• trace_curves (bool, optional) – Compute the trace curves, i.e. the zeros of the Jacobian’s trace. Default is False.

• analyse_blocks (bool, optional) – Analyse the diagonal blocks of the Jacobian rather than the Jacobian itself. This is relevant if the Jacobian is block-diagonal.

• kwargs_bif (dict, optional) – Passed to bifurcation_curves() Default is dict().

substitution_solve_dof(solve_dof)

Set every oscillator amplitude to 0 except the one to solve for.

Notes

If one wants to solve for \(a_i\), then the system is evaluated for \(a_j=0, \; \forall j \neq i\).

trace_curves(trJ, sol_type='forced', var_a=False, var_sig=True, solver=<function solve_poly2>)

Compute curves capturing the variations of the trace of the Jacobian. For 2 by 2 Jacobians, a negative trace with a positive determinant indicates a stable response.

Parameters:

• trJ (sympy.Expr) – The trace of the matrix.

• sol_type (str, optional) – The solution type, among {“forced”, “bbc”, “LC”} Default is “forced”.

• var_a (bool, optional) – Consider the \(i^{\textrm{th}}\) oscillator’s amplitude \(a_i\) as the variable and find the trace curve as an expression for \(a_i\). Default is False.

• var_sig (bool, optional) – Consider the detuning \(\sigma\) as the variable and find the trace curve as an expression for \(\sigma\). Default is True.

• solver (function, optional) – The solver to use to compute the trace curves. Available are solver called as solve(expr, x), which solve expr=0 for x. solve() can be used but is sometimes slow. Default is solve_poly2().

Returns:

• tr_a (list) – The trace curves for \(a_i^2\).

• tr_sig (list) – The trace curves for \(\sigma\).

class oscilate.MMS.steady_state.Substitutions_SS(mms)

Bases: object

Substitutions used in the steady state evaluations.

visualisation

Started on Tue Feb 15 17:25:59 2022

@author: Vincent MAHE

Analyse systems of coupled nonlinear equations using the Method of Multiple Scales (MMS). This sub-module evaluates the symbolic expressions for given numerical parameters and allows to plot the resulting numerical results.

class oscilate.MMS.visualisation.Amplitude_response_curve(mms, ss, dyn, param)

Bases: object

Evaluate the amplitude response curves (ARC) for given numerical parameters. ‘Amplitude’ refers to the forcing amplitude. This transforms the sympy expressions to numpy arrays. They can then be plotted.

Parameters:

• mms (Multiple_scales_system) – The MMS object.

• ss (Steady_state) – The MMS results evaluated at steady state.

• dyn (Dynamical_system) – The initial dynamical system.

• param (list[tuple]) –

  A list whose values are tuples with 2 elements:

  1. The sympy symbol of a parameter,

  2. The numerical value(s) taken by that parameter.

plot(**kwargs)

Plots the amplitude-response curves (ARC), both forcing amplitude-amplitude and forcing amplitude-phase.

Parameters:

ss (Steady_state, optional) – Steady state object. Used to name the axis labels.

Returns:

• fig1 (Figure) – The amplitude plot \(a(F)\).

• fig2 (Figure) – The phase plot \(\beta(F)\).

class oscilate.MMS.visualisation.Backbone_curve(mms, ss, dyn, param)

Bases: object

Evaluate the backbone curve for given numerical parameters. This transforms the sympy expressions to numpy arrays. They can then be plotted.

Parameters:

• mms (Multiple_scales_system) – The MMS object.

• ss (Steady_state) – The MMS results evaluated at steady state.

• dyn (Dynamical_system) – The initial dynamical system.

• param (list[tuple]) –

  A list whose values are tuples with 2 elements:

  1. The sympy symbol of a parameter,

  2. The numerical value(s) taken by that parameter.

plot(**kwargs)

Plots the backbone curve.

Parameters:

ss (Steady_state, optional) – Steady state object. Used to name the axis labels.

Returns:

fig – The amplitude plot \(a(\omega)\).

Return type:

Figure

class oscilate.MMS.visualisation.Frequency_response_curve(mms, ss, dyn, param, bif=True)

Bases: object

Evaluate the frequency response curves (FRC) and bifurcation curves (if computed) for given numerical parameters. This transforms the sympy expressions to numpy arrays. They can then be plotted.

Parameters:

• mms (Multiple_scales_system) – The MMS object.

• ss (Steady_state) – The MMS results evaluated at steady state.

• dyn (Dynamical_system) – The initial dynamical system.

• param (list[tuple]) –

  A list whose values are tuples with 2 elements:

  1. The sympy symbol of a parameter,

  2. The numerical value(s) taken by that parameter.

• bif (bool, optional) – Evaluate the bifurcation curves. Default is True.

plot(bbc=None, **kwargs)

Plots the frequency response curves (FRC), both frequency-amplitude and frequency-phase. Also includes the stability information if computed, and the backbone curve if given.

Parameters:

• bbc (Backbone_curve, optional) – The evaluated backbone curve. Default is None.

• ss (Steady_state, optional) – Steady state object. Used to name the axis labels.

Returns:

• fig1 (Figure) – The amplitude plot \(a(\omega)\).

• fig2 (Figure) – The phase plot \(\beta(\omega)\).

class oscilate.MMS.visualisation.Limit_cycle(mms, ss, param, Npts=1000)

Bases: object

Evaluate the limit cycle for given numerical parameters. This transforms the sympy expressions to numpy arrays. They can then be plotted in the phase portrait and as time signals.

Parameters:

• mms (Multiple_scales_system) – The MMS object.

• ss (Steady_state) – The SS object.

• param (list[tuple]) –

  A list whose values are tuples with 2 elements:

  1. The sympy symbol of a parameter,

  2. The numerical value(s) taken by that parameter.

• Npts (int, optional) – Number of time points. Default is 1000.

plot_PP(c='tab:blue', lw=2, **kwargs)

Plots the transient trajectory in the phase portrait.

Parameters:

• c (str, optional) – The color of the line.

• lw (float, optional) – The linewidth Default is 2.

Returns:

fig – The phase portrait plot.

Return type:

Figure

plot_time(c='tab:blue', **kwargs)

Plots the transient time response.

Parameters:

c (str, optional) – The color of the line.

Returns:

fig – The time plot.

Return type:

Figure

class oscilate.MMS.visualisation.Transient_response(mms, param)

Bases: object

Evaluate the transient response for given numerical parameters. This transforms the sympy expressions to numpy arrays. They can then be plotted in the phase portrait and as time signals.

Parameters:

• mms (Multiple_scales_system) – The MMS object.

• param (list[tuple]) –

  A list whose values are tuples with 2 elements:

  1. The sympy symbol of a parameter,

  2. The numerical value(s) taken by that parameter.

plot_PP(c='tab:blue', **kwargs)

Plots the transient trajectory in the phase portrait.

Parameters:

c (str, optional) – The color of the line.

Returns:

fig – The phase portrait plot.

Return type:

Figure

plot_time(c='tab:blue', **kwargs)

Plots the transient time response.

Parameters:

c (str, optional) – The color of the line.

Returns:

fig – The time plot.

Return type:

Figure

oscilate.MMS.visualisation.numpise_F_ARC(mms, ss, param)

Numpise the forced response’s forcing amplitude \(F\).

Parameters:

• mms (Multiple_scales_system)

• ss (Steady_state)

• param (dict) – See sympy_to_numpy().

Returns:

F – Numpised forced response’s forcing amplitude.

Return type:

numpy.ndarray

oscilate.MMS.visualisation.numpise_LC(mms, ss, param)

Numpise the limit cycle solution.

Parameters:

• mms (Multiple_scales_system)

• ss (Steady_state_system)

• param (dict) – See sympy_to_numpy().

Returns:

• a (float) – Numpised LC amplitude.

• beta (float) – Numpised LC initial phase.

• omega (float) – Numpised LC frequency.

oscilate.MMS.visualisation.numpise_LC_trajectory(mms, ss, param)

Numpise the oscillator’s LC trajectory

Parameters:

• mms (Multiple_scales_system)

• ss (Steady_state_system)

• param (dict) – See sympy_to_numpy().

Returns:

• x (numpy.ndarray) – Numpised LC motion.

• dxdt (numpy.ndarray) – Numpised LC velocity.

oscilate.MMS.visualisation.numpise_omegaMMS(mms, param)

Numpise the frequency around which a solution is sought.

Parameters:

• mms (Multiple_scales_system)

• param (dict) – See sympy_to_numpy().

Returns:

omegaMMS – Numpised MMS frequency.

Return type:

float

oscilate.MMS.visualisation.numpise_omega_FRC(mms, ss, param)

Numpise the forced response’s frequency \(\omega\).

Parameters:

• mms (Multiple_scales_system)

• ss (Steady_state)

• param (dict) – See sympy_to_numpy().

Returns:

omega – Numpised forced response’s frequency.

Return type:

numpy.ndarray

oscilate.MMS.visualisation.numpise_omega_bbc(mms, ss, param)

Numpise the backbone curve’s frequency \(\omega_{\textrm{bbc}}\).

Parameters:

• mms (Multiple_scales_system)

• ss (Steady_state)

• param (dict) – See sympy_to_numpy().

Returns:

omega_bbc – Numpised backbone curve’s frequency.

Return type:

numpy.ndarray

oscilate.MMS.visualisation.numpise_omega_bif(mms, ss, param)

Numpise the bifurcation curves’ frequency \(\omega_{\textrm{bif}}\).

Parameters:

• mms (Multiple_scales_system)

• ss (Steady_state)

• param (dict) – See sympy_to_numpy().

Returns:

omega_bif – Numpised bifurcation curves’ frequency.

Return type:

list of numpy.ndarray

oscilate.MMS.visualisation.numpise_phase(mms, ss, dyn, param, omega, F)

Numpise the phase \(\beta_i\).

Parameters:

• mms (Multiple_scales_system)

• ss (Steady_state)

• dyn (Dynamical_system)

• param (dict) – See sympy_to_numpy().

• omega (numpy.ndarray) – The frequency array.

• F (numpy.ndarray) – The forcing amplitude array.

Returns:

phase – Numpised phase.

Return type:

numpy.ndarray

oscilate.MMS.visualisation.numpise_transient_slow_time(mms, param)

Numpise the slow time transient response.

Parameters:

• mms (Multiple_scales_system)

• param (dict) – See sympy_to_numpy().

Returns:

• a (numpy.ndarray) – Numpised transient amplitude.

• psi (numpy.ndarray) – Numpised transient absolute phase.

• dadt (numpy.ndarray) – Numpised time derivative of the transient amplitude.

• dpsidt (numpy.ndarray) – Numpised time derivative of the transient absolute phase.

oscilate.MMS.visualisation.numpise_transient_trajectory(mms, param)

Numpise the transient oscillator’s trajectory

Parameters:

• mms (Multiple_scales_system)

• param (dict) – See sympy_to_numpy().

Returns:

• x (numpy.ndarray) – Numpised transient motion.

• dxdt (numpy.ndarray) – Numpised transient velocity.

oscilate.MMS.visualisation.numpise_xmax_bbc(mms, ss, param)

Numpise the peak oscillator’s amplitude \(x_{\textrm{max}}\) on the backbone curve.

Parameters:

• mms (Multiple_scales_system)

• ss (Steady_state)

• param (dict) – See sympy_to_numpy().

Returns:

xmax – Numpised peak amplitude on the backbone curve.

Return type:

numpy.ndarray

sympy_functions Module

Started on Wed Apr 9 13:39:41 2025

@author: Vincent MAHE

Sympy functions useful to the use of MMS functions.

oscilate.sympy_functions.check_solvability(poly, x)

Check the solvability of a polynomial \(p(x)\).

Parameters:

• poly (sympy.Expr) – The polynomial considered.

• x (sympy.Symbol) – The variable to solve for.

Returns:

bool – True is solvable, False otherwise.

Return type:

bool,

oscilate.sympy_functions.get_block_diagonal_indices(matrix, block_sizes)

Generate a list of \((i, j)\) indices for all elements in the diagonal blocks of a block-diagonal matrix.

Parameters:

• matrix (sympy.Matrix) – The matrix to check for block diagonality.

• block_sizes (int or list of int) – Size(s) of the diagonal blocks.

Returns:

indices – A list of tuples (i, j) representing the indices of elements in the diagonal blocks.

Return type:

list

oscilate.sympy_functions.get_exponent(expr, x)

Get the exponent of \(x\) in an expression of the type \(\lambda x^n\) where \(\lambda\) is a constant while \(n\) is an integer or rational.

Parameters:

• expr (sympy.Expr) – The expression in which one wants to identify the exponent of x.

• x (sympy.Symbol) – The variable whose exponent is to be known.

oscilate.sympy_functions.is_block_diagonal(matrix, block_sizes)

Check if a matrix is block-diagonal given block sizes.

Parameters:

• matrix (sympy.Matrix) – The matrix to check for block diagonality.

• block_sizes (int or list of int) – Size(s) of the diagonal blocks.

Returns:

bool

Return type:

True if the matrix is block-diagonal, False otherwise.

oscilate.sympy_functions.polynomial_terms(poly, x)

Identify the terms of a polynomial.

Parameters:

• poly (sympy.Expr) – The polynomial considered.

• x (sympy.Symbol) – The variable to solve for.

Returns:

dic_x – The polynomial terms.

Return type:

dict

Notes

If the expression given for poly is of the form

\[p(x) = q(x) x^{-n},\]

where the powers of \(x\) in \(q(x)\) are all superior or equal to \(0\), then an auxiliary polynomial

\[P(x) = \dfrac{p(x)}{x^{-n}}\]

is constructed. It is the terms of that positive powers polynomial that are returned.

oscilate.sympy_functions.solve_poly2(poly, x)

Finds the roots of a polynomial of degree 2.

Parameters:

• poly (sympy.Expr) – polynomial whose roots are to be computed

• x (sympy.Symbol) – Variable of the polynomial

Returns:

x_sol – list containing the two roots of the polynomial

Return type:

list of sympy.Expr

Notes

The polynomial of degree 2 takes the form

\[p(x) = a x^2 + bx + c.\]

Note that \(b\) can be null but not \(a\) nor \(c\).

It is a workaround to using solve() or solveset(). These two work but can be very long when coefficients \(a,\; b,\; c\) are expressions involving many parameters. Note that solve() is significantly slower than solveset().

oscilate.sympy_functions.sub_deep(expr, sub)

Performs deep substitutions of an expression.

Parameters:

• expr (sympy.Expr) – Expression on which substitutions are to be performed.

• sub (list of tuples) – The substitutions to perform.

Returns:

expr_sub – The expression with substitutions performed.

Return type:

sympy.Expr

Notes

Deep substitutions are needed when a substitution involves terms that can still be substituted. For instance, one wants to substitute \(a_1\) and \(a_2\) by expressions, but \(a_1\) is actually a function of \(a_2\), so at least 2 substitutions are required.

oscilate.sympy_functions.sympy_to_numpy(expr_sy, param)

Transform a sympy expression into a numpy array.

Parameters:

• expr_sy (sympy.Expr) – A sympy expression.

• param (dict) –

  A dictionnary whose values are tuples with 2 elements:

  1. The sympy symbol of a parameter

  2. The numerical value(s) taken by that parameter

Returns:

expr_np – The numerical values taken by the sympy expression evaluated.

Return type:

numpy.ndarray