AsymptoticStability
- class sharpy.postproc.asymptoticstability.AsymptoticStability[source]
Calculates the asymptotic stability properties of the linearised system by computing the corresponding eigenvalues and eigenvectors.
The output of this solver is written to the
stability
directory in the case output.The stability of the systems specified in
target_systems
is performed. If the system has been previously scaled, areference_velocity
should be provided to compute the stability at such point.The eigenvalues can be truncated, keeping a minimum
num_evals
(sorted by decreasing real part) or by limiting the higher frequency modes throughfrequency_cutoff
.Results can be saved to file using
export_eigenvalues
. The settingdisplay_root_locus
shows a simple Argand diagram where the continuous time eigenvalues are displayed.The eigenvectors can be displayed (in .vtu format for use in Paraview) with the
modes_to_plot
setting, whereby the user specifies the mode indices that are to be plotted (sorted with in the same way as the eigenvalue table, i.e. by decreasing real part). A snapshot of the eigenvector is produced for the 0, 45, 90 and 135 degree phases. This feature currently supports the flexible structural modes and does not show the rigid body contribution. The output is written to thestability/modes
folder and includes the structure at the reference linearisation state.The settings that this solver accepts are given by a dictionary, with the following key-value pairs:
Name
Type
Description
Default
print_info
bool
Print information and table of eigenvalues
False
reference_velocity
float
Reference velocity at which to compute eigenvalues for scaled systems
1.0
frequency_cutoff
float
Truncate higher frequency modes. If zero none are truncated
0
export_eigenvalues
bool
Save eigenvalues and eigenvectors to file.
False
output_file_format
str
Eigenvalue/eigenvector output file format. HDF5 or text (.dat) files.
dat
display_root_locus
bool
Show plot with eigenvalues on Argand diagram
False
velocity_analysis
list(float)
List containing min, max and number of velocities to analyse the system
[]
target_system
list(str)
System or systems for which to find frequency response.
['aeroelastic']
num_evals
int
Number of eigenvalues to retain.
200
modes_to_plot
list(int)
List of mode numbers to plot. Plots the 0, 45, 90 and 135degree phases.
[]
- compute_eigenvalues(ss, system_name_list=None, not_scaled=True)[source]
Computes the eigenvalues and eigenvectors of the state-space
- Parameters:
ss (libss.StateSpace or list([libss.StateSpace]) – State-space or list of state-spaces
system_name_list (list([str]) – Names of systems in the case multiple systems are required
not_scaled (bool) – Flag to indicate whether the systems are assembled in non-dimensional time
- convert_to_continuoustime(dt, discrete_time_eigenvalues, not_scaled=False)[source]
Convert eigenvalues to discrete time. The
not_scaled
argument can be used to bypass the search from within SHARPy of scaling factors. For instance, when the state-space of choice is not part of a standard SHARPy case but rather an interpolated ROM etc.The eigenvalues are converted to continuous time using
\[\lambda_{ct} = \frac{\log (\lambda_{dt})}{\Delta t}\]If the system is scaled, the dimensional time step is retrieved as
\[\Delta t_{dim} = \bar{\Delta t} \frac{l_{ref}}{U_{\infty, actual}}\]where \(l_{ref}\) is the reference length and \(U_{\infty, actual}\) is the free stream velocity at which to calculate the eigenvalues.
- Parameters:
dt (float) – Discrete time increment.
discrete_time_eigenvalues (np.ndarray) – Array of discrete time eigenvalues.
not_scaled (bool) – Treat the system as not scaled. No Scaling Factors will be searched in SHARPy.
- static display_root_locus(eigenvalues)[source]
Displays root locus diagrams.
Returns the
fig
andax
handles for further editing.- Returns:
ax:
- Return type:
fig
- export_eigenvalues(num_evals, eigenvalues, eigenvectors, filename=None)[source]
Saves a
num_evals
number of eigenvalues and eigenvectors to file.The files are saved in the output directory and include:
{system_name}_eigenvalues.dat
: Array of eigenvalues of shape(num_evals, 2)
where the first column corresponds to the real part and the second column to the imaginary part.{system_name}_stability.h5
: An.h5
file containing the desired number of eigenvalues and eigenvectors of the chosen systems.
The units of the eigenvalues are
rad/s
.- Parameters:
num_evals (int) – Number of eigenvalues to save.
eigenvalues (np.ndarray) – Eigenvalue array.
eigenvectors (np.ndarray) – Matrix of eigenvectors.
filename (str (optional)) – Optional prefix of the output filenames.
See also
Loading and saving complex arrays: https://stackoverflow.com/questions/6494102/how-to-save-and-load-an-array-of-complex-numbers-using-numpy-savetxt/6522396
- plot_modes(eigenvectors)[source]
Plot the aeroelastic mode shapes for the first
n_modes_to_plot
Plots the 0, 45, 90 and 135 degrees phase of the mode. Also plots the reference at the linearisation state.
\[x_{out} = Re(\Phi_i e^{i\theta})\]for \(\theta \in \{0, \pi/4, \pi/2, 3\pi/4}\).
- run(**kwargs)[source]
Computes the eigenvalues and eigenvectors
- Returns:
Eigenvalues sorted and frequency truncated eigenvectors (np.ndarray): Corresponding mode shapes
- Return type:
eigenvalues (np.ndarray)
- static sort_eigenvalues(eigenvalues, eigenvectors, frequency_cutoff=0, number_of_eigenvalues=None)[source]
Sort continuous-time eigenvalues by order of magnitude.
The conjugate of complex eigenvalues is removed, then if specified, high frequency modes are truncated. Finally, the eigenvalues are sorted by largest to smallest real part.
- Parameters:
eigenvalues (np.ndarray) – Continuous-time eigenvalues
eigenvectors (np.ndarray) – Corresponding right eigenvectors
frequency_cutoff (float) – Cutoff frequency for truncation
[rad/s]
number_of_eigenvalues (int (optional)) – Number of eigenvalues to retain
- Returns:
eigenvalues and eigenvectors
- Return type:
tuple(np.array, np.array)
- velocity_analysis()[source]
Velocity analysis for scaled systems.
Runs the stability analysis for different velocities for aeroelastic systems that have been previously scaled.
For every velocity, the linear system is updated. This involves updating the structural matrices and the coupling matrix. The eigenvalues saved are in continuous time.
It saves the results to a
.dat
file where the first column corresponds to the free stream velocity and the second and third columns to the real and imaginary parts of the eigenvalues.