class sharpy.linear.src.linuvlm.DynamicBlock(tsdata, dt=None, dynamic_settings=None, integr_order=2, RemovePredictor=True, ScalingDict=None, UseSparse=True, for_vel=numpy.zeros)[source]

Class for dynamic linearised UVLM solution. Linearisation around steady-state are only supported.

The class is a low-memory implementation of Dynamic, and inherits most of the methods contained there. State-space models are allocated in list-block form (as per numpy.block) to minimise memory usage. This class provides lower memory / computational time assembly, frequency response and frequency limited balancing.


  • tsdata: aero timestep data from SHARPy solution

  • dt: time-step

  • integr_order=2: integration order for UVLM unsteady aerodynamic force

  • RemovePredictor=True: if true, the state-space model is modified so as to accept in input perturbations, u, evaluated at time-step n rather than n+1.

  • ScalingDict=None: disctionary containing fundamental reference units

    >>> {'length':  reference_length,
         'speed':   reference_speed,
         'density': reference density}

    used to derive scaling quantities for the state-space model variables. The scaling factors are stores in self.ScalingFact.

    Note that while time, circulation, angular speeds) are scaled accordingly, FORCES ARE NOT. These scale by qinf*b**2, where b is the reference length and qinf is the dinamic pressure.

  • UseSparse=False: builds the A and B matrices in sparse form. C and D are dense, hence the sparce format is not used.

- nondimss

normalises a dimensional state-space model based on the scaling factors in self.ScalingFact.

- dimss

inverse of nondimss.

- assemble_ss

builds state-space model. See function for more details.

- assemble_ss_profiling

generate profiling report of the assembly and saves it into self.prof_out. To read the report:

>>> import pstats
- freqresp

ad-hoc method for fast frequency response (only implemented) for remove_predictor=False

To do: upgrade to linearise around unsteady snapshot (adjoint)


Produces block-form of state-space model

\[\begin{split}\mathbf{x}_{n+1} &= \mathbf{A}\,\mathbf{x}_n + \mathbf{B} \mathbf{u}_{n+1} \\ \mathbf{y}_n &= \mathbf{C}\,\mathbf{x}_n + \mathbf{D} \mathbf{u}_n\end{split}\]

where the state, inputs and outputs are:

\[\mathbf{x}_n = \{ \delta \mathbf{\Gamma}_n,\, \delta \mathbf{\Gamma_{w_n}},\, \Delta t\,\delta\mathbf{\Gamma}'_n,\, \delta\mathbf{\Gamma}_{n-1} \}\]
\[\mathbf{u}_n = \{ \delta\mathbf{\zeta}_n,\, \delta\mathbf{\zeta}'_n,\, \delta\mathbf{u}_{ext,n} \}\]
\[\mathbf{y} = \{\delta\mathbf{f}\}\]

with \(\mathbf{\Gamma}\) being the vector of vortex circulations, \(\mathbf{\zeta}\) the vector of vortex lattice coordinates and \(\mathbf{f}\) the vector of aerodynamic forces and moments. Note that \((\bullet)'\) denotes a derivative with respect to time.

Note that the input is atypically defined at time n+1, therefore by default self.remove_predictor = True and the predictor term u_{n+1} is eliminated through the change of state[1]:

\[\begin{split}\mathbf{h}_n &= \mathbf{x}_n - \mathbf{B}\,\mathbf{u}_n \\\end{split}\]

such that:

\[\begin{split}\mathbf{h}_{n+1} &= \mathbf{A}\,\mathbf{h}_n + \mathbf{A\,B}\,\mathbf{u}_n \\ \mathbf{y}_n &= \mathbf{C\,h}_n + (\mathbf{C\,B}+\mathbf{D})\,\mathbf{u}_n\end{split}\]

which only modifies the equivalent \(\mathbf{B}\) and \(\mathbf{D}\) matrices.


[1] Franklin, GF and Powell, JD. Digital Control of Dynamic Systems, Addison-Wesley Publishing Company, 1980

To do: - remove all calls to scipy.linalg.block_diag

balfreq(DictBalFreq, wake_prop_settings=None)[source]

Low-rank method for frequency limited balancing. The Observability ad controllability Gramians over the frequencies kv are solved in factorised form. Balancd modes are then obtained with a square-root method.

Details: Observability and controllability Gramians are solved in factorised form through explicit integration. The number of integration points determines both the accuracy and the maximum size of the balanced model.

Stability over all (Nb) balanced states is achieved if:
  1. one of the Gramian is integrated through the full Nyquist range

  2. the integration points are enough.

Note, however, that even when stability is not achieved over the full balanced states, stability of the balanced truncated model with Ns<=Nb states is normally observed even when a low number of integration points is used. Two integration methods (trapezoidal rule on uniform grid and Gauss-Legendre quadrature) are provided.


  • DictBalFreq: dictionary specifying integration method with keys:

    • ‘frequency’: defines limit frequencies for balancing. The balanced

    model will be accurate in the range [0,F], where F is the value of this key. Note that F units must be consistent with the units specified in the self.ScalingFacts dictionary.

    • ‘method_low’: [‘gauss’,’trapz’] specifies whether to use gauss

    quadrature or trapezoidal rule in the low-frequency range [0,F]

    • ‘options_low’: options to use for integration in the low-frequencies.

    These depend on the integration scheme (See below).

    • ‘method_high’: method to use for integration in the range [F,F_N],

    where F_N is the Nyquist frequency. See ‘method_low’.

    • ‘options_high’: options to use for integration in the high-frequencies.

    • ‘check_stability’: if True, the balanced model is truncated to

    eliminate unstable modes - if any is found. Note that very accurate balanced model can still be obtained, even if high order modes are unstable. Note that this option is overridden if “”

    • ‘get_frequency_response’: if True, the function also returns the

    frequency response evaluated at the low-frequency range integration points. If True, this option also allows to automatically tune the balanced model.

Future options:

  • ‘truncation_tolerance’: if ‘get_frequency_response’ is True, allows

to truncatethe balanced model so as to achieved a prescribed tolerance in the low-frequwncy range.

  • Ncpu: for parallel run

The following integration schemes are available:
  • ‘trapz’: performs integration over equally spaced points using

trapezoidal rule. It accepts options dictionaries with keys:
  • ‘points’: number of integration points to use (including

domain boundary)

  • ‘gauss’ performs gauss-lobotto quadrature. The domain can be

partitioned in Npart sub-domain in which the gauss-lobotto quadrature of order Ord can be applied. A total number of Npart*Ord points is required. It accepts options dictionaries of the form:

  • ‘partitions’: number of partitions

  • ‘order’: quadrature order.

Example: The following dictionary

DictBalFreq={ ‘frequency’: 1.2,

‘method_low’: ‘trapz’, ‘options_low’: {‘points’: 12}, ‘method_high’: ‘gauss’, ‘options_high’: {‘partitions’: 2, ‘order’: 8}, ‘check_stability’: True }

balances the state-space model self.SS in the frequency range [0, 1.2] using

  1. 12 equally-spaced points integration of the Gramians in

the low-frequency range [0,1.2] and

(b) a 2 Gauss-Lobotto 8-th order quadratures of the controllability Gramian in the high-frequency range.

A total number of 28 integration points will be required, which will result into a balanced model with number of states

min{ 2*28* number_inputs, 2*28* number_outputs }

The model is finally truncated so as to retain only the first Ns stable modes.

freqresp(kv, wake_prop_settings=None)[source]

Ad-hoc method for fast UVLM frequency response over the frequencies kv. The method, only requires inversion of a K x K matrix at each frequency as the equation for propagation of wake circulation are solved exactly. The algorithm implemented here can be used also upon projection of the state-space model.

Note: This method is very similar to the “minsize” solution option is the steady_solve.


Scale state-space model based of self.ScalingFacts.

solve_step(x_n, u_n, u_n1=None, transform_state=False)[source]

Solve step.

If the predictor term has not been removed (remove_predictor = False) then the system is solved as:

\[\begin{split}\mathbf{x}^{n+1} &= \mathbf{A\,x}^n + \mathbf{B\,u}^n \\ \mathbf{y}^{n+1} &= \mathbf{C\,x}^{n+1} + \mathbf{D\,u}^n\end{split}\]

Else, if remove_predictor = True, the state is modified as

\[\mathbf{h}^n = \mathbf{x}^n - \mathbf{B\,u}^n\]

And the system solved by:

\[\begin{split}\mathbf{h}^{n+1} &= \mathbf{A\,h}^n + \mathbf{B_{mod}\,u}^{n} \\ \mathbf{y}^{n+1} &= \mathbf{C\,h}^{n+1} + \mathbf{D_{mod}\,u}^{n+1}\end{split}\]

Finally, the original state is recovered using the reverse transformation:

\[\mathbf{x}^{n+1} = \mathbf{h}^{n+1} + \mathbf{B\,u}^{n+1}\]

where the modifications to the \(\mathbf{B}_{mod}\) and \(\mathbf{D}_{mod}\) are detailed in Dynamic.assemble_ss().


Although the original equations include the term \(\mathbf{u}_{n+1}\), it is a reasonable approximation to take \(\mathbf{u}_{n+1}\approx\mathbf{u}_n\) given a sufficiently small time step, hence if the input at time n+1 is not parsed, it is estimated from \(u^n\).

  • x_n (np.array) – State vector at the current time step \(\mathbf{x}^n\)

  • u_n (np.array) – Input vector at time step \(\mathbf{u}^n\)

  • u_n1 (np.array) – Input vector at time step \(\mathbf{u}^{n+1}\)

  • transform_state (bool) – When the predictor term is removed, if true it will transform the state vector. If false it will be assumed that the state vector that is parsed is already transformed i.e. it is \(\mathbf{h}\).


Updated state and output vector packed in a tuple \((\mathbf{x}^{n+1},\,\mathbf{y}^{n+1})\)

Return type



Because in BlockDynamics the predictor is never removed when building ‘self.SS’, the implementation change with respect to Dynamic. However, formulas are consistent.