{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# Asymptotic Stability of a Flying Wing in Cruise Trimmed Conditions\n", "\n", "A Horten flying wing is analysed. The nonlinear trim condition is found and the system is linearised. The eigenvalues of the linearised system are then used to evaluate the stability at the cruise trimmed flight conditions." ] }, { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [], "source": [ "# required packages\n", "import sharpy.utils.algebra as algebra\n", "import sharpy.sharpy_main\n", "from sharpy.cases.hangar.richards_wing import Baseline\n", "import numpy as np\n", "import configobj\n", "import matplotlib.pyplot as plt" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Flight Conditions\n", "\n", "Initial flight conditions. The values for angle of attack ``alpha``, control surface deflection ``cs_deflection`` and ``thrust`` are only initial values. The values required for trim will be calculated by the ``StaticTrim`` routine" ] }, { "cell_type": "code", "execution_count": 2, "metadata": {}, "outputs": [], "source": [ "u_inf = 28\n", "alpha_deg = 4.5135\n", "cs_deflection = 0.1814\n", "thrust = 5.5129" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Discretisation" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "M = 4 # chordwise panels\n", "N = 11 # spanwise panels\n", "Msf = 5 # wake length in chord numbers" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Create Horten Wing" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "0\n", "Section Mass: 11.88 \n", "Linear Mass: 11.88\n", "Section Ixx: 1.8777\n", "Section Iyy: 1.0137\n", "Section Izz: 2.5496\n", "Linear Ixx: 1.88\n", "1\n", "Section Mass: 10.99 \n", "Linear Mass: 10.99\n", "Section Ixx: 1.4694\n", "Section Iyy: 0.9345\n", "Section Izz: 2.1501\n", "Linear Ixx: 1.74\n", "2\n", "Section Mass: 10.10 \n", "Linear Mass: 10.10\n", "Section Ixx: 1.1257\n", "Section Iyy: 0.8561\n", "Section Izz: 1.7993\n", "Linear Ixx: 1.60\n", "3\n", "Section Mass: 9.21 \n", "Linear Mass: 9.21\n", "Section Ixx: 0.8410\n", "Section Iyy: 0.7783\n", "Section Izz: 1.4933\n", "Linear Ixx: 1.46\n", "4\n", "Section Mass: 8.32 \n", "Linear Mass: 8.32\n", "Section Ixx: 0.6096\n", "Section Iyy: 0.7011\n", "Section Izz: 1.2280\n", "Linear Ixx: 1.31\n", "5\n", "Section Mass: 7.43 \n", "Linear Mass: 7.43\n", "Section Ixx: 0.4260\n", "Section Iyy: 0.6246\n", "Section Izz: 0.9996\n", "Linear Ixx: 1.17\n", "6\n", "Section Mass: 6.54 \n", "Linear Mass: 6.54\n", "Section Ixx: 0.2845\n", "Section Iyy: 0.5485\n", "Section Izz: 0.8040\n", "Linear Ixx: 1.03\n", "7\n", "Section Mass: 5.64 \n", "Linear Mass: 5.64\n", "Section Ixx: 0.1796\n", "Section Iyy: 0.4728\n", "Section Izz: 0.6374\n", "Linear Ixx: 0.89\n", "8\n", "Section Mass: 4.75 \n", "Linear Mass: 4.75\n", "Section Ixx: 0.1055\n", "Section Iyy: 0.3975\n", "Section Izz: 0.4959\n", "Linear Ixx: 0.75\n", "9\n", "Section Mass: 3.86 \n", "Linear Mass: 3.86\n", "Section Ixx: 0.0567\n", "Section Iyy: 0.3226\n", "Section Izz: 0.3753\n", "Linear Ixx: 0.61\n", "10\n", "Section Mass: 2.97 \n", "Linear Mass: 2.97\n", "Section Ixx: 0.0275\n", "Section Iyy: 0.2479\n", "Section Izz: 0.2719\n", "Linear Ixx: 0.47\n" ] } ], "source": [ "ws = Baseline(M=M,\n", " N=N,\n", " Mstarfactor=Msf,\n", " u_inf=u_inf,\n", " rho=1.02,\n", " alpha_deg=alpha_deg,\n", " roll_deg=0,\n", " cs_deflection_deg=cs_deflection,\n", " thrust=thrust,\n", " physical_time=20,\n", " case_name='horten',\n", " case_name_format=4,\n", " case_remarks='M%gN%gMsf%g' % (M, N, Msf))\n", "\n", "ws.set_properties()\n", "ws.initialise()\n", "ws.clean_test_files()\n", "\n", "ws.update_mass_stiffness(sigma=1., sigma_mass=2.5)\n", "ws.update_fem_prop()\n", "ws.generate_fem_file()\n", "ws.update_aero_properties()\n", "ws.generate_aero_file()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Simulation Information\n", "\n", "The ``flow`` setting tells SHARPy which solvers to run and in which order. You may be stranged by the presence of the ``DynamicCoupled`` solver but it is necessary to give an initial speed to the structure. This will allow proper linearisation of the structural and rigid body equations." ] }, { "cell_type": "code", "execution_count": 5, "metadata": {}, "outputs": [], "source": [ "flow = ['BeamLoader',\n", " 'AerogridLoader',\n", " 'StaticTrim',\n", " 'BeamPlot',\n", " 'AerogridPlot',\n", " 'AeroForcesCalculator',\n", " 'DynamicCoupled',\n", " 'Modal',\n", " 'LinearAssembler',\n", " 'AsymptoticStability',\n", " ]" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### SHARPy Settings" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "settings = dict()\n", "settings['SHARPy'] = {'case': ws.case_name,\n", " 'route': ws.case_route,\n", " 'flow': flow,\n", " 'write_screen': 'on',\n", " 'write_log': 'on',\n", " 'log_folder': './output/',\n", " 'log_file': ws.case_name + '.log'}\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Loaders" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "settings['BeamLoader'] = {'unsteady': 'off',\n", " 'orientation': algebra.euler2quat(np.array([ws.roll,\n", " ws.alpha,\n", " ws.beta]))}\n", "settings['AerogridLoader'] = {'unsteady': 'off',\n", " 'aligned_grid': 'on',\n", " 'mstar': int(ws.M * ws.Mstarfactor),\n", " 'freestream_dir': ['1', '0', '0'],\n", " 'control_surface_deflection': [''],\n", " 'wake_shape_generator': 'StraightWake',\n", " 'wake_shape_generator_input': {'u_inf': ws.u_inf,\n", " 'u_inf_direction': ['1', '0', '0'],\n", " 'dt': ws.dt}}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### StaticCoupled Solver" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "settings['StaticCoupled'] = {'print_info': 'on',\n", " 'structural_solver': 'NonLinearStatic',\n", " 'structural_solver_settings': {'print_info': 'off',\n", " 'max_iterations': 200,\n", " 'num_load_steps': 1,\n", " 'delta_curved': 1e-5,\n", " 'min_delta': ws.tolerance,\n", " 'gravity_on': 'on',\n", " 'gravity': 9.81},\n", " 'aero_solver': 'StaticUvlm',\n", " 'aero_solver_settings': {'print_info': 'on',\n", " 'horseshoe': ws.horseshoe,\n", " 'num_cores': 4,\n", " 'n_rollup': int(0),\n", " 'rollup_dt': ws.dt,\n", " 'rollup_aic_refresh': 1,\n", " 'rollup_tolerance': 1e-4,\n", " 'velocity_field_generator': 'SteadyVelocityField',\n", " 'velocity_field_input': {'u_inf': ws.u_inf,\n", " 'u_inf_direction': [1., 0, 0]},\n", " 'rho': ws.rho},\n", " 'max_iter': 200,\n", " 'n_load_steps': 1,\n", " 'tolerance': ws.tolerance,\n", " 'relaxation_factor': 0.2}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Trim solver" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [], "source": [ "settings['StaticTrim'] = {'solver': 'StaticCoupled',\n", " 'solver_settings': settings['StaticCoupled'],\n", " 'thrust_nodes': ws.thrust_nodes,\n", " 'initial_alpha': ws.alpha,\n", " 'initial_deflection': ws.cs_deflection,\n", " 'initial_thrust': ws.thrust,\n", " 'max_iter': 200,\n", " 'fz_tolerance': 1e-2,\n", " 'fx_tolerance': 1e-2,\n", " 'm_tolerance': 1e-2}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Nonlinear Equilibrium Post-process" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "settings['AerogridPlot'] = {\n", " 'include_rbm': 'off',\n", " 'include_applied_forces': 'on',\n", " 'minus_m_star': 0,\n", " 'u_inf': ws.u_inf\n", " }\n", "settings['AeroForcesCalculator'] = {\n", " 'write_text_file': 'off',\n", " 'text_file_name': ws.case_name + '_aeroforces.csv',\n", " 'screen_output': 'on',\n", " 'coefficients': True,\n", " 'q_ref': 0.5 * ws.rho * ws.u_inf ** 2,\n", " 'S_ref': 12.809,\n", " }\n", "\n", "settings['BeamPlot'] = {\n", " 'include_rbm': 'on',\n", " 'include_applied_forces': 'on',\n", " 'include_FoR': 'on'}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### DynamicCoupled Solver\n", "\n", "As mentioned before, a single time step of ``DynamicCoupled`` is required to give the structure the velocity required for the linearisation of the rigid body equations to be correct. Hence `n_time_steps = 1`" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [], "source": [ "struct_solver_settings = {'print_info': 'off',\n", " 'initial_velocity_direction': [-1., 0., 0.],\n", " 'max_iterations': 950,\n", " 'delta_curved': 1e-6,\n", " 'min_delta': ws.tolerance,\n", " 'newmark_damp': 5e-3,\n", " 'gravity_on': True,\n", " 'gravity': 9.81,\n", " 'num_steps': ws.n_tstep,\n", " 'dt': ws.dt,\n", " 'initial_velocity': ws.u_inf * 1}\n", "\n", "step_uvlm_settings = {'print_info': 'on',\n", " 'num_cores': 4,\n", " 'convection_scheme': ws.wake_type,\n", " 'velocity_field_generator': 'SteadyVelocityField',\n", " 'velocity_field_input': {'u_inf': ws.u_inf * 0,\n", " 'u_inf_direction': [1., 0., 0.]},\n", " 'rho': ws.rho,\n", " 'n_time_steps': ws.n_tstep,\n", " 'dt': ws.dt,\n", " 'gamma_dot_filtering': 3}\n", "\n", "settings['DynamicCoupled'] = {'print_info': 'on',\n", " 'structural_solver': 'NonLinearDynamicCoupledStep',\n", " 'structural_solver_settings': struct_solver_settings,\n", " 'aero_solver': 'StepUvlm',\n", " 'aero_solver_settings': step_uvlm_settings,\n", " 'fsi_substeps': 200,\n", " 'fsi_tolerance': ws.fsi_tolerance,\n", " 'relaxation_factor': ws.relaxation_factor,\n", " 'minimum_steps': 1,\n", " 'relaxation_steps': 150,\n", " 'final_relaxation_factor': 0.5,\n", " 'n_time_steps': 1,\n", " 'dt': ws.dt,\n", " 'include_unsteady_force_contribution': 'off',\n", " }" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Modal Solver Settings" ] }, { "cell_type": "code", "execution_count": 12, "metadata": {}, "outputs": [], "source": [ "settings['Modal'] = {'print_info': True,\n", " 'use_undamped_modes': True,\n", " 'NumLambda': 30,\n", " 'rigid_body_modes': True,\n", " 'write_modes_vtk': 'on',\n", " 'continuous_eigenvalues': 'off',\n", " 'dt': ws.dt,\n", " 'plot_eigenvalues': False,\n", " 'rigid_modes_cg': False}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Linear Assembler Settings\n", "\n", "Note that for the assembly of the linear system, we replace the parametrisation of the orientation with Euler angles instead of quaternions." ] }, { "cell_type": "code", "execution_count": 13, "metadata": {}, "outputs": [], "source": [ "settings['LinearAssembler'] = {'linear_system': 'LinearAeroelastic',\n", " 'linear_system_settings': {\n", " 'beam_settings': {'modal_projection': 'off',\n", " 'inout_coords': 'modes',\n", " 'discrete_time': True,\n", " 'newmark_damp': 0.5e-2,\n", " 'discr_method': 'newmark',\n", " 'dt': ws.dt,\n", " 'proj_modes': 'undamped',\n", " 'num_modes': 9,\n", " 'print_info': 'on',\n", " 'gravity': 'on',\n", " 'remove_dofs': []},\n", " 'aero_settings': {'dt': ws.dt,\n", " 'integr_order': 2,\n", " 'density': ws.rho,\n", " 'remove_predictor': 'off',\n", " 'use_sparse': 'off',\n", " 'remove_inputs': ['u_gust']},\n", " 'track_body': 'on',\n", " 'use_euler': 'on',\n", " }}" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Asymptotic Stability Post-processor" ] }, { "cell_type": "code", "execution_count": 14, "metadata": {}, "outputs": [], "source": [ "settings['AsymptoticStability'] = {\n", " 'print_info': 'on',\n", " 'frequency_cutoff': 0,\n", " 'export_eigenvalues': 'on',\n", " 'num_evals': 100,\n", " }" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Write solver file" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [], "source": [ "config = configobj.ConfigObj()\n", "np.set_printoptions(precision=16)\n", "file_name = ws.case_route + '/' + ws.case_name + '.sharpy'\n", "config.filename = file_name\n", "for k, v in settings.items():\n", " config[k] = v\n", "config.write()" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Run Simulation" ] }, { "cell_type": "code", "execution_count": 16, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "--------------------------------------------------------------------------------\u001b[0m\n", " ###### ## ## ### ######## ######## ## ##\u001b[0m\n", " ## ## ## ## ## ## ## ## ## ## ## ##\u001b[0m\n", " ## ## ## ## ## ## ## ## ## ####\u001b[0m\n", " ###### ######### ## ## ######## ######## ##\u001b[0m\n", " ## ## ## ######### ## ## ## ##\u001b[0m\n", " ## ## ## ## ## ## ## ## ## ##\u001b[0m\n", " ###### ## ## ## ## ## ## ## ##\u001b[0m\n", "--------------------------------------------------------------------------------\u001b[0m\n", "Aeroelastics Lab, Aeronautics Department.\u001b[0m\n", " Copyright (c), Imperial College London.\u001b[0m\n", " All rights reserved.\u001b[0m\n", " License available at https://github.com/imperialcollegelondon/sharpy\u001b[0m\n", "\u001b[36mRunning SHARPy from /home/ng213/2TB/pazy_code/pazy-sharpy/lib/sharpy/docs/source/content/example_notebooks\u001b[0m\n", "\u001b[36mSHARPy being run is in /home/ng213/2TB/pazy_code/pazy-sharpy/lib/sharpy\u001b[0m\n", "\u001b[36mThe branch being run is dev_setting_error\u001b[0m\n", "\u001b[36mThe version and commit hash are: v1.2.1-344-g0239644-0239644\u001b[0m\n", "SHARPy output folder set\u001b[0m\n", "\u001b[34m\t./output//horten_u_inf2800_M4N11Msf5/\u001b[0m\n", "\u001b[36mGenerating an instance of BeamLoader\u001b[0m\n", "Variable for_pos has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: [0.0, 0, 0]\u001b[0m\n", "\u001b[36mGenerating an instance of AerogridLoader\u001b[0m\n", "Variable control_surface_deflection_generator_settings has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: {}\u001b[0m\n", "Variable dx1 has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: -1.0\u001b[0m\n", "Variable ndx1 has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 1\u001b[0m\n", "Variable r has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 1.0\u001b[0m\n", "Variable dxmax has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: -1.0\u001b[0m\n", "\u001b[34mThe aerodynamic grid contains 4 surfaces\u001b[0m\n", "\u001b[34m Surface 0, M=4, N=2\u001b[0m\n", " Wake 0, M=20, N=2\u001b[0m\n", "\u001b[34m Surface 1, M=4, N=22\u001b[0m\n", " Wake 1, M=20, N=22\u001b[0m\n", "\u001b[34m Surface 2, M=4, N=2\u001b[0m\n", " Wake 2, M=20, N=2\u001b[0m\n", "\u001b[34m Surface 3, M=4, N=22\u001b[0m\n", " Wake 3, M=20, N=22\u001b[0m\n", " In total: 192 bound panels\u001b[0m\n", " In total: 960 wake panels\u001b[0m\n", " Total number of panels = 1152\u001b[0m\n", "\u001b[36mGenerating an instance of StaticTrim\u001b[0m\n", "Variable print_info has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: True\u001b[0m\n", "Variable tail_cs_index has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0\u001b[0m\n", "Variable initial_angle_eps has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0.05\u001b[0m\n", "Variable initial_thrust_eps has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 2.0\u001b[0m\n", "Variable relaxation_factor has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0.2\u001b[0m\n", "Variable save_info has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "\u001b[36mGenerating an instance of StaticCoupled\u001b[0m\n", "Variable correct_forces_method has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: \u001b[0m\n", "Variable runtime_generators has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: {}\u001b[0m\n", "\u001b[36mGenerating an instance of NonLinearStatic\u001b[0m\n", "Variable newmark_damp has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0.0001\u001b[0m\n", "Variable gravity_dir has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: [0.0, 0.0, 1.0]\u001b[0m\n", "Variable relaxation_factor has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0.3\u001b[0m\n", "Variable dt has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0.01\u001b[0m\n", "Variable num_steps has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 500\u001b[0m\n", "Variable initial_position has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: [0. 0. 0.]\u001b[0m\n", "\u001b[36mGenerating an instance of StaticUvlm\u001b[0m\n", "Variable iterative_solver has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "Variable iterative_tol has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0.0001\u001b[0m\n", "Variable iterative_precond has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "Variable cfl1 has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: True\u001b[0m\n", "Variable vortex_radius has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 1e-06\u001b[0m\n", "Variable vortex_radius_wake_ind has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 1e-06\u001b[0m\n", "Variable rbm_vel_g has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: [0.0, 0.0, 0.0, 0.0, 0.0, 0.0]\u001b[0m\n", "Variable centre_rot_g has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: [0.0, 0.0, 0.0]\u001b[0m\n", "\u001b[0m\n", "\u001b[0m\n", "\u001b[0m\n", "|=====|=====|============|==========|==========|==========|==========|==========|==========|\u001b[0m\n", "|iter |step | log10(res) | Fx | Fy | Fz | Mx | My | Mz |\u001b[0m\n", "|=====|=====|============|==========|==========|==========|==========|==========|==========|\u001b[0m\n", "\u001b[0m\n", "\u001b[0m\n", "\u001b[0m\n", "|==========|==========|==========|==========|==========|==========|==========|==========|==========|==========|\u001b[0m\n", "| iter |alpha[deg]|elev[deg] | thrust | Fx | Fy | Fz | Mx | My | Mz |\u001b[0m\n", "|==========|==========|==========|==========|==========|==========|==========|==========|==========|==========|\u001b[0m\n", "| 0 | 0 | 0.00000 | -0.1051 | -0.0000 | 0.0598 | -0.0000 | 1.0837 | -0.0000 |\u001b[0m\n", "| 1 | 0 | -7.62384 | -0.1284 | -0.0000 | 0.1276 | -0.0000 | 0.0045 | -0.0000 |\u001b[0m\n", "| 2 | 0 | -8.33392 | -0.1190 | -0.0000 | 0.0397 | -0.0000 | -0.0774 | -0.0000 |\u001b[0m\n", "| 3 | 0 | -9.30379 | -0.1133 | 0.0000 | 0.0011 | -0.0000 | -0.0070 | 0.0000 |\u001b[0m\n", "| 4 | 0 | -10.71602 | -0.1136 | -0.0000 | 0.0032 | 0.0000 | -0.0100 | -0.0000 |\u001b[0m\n", "| 5 | 0 | -10.88827 | -0.1138 | -0.0000 | 0.0043 | 0.0000 | -0.0119 | -0.0000 |\u001b[0m\n", "| 6 | 0 | -11.66331 | -0.1138 | -0.0000 | 0.0042 | 0.0000 | -0.0116 | -0.0000 |\u001b[0m\n", "| 7 | 0 | -12.88496 | -0.1138 | -0.0000 | 0.0041 | 0.0000 | -0.0116 | 0.0000 |\u001b[0m\n", "| 0 | 4.5135 | 0.1814 | 5.5129 | -0.1138 | -0.0000 | 0.0041 | 0.0000 | -0.0116 | 0.0000 |\u001b[0m\n", "| 0 | 0 | 0.00000 |-116.9870 | -0.0000 | 994.8063 | -0.0000 |-882.4104 | -0.0000 |\u001b[0m\n", "| 1 | 0 | -5.79178 | -72.3841 | 0.0000 | 944.6140 | 0.0000 |-802.5912 | 0.0000 |\u001b[0m\n", "| 2 | 0 | -6.63730 | -62.4378 | 0.0000 | 937.6662 | 0.0000 |-792.1259 | -0.0000 |\u001b[0m\n", "| 3 | 0 | -7.22937 | -62.8923 | -0.0000 | 939.7866 | -0.0000 |-795.7093 | -0.0000 |\u001b[0m\n", "| 4 | 0 | -8.65323 | -62.8757 | -0.0000 | 939.7100 | -0.0000 |-795.5764 | -0.0000 |\u001b[0m\n", "| 5 | 0 | -8.81438 | -62.8640 | -0.0000 | 939.6554 | -0.0000 |-795.4837 | -0.0000 |\u001b[0m\n", "| 6 | 0 | -9.59386 | -62.8660 | 0.0000 | 939.6644 | -0.0000 |-795.4991 | 0.0000 |\u001b[0m\n", "| 7 | 0 | -10.80422 | -62.8661 | 0.0000 | 939.6650 | -0.0000 |-795.5000 | 0.0000 |\u001b[0m\n", "| 8 | 0 | -11.01365 | -62.8660 | -0.0000 | 939.6647 | -0.0000 |-795.4994 | -0.0000 |\u001b[0m\n", "| 9 | 0 | -12.15197 | -62.8660 | 0.0000 | 939.6647 | -0.0000 |-795.4995 | 0.0000 |\u001b[0m\n", "| 0 | 7.3783 | -2.6834 | 5.5129 | -62.8660 | 0.0000 | 939.6647 | -0.0000 |-795.4995 | 0.0000 |\u001b[0m\n", "| 0 | 0 | 0.00000 | -32.9132 | -0.0000 | 371.4715 | -0.0000 |-902.7953 | -0.0000 |\u001b[0m\n", "| 1 | 0 | -5.48409 | -8.7241 | -0.0000 | 298.8484 | 0.0000 |-777.2938 | 0.0000 |\u001b[0m\n", "| 2 | 0 | -6.39387 | -4.0957 | -0.0000 | 289.8092 | -0.0000 |-761.4855 | -0.0000 |\u001b[0m\n", "| 3 | 0 | -6.85613 | -4.6263 | -0.0000 | 293.2407 | -0.0000 |-767.3376 | 0.0000 |\u001b[0m\n", "| 4 | 0 | -8.25962 | -4.6052 | -0.0000 | 293.1048 | -0.0000 |-767.1066 | 0.0000 |\u001b[0m\n", "| 5 | 0 | -8.44065 | -4.5914 | -0.0000 | 293.0156 | 0.0000 |-766.9545 | 0.0000 |\u001b[0m\n", "| 6 | 0 | -9.20968 | -4.5937 | -0.0000 | 293.0308 | -0.0000 |-766.9804 | -0.0000 |\u001b[0m\n", "| 7 | 0 | -10.44736 | -4.5939 | -0.0000 | 293.0316 | -0.0000 |-766.9819 | -0.0000 |\u001b[0m\n", "| 8 | 0 | -10.63000 | -4.5938 | -0.0000 | 293.0311 | 0.0000 |-766.9809 | -0.0000 |\u001b[0m\n", "| 9 | 0 | -11.74670 | -4.5938 | 0.0000 | 293.0311 | 0.0000 |-766.9810 | 0.0000 |\u001b[0m\n", "| 10 | 0 | -12.29943 | -4.5938 | -0.0000 | 293.0311 | 0.0000 |-766.9810 | 0.0000 |\u001b[0m\n", "| 0 | 4.5135 | 3.0462 | 5.5129 | -4.5938 | -0.0000 | 293.0311 | 0.0000 |-766.9810 | 0.0000 |\u001b[0m\n", "| 0 | 0 | 0.00000 | -4.1051 | -0.0000 | 0.0598 | -0.0000 | 1.0834 | -0.0000 |\u001b[0m\n", "| 1 | 0 | -7.62384 | -4.1284 | -0.0000 | 0.1276 | -0.0000 | 0.0042 | -0.0000 |\u001b[0m\n", "| 2 | 0 | -8.33392 | -4.1190 | -0.0000 | 0.0397 | 0.0000 | -0.0778 | -0.0000 |\u001b[0m\n", "| 3 | 0 | -9.30379 | -4.1133 | -0.0000 | 0.0011 | 0.0000 | -0.0074 | -0.0000 |\u001b[0m\n", "| 4 | 0 | -10.71602 | -4.1136 | 0.0000 | 0.0032 | -0.0000 | -0.0104 | 0.0000 |\u001b[0m\n", "| 5 | 0 | -10.88827 | -4.1138 | 0.0000 | 0.0043 | 0.0000 | -0.0123 | 0.0000 |\u001b[0m\n", "| 6 | 0 | -11.66331 | -4.1138 | -0.0000 | 0.0042 | 0.0000 | -0.0119 | -0.0000 |\u001b[0m\n", "| 7 | 0 | -12.88496 | -4.1138 | -0.0000 | 0.0041 | 0.0000 | -0.0119 | -0.0000 |\u001b[0m\n", "| 0 | 4.5135 | 0.1814 | 7.5129 | -4.1138 | -0.0000 | 0.0041 | 0.0000 | -0.0119 | -0.0000 |\u001b[0m\n", "| 0 | 0 | 0.00000 | 0.0095 | -0.0000 | 0.0498 | 0.0000 | 1.1013 | -0.0000 |\u001b[0m\n", "| 1 | 0 | -7.62357 | -0.0140 | 0.0000 | 0.1189 | 0.0000 | 0.0198 | 0.0000 |\u001b[0m\n", "| 2 | 0 | -8.33375 | -0.0046 | -0.0000 | 0.0312 | -0.0000 | -0.0624 | 0.0000 |\u001b[0m\n", "| 3 | 0 | -9.30318 | 0.0010 | -0.0000 | -0.0075 | 0.0000 | 0.0081 | -0.0000 |\u001b[0m\n", "| 4 | 0 | -10.71542 | 0.0007 | -0.0000 | -0.0054 | 0.0000 | 0.0051 | 0.0000 |\u001b[0m\n", "| 5 | 0 | -10.88766 | 0.0006 | -0.0000 | -0.0043 | 0.0000 | 0.0032 | -0.0000 |\u001b[0m\n", "| 6 | 0 | -11.66271 | 0.0006 | 0.0000 | -0.0044 | -0.0000 | 0.0035 | 0.0000 |\u001b[0m\n", "| 7 | 0 | -12.88441 | 0.0006 | -0.0000 | -0.0045 | 0.0000 | 0.0035 | -0.0000 |\u001b[0m\n", "| 1 | 4.5135 | 0.1814 | 5.4560 | 0.0006 | -0.0000 | -0.0045 | 0.0000 | 0.0035 | -0.0000 |\u001b[0m\n", "\u001b[36mGenerating an instance of BeamPlot\u001b[0m\n", "Variable include_applied_moments has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: True\u001b[0m\n", "Variable name_prefix has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: \u001b[0m\n", "Variable output_rbm has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: True\u001b[0m\n", "\u001b[34m...Finished\u001b[0m\n", "\u001b[36mGenerating an instance of AerogridPlot\u001b[0m\n", "Variable include_forward_motion has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "Variable include_unsteady_applied_forces has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "Variable name_prefix has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: \u001b[0m\n", "Variable dt has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0.0\u001b[0m\n", "Variable include_velocities has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "Variable include_incidence_angle has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "Variable num_cores has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 1\u001b[0m\n", "Variable vortex_radius has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 1e-06\u001b[0m\n", "\u001b[34m...Finished\u001b[0m\n", "\u001b[36mGenerating an instance of AeroForcesCalculator\u001b[0m\n", "--------------------------------------------------------------------------------\u001b[0m\n", "\u001b[34mtstep | fx_g | fy_g | fz_g | Cfx_g | Cfy_g | Cfz_g \u001b[0m\n", "\u001b[34m 0 | 1.088e+01 | -4.476e-13 | 1.835e+03 | 2.124e-03 | -8.740e-17 | 3.583e-01\u001b[0m\n", "\u001b[34m...Finished\u001b[0m\n", "\u001b[36mGenerating an instance of DynamicCoupled\u001b[0m\n", "Variable structural_substeps has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0\u001b[0m\n", "Variable dynamic_relaxation has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "Variable postprocessors has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: []\u001b[0m\n", "Variable postprocessors_settings has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: {}\u001b[0m\n", "Variable controller_id has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: {}\u001b[0m\n", "Variable controller_settings has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: {}\u001b[0m\n", "Variable cleanup_previous_solution has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "Variable steps_without_unsteady_force has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0\u001b[0m\n", "Variable pseudosteps_ramp_unsteady_force has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0\u001b[0m\n", "Variable correct_forces_method has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: \u001b[0m\n", "Variable network_settings has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: {}\u001b[0m\n", "Variable runtime_generators has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: {}\u001b[0m\n", "\u001b[36mGenerating an instance of NonLinearDynamicCoupledStep\u001b[0m\n", "Variable num_load_steps has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 1\u001b[0m\n", "Variable gravity_dir has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: [0.0, 0.0, 1.0]\u001b[0m\n", "Variable relaxation_factor has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0.3\u001b[0m\n", "Variable balancing has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "\u001b[36mGenerating an instance of StepUvlm\u001b[0m\n", "Variable iterative_solver has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "Variable iterative_tol has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0.0001\u001b[0m\n", "Variable iterative_precond has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "Variable cfl1 has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: True\u001b[0m\n", "Variable vortex_radius has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 1e-06\u001b[0m\n", "Variable vortex_radius_wake_ind has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 1e-06\u001b[0m\n", "Variable interp_coords has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0\u001b[0m\n", "Variable filter_method has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0\u001b[0m\n", "Variable interp_method has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0\u001b[0m\n", "Variable yaw_slerp has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0.0\u001b[0m\n", "Variable centre_rot has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: [0.0, 0.0, 0.0]\u001b[0m\n", "Variable quasi_steady has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "\u001b[0m\n", "\u001b[0m\n", "\u001b[0m\n", "|=======|========|======|==============|==============|==============|==============|==============|\u001b[0m\n", "| ts | t | iter | struc ratio | iter time | residual vel | FoR_vel(x) | FoR_vel(z) |\u001b[0m\n", "|=======|========|======|==============|==============|==============|==============|==============|\u001b[0m\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "/home/ng213/2TB/pazy_code/pazy-sharpy/lib/sharpy/sharpy/aero/utils/uvlmlib.py:264: RuntimeWarning: invalid value encountered in true_divide\n", " flightconditions.uinf_direction = np.ctypeslib.as_ctypes(ts_info.u_ext[0][:, 0, 0]/flightconditions.uinf)\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "| 1 | 0.0089 | 3 | 0.877319 | 0.593232 | -10.598250 |-2.791317e+01 |-2.203426e+00 |\u001b[0m\n", "\u001b[34m...Finished\u001b[0m\n", "\u001b[36mGenerating an instance of Modal\u001b[0m\n", "Variable keep_linear_matrices has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: True\u001b[0m\n", "Variable write_dat has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: True\u001b[0m\n", "Variable delta_curved has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0.01\u001b[0m\n", "Variable max_rotation_deg has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 15.0\u001b[0m\n", "Variable max_displacement has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 0.15\u001b[0m\n", "Variable use_custom_timestep has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: -1\u001b[0m\n", "Structural eigenvalues\u001b[0m\n", "\u001b[0m\n", "\u001b[0m\n", "\u001b[0m\n", "|==============|==============|==============|==============|==============|==============|==============|\u001b[0m\n", "| mode | eval_real | eval_imag | freq_n (Hz) | freq_d (Hz) | damping | period (s) |\u001b[0m\n", "|==============|==============|==============|==============|==============|==============|==============|\u001b[0m\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "/home/ng213/2TB/pazy_code/pazy-sharpy/lib/sharpy/sharpy/solvers/modal.py:265: UserWarning: Projecting a system with damping on undamped modal shapes\n", " warnings.warn('Projecting a system with damping on undamped modal shapes')\n" ] }, { "name": "stdout", "output_type": "stream", "text": [ "| 0 | -0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 1 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 2 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 3 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 4 | -0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 5 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 6 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 7 | -0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 8 | -0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 9 | -0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 10 | 0.000000 | 28.293939 | 4.503120 | 4.503120 | -0.000000 | 0.222068 |\u001b[0m\n", "| 11 | 0.000000 | 29.271318 | 4.658675 | 4.658675 | -0.000000 | 0.214653 |\u001b[0m\n", "| 12 | 0.000000 | 54.780234 | 8.718545 | 8.718545 | -0.000000 | 0.114698 |\u001b[0m\n", "| 13 | 0.000000 | 58.999779 | 9.390106 | 9.390106 | -0.000000 | 0.106495 |\u001b[0m\n", "| 14 | 0.000000 | 70.520741 | 11.223724 | 11.223724 | -0.000000 | 0.089097 |\u001b[0m\n", "| 15 | 0.000000 | 76.917111 | 12.241738 | 12.241738 | -0.000000 | 0.081688 |\u001b[0m\n", "| 16 | 0.000000 | 87.324076 | 13.898058 | 13.898058 | -0.000000 | 0.071952 |\u001b[0m\n", "| 17 | 0.000000 | 108.035577 | 17.194396 | 17.194396 | -0.000000 | 0.058158 |\u001b[0m\n", "| 18 | 0.000000 | 119.692139 | 19.049596 | 19.049596 | -0.000000 | 0.052495 |\u001b[0m\n", "| 19 | 0.000000 | 133.495187 | 21.246419 | 21.246419 | -0.000000 | 0.047067 |\u001b[0m\n", "| 20 | 0.000000 | 134.444788 | 21.397553 | 21.397553 | -0.000000 | 0.046734 |\u001b[0m\n", "| 21 | 0.000000 | 151.060442 | 24.042016 | 24.042016 | -0.000000 | 0.041594 |\u001b[0m\n", "| 22 | 0.000000 | 159.369020 | 25.364367 | 25.364367 | -0.000000 | 0.039425 |\u001b[0m\n", "| 23 | 0.000000 | 171.256102 | 27.256255 | 27.256255 | -0.000000 | 0.036689 |\u001b[0m\n", "| 24 | 0.000000 | 173.895881 | 27.676389 | 27.676389 | -0.000000 | 0.036132 |\u001b[0m\n", "| 25 | 0.000000 | 199.016557 | 31.674469 | 31.674469 | -0.000000 | 0.031571 |\u001b[0m\n", "| 26 | 0.000000 | 205.412581 | 32.692428 | 32.692428 | -0.000000 | 0.030588 |\u001b[0m\n", "| 27 | 0.000000 | 205.419531 | 32.693534 | 32.693534 | -0.000000 | 0.030587 |\u001b[0m\n", "| 28 | 0.000000 | 223.563796 | 35.581283 | 35.581283 | -0.000000 | 0.028105 |\u001b[0m\n", "| 29 | 0.000000 | 227.924750 | 36.275351 | 36.275351 | -0.000000 | 0.027567 |\u001b[0m\n", "\u001b[36mGenerating an instance of LinearAssembler\u001b[0m\n", "Variable linearisation_tstep has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: -1\u001b[0m\n", "Variable modal_tstep has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: -1\u001b[0m\n", "Variable inout_coordinates has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: \u001b[0m\n", "Variable retain_inputs has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: []\u001b[0m\n", "Variable retain_outputs has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: []\u001b[0m\n", "\u001b[36mGenerating an instance of LinearAeroelastic\u001b[0m\n", "Variable uvlm_filename has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: \u001b[0m\n", "\u001b[36mGenerating an instance of LinearUVLM\u001b[0m\n", "Variable ScalingDict has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: {}\u001b[0m\n", "Variable gust_assembler has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: \u001b[0m\n", "Variable rom_method has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: []\u001b[0m\n", "Variable rom_method_settings has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: {}\u001b[0m\n", "Variable vortex_radius has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 1e-06\u001b[0m\n", "Variable cfl1 has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: True\u001b[0m\n", "Variable length has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 1.0\u001b[0m\n", "Variable speed has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 1.0\u001b[0m\n", "Variable density has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 1.0\u001b[0m\n", "Variable velocity_field_generator has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: SteadyVelocityField\u001b[0m\n", "Variable velocity_field_input has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: {}\u001b[0m\n", "Variable physical_model has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: True\u001b[0m\n", "Variable track_body has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "Variable track_body_number has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: -1\u001b[0m\n", "Initialising Static linear UVLM solver class...\u001b[0m\n", "\t\t\t...done in 0.39 sec\u001b[0m\n", "State-space realisation of UVLM equations started...\u001b[0m\n", "\u001b[34mComputing wake propagation matrix with CFL1=True\u001b[0m\n", "\u001b[34m\tstate-space model produced in form:\u001b[0m\n", "\u001b[34m\tx_{n+1} = A x_{n} + Bp u_{n+1}\u001b[0m\n", "\t\t\t...done in 2.43 sec\u001b[0m\n", "\u001b[36mGenerating an instance of LinearBeam\u001b[0m\n", "Variable remove_sym_modes has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "Warning, projecting system with damping onto undamped modes\u001b[0m\n", "\u001b[0m\n", "Linearising gravity terms...\u001b[0m\n", "\u001b[34m\tM = 187.12 kg\u001b[0m\n", "\u001b[34m\tX_CG A -> 1.19 0.00 0.01\u001b[0m\n", "\u001b[36mNode 1 \t-> B 0.000 -0.089 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 0.089 0.206 0.000\u001b[0m\n", "\u001b[36m\t\t\t-> G 0.089 0.206 -0.007\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 2.6141\u001b[0m\n", "\u001b[36mNode 2 \t-> B -0.010 -0.019 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 0.019 0.403 0.000\u001b[0m\n", "\u001b[36m\t\t\t-> G 0.019 0.403 -0.001\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 7.3672\u001b[0m\n", "\u001b[36mNode 3 \t-> B -0.019 -0.087 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 0.234 0.800 0.000\u001b[0m\n", "\u001b[36m\t\t\t-> G 0.234 0.800 -0.018\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 5.8780\u001b[0m\n", "\u001b[36mNode 4 \t-> B -0.019 -0.084 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 0.390 1.238 0.001\u001b[0m\n", "\u001b[36m\t\t\t-> G 0.389 1.238 -0.030\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 2.8288\u001b[0m\n", "\u001b[36mNode 5 \t-> B -0.018 -0.081 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 0.546 1.676 0.001\u001b[0m\n", "\u001b[36m\t\t\t-> G 0.544 1.676 -0.041\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 5.4372\u001b[0m\n", "\u001b[36mNode 6 \t-> B -0.017 -0.078 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 0.702 2.113 0.002\u001b[0m\n", "\u001b[36m\t\t\t-> G 0.700 2.113 -0.053\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 2.6084\u001b[0m\n", "\u001b[36mNode 7 \t-> B -0.016 -0.074 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 0.857 2.551 0.003\u001b[0m\n", "\u001b[36m\t\t\t-> G 0.855 2.551 -0.064\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 4.9963\u001b[0m\n", "\u001b[36mNode 8 \t-> B -0.016 -0.071 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 1.013 2.988 0.004\u001b[0m\n", "\u001b[36m\t\t\t-> G 1.010 2.988 -0.076\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 2.3879\u001b[0m\n", "\u001b[36mNode 9 \t-> B -0.015 -0.068 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 1.169 3.426 0.005\u001b[0m\n", "\u001b[36m\t\t\t-> G 1.166 3.426 -0.087\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 4.5555\u001b[0m\n", "\u001b[36mNode 10 \t-> B -0.014 -0.065 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 1.325 3.863 0.006\u001b[0m\n", "\u001b[36m\t\t\t-> G 1.321 3.863 -0.098\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 2.1675\u001b[0m\n", "\u001b[36mNode 11 \t-> B -0.013 -0.061 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 1.480 4.301 0.007\u001b[0m\n", "\u001b[36m\t\t\t-> G 1.476 4.301 -0.109\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 4.1146\u001b[0m\n", "\u001b[36mNode 12 \t-> B -0.013 -0.058 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 1.636 4.739 0.009\u001b[0m\n", "\u001b[36m\t\t\t-> G 1.632 4.739 -0.120\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 1.9471\u001b[0m\n", "\u001b[36mNode 13 \t-> B -0.012 -0.055 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 1.792 5.176 0.010\u001b[0m\n", "\u001b[36m\t\t\t-> G 1.787 5.176 -0.131\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 3.6738\u001b[0m\n", "\u001b[36mNode 14 \t-> B -0.011 -0.052 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 1.948 5.614 0.011\u001b[0m\n", "\u001b[36m\t\t\t-> G 1.943 5.614 -0.142\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 1.7267\u001b[0m\n", "\u001b[36mNode 15 \t-> B -0.011 -0.048 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 2.104 6.052 0.012\u001b[0m\n", "\u001b[36m\t\t\t-> G 2.098 6.052 -0.153\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 3.2329\u001b[0m\n", "\u001b[36mNode 16 \t-> B -0.010 -0.045 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 2.260 6.489 0.014\u001b[0m\n", "\u001b[36m\t\t\t-> G 2.254 6.489 -0.164\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 1.5062\u001b[0m\n", "\u001b[36mNode 17 \t-> B -0.009 -0.042 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 2.415 6.927 0.015\u001b[0m\n", "\u001b[36m\t\t\t-> G 2.409 6.927 -0.175\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 2.7921\u001b[0m\n", "\u001b[36mNode 18 \t-> B -0.008 -0.039 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 2.571 7.364 0.016\u001b[0m\n", "\u001b[36m\t\t\t-> G 2.564 7.364 -0.186\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 1.2858\u001b[0m\n", "\u001b[36mNode 19 \t-> B -0.008 -0.035 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 2.727 7.802 0.017\u001b[0m\n", "\u001b[36m\t\t\t-> G 2.720 7.802 -0.197\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 2.3512\u001b[0m\n", "\u001b[36mNode 20 \t-> B -0.007 -0.032 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 2.883 8.239 0.019\u001b[0m\n", "\u001b[36m\t\t\t-> G 2.875 8.239 -0.208\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 1.0654\u001b[0m\n", "\u001b[36mNode 21 \t-> B -0.006 -0.028 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 3.038 8.677 0.020\u001b[0m\n", "\u001b[36m\t\t\t-> G 3.030 8.677 -0.219\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 1.9104\u001b[0m\n", "\u001b[36mNode 22 \t-> B -0.006 -0.026 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 3.194 9.114 0.021\u001b[0m\n", "\u001b[36m\t\t\t-> G 3.186 9.114 -0.230\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 0.8450\u001b[0m\n", "\u001b[36mNode 23 \t-> B -0.005 -0.022 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 3.350 9.552 0.023\u001b[0m\n", "\u001b[36m\t\t\t-> G 3.341 9.552 -0.241\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 1.4695\u001b[0m\n", "\u001b[36mNode 24 \t-> B -0.005 -0.022 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 3.508 9.988 0.024\u001b[0m\n", "\u001b[36m\t\t\t-> G 3.499 9.988 -0.252\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 0.3674\u001b[0m\n", "\u001b[36mNode 25 \t-> B 0.000 0.089 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 0.089 -0.206 0.000\u001b[0m\n", "\u001b[36m\t\t\t-> G 0.089 -0.206 -0.007\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 2.6141\u001b[0m\n", "\u001b[36mNode 26 \t-> B -0.010 0.019 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 0.019 -0.403 0.000\u001b[0m\n", "\u001b[36m\t\t\t-> G 0.019 -0.403 -0.001\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 7.3672\u001b[0m\n", "\u001b[36mNode 27 \t-> B -0.019 0.087 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 0.234 -0.800 0.000\u001b[0m\n", "\u001b[36m\t\t\t-> G 0.234 -0.800 -0.018\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 5.8780\u001b[0m\n", "\u001b[36mNode 28 \t-> B -0.019 0.084 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 0.390 -1.238 0.001\u001b[0m\n", "\u001b[36m\t\t\t-> G 0.389 -1.238 -0.030\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 2.8288\u001b[0m\n", "\u001b[36mNode 29 \t-> B -0.018 0.081 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 0.546 -1.676 0.001\u001b[0m\n", "\u001b[36m\t\t\t-> G 0.544 -1.676 -0.041\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 5.4372\u001b[0m\n", "\u001b[36mNode 30 \t-> B -0.017 0.078 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 0.702 -2.113 0.002\u001b[0m\n", "\u001b[36m\t\t\t-> G 0.700 -2.113 -0.053\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 2.6084\u001b[0m\n", "\u001b[36mNode 31 \t-> B -0.016 0.074 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 0.857 -2.551 0.003\u001b[0m\n", "\u001b[36m\t\t\t-> G 0.855 -2.551 -0.064\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 4.9963\u001b[0m\n", "\u001b[36mNode 32 \t-> B -0.016 0.071 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 1.013 -2.988 0.004\u001b[0m\n", "\u001b[36m\t\t\t-> G 1.010 -2.988 -0.076\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 2.3879\u001b[0m\n", "\u001b[36mNode 33 \t-> B -0.015 0.068 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 1.169 -3.426 0.005\u001b[0m\n", "\u001b[36m\t\t\t-> G 1.166 -3.426 -0.087\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 4.5555\u001b[0m\n", "\u001b[36mNode 34 \t-> B -0.014 0.065 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 1.325 -3.863 0.006\u001b[0m\n", "\u001b[36m\t\t\t-> G 1.321 -3.863 -0.098\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 2.1675\u001b[0m\n", "\u001b[36mNode 35 \t-> B -0.013 0.061 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 1.480 -4.301 0.007\u001b[0m\n", "\u001b[36m\t\t\t-> G 1.476 -4.301 -0.109\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 4.1146\u001b[0m\n", "\u001b[36mNode 36 \t-> B -0.013 0.058 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 1.636 -4.739 0.009\u001b[0m\n", "\u001b[36m\t\t\t-> G 1.632 -4.739 -0.120\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 1.9471\u001b[0m\n", "\u001b[36mNode 37 \t-> B -0.012 0.055 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 1.792 -5.176 0.010\u001b[0m\n", "\u001b[36m\t\t\t-> G 1.787 -5.176 -0.131\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 3.6738\u001b[0m\n", "\u001b[36mNode 38 \t-> B -0.011 0.052 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 1.948 -5.614 0.011\u001b[0m\n", "\u001b[36m\t\t\t-> G 1.943 -5.614 -0.142\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 1.7267\u001b[0m\n", "\u001b[36mNode 39 \t-> B -0.011 0.048 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 2.104 -6.052 0.012\u001b[0m\n", "\u001b[36m\t\t\t-> G 2.098 -6.052 -0.153\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 3.2329\u001b[0m\n", "\u001b[36mNode 40 \t-> B -0.010 0.045 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 2.260 -6.489 0.014\u001b[0m\n", "\u001b[36m\t\t\t-> G 2.254 -6.489 -0.164\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 1.5062\u001b[0m\n", "\u001b[36mNode 41 \t-> B -0.009 0.042 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 2.415 -6.927 0.015\u001b[0m\n", "\u001b[36m\t\t\t-> G 2.409 -6.927 -0.175\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 2.7921\u001b[0m\n", "\u001b[36mNode 42 \t-> B -0.008 0.039 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 2.571 -7.364 0.016\u001b[0m\n", "\u001b[36m\t\t\t-> G 2.564 -7.364 -0.186\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 1.2858\u001b[0m\n", "\u001b[36mNode 43 \t-> B -0.008 0.035 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 2.727 -7.802 0.017\u001b[0m\n", "\u001b[36m\t\t\t-> G 2.720 -7.802 -0.197\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 2.3512\u001b[0m\n", "\u001b[36mNode 44 \t-> B -0.007 0.032 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 2.883 -8.239 0.019\u001b[0m\n", "\u001b[36m\t\t\t-> G 2.875 -8.239 -0.208\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 1.0654\u001b[0m\n", "\u001b[36mNode 45 \t-> B -0.006 0.028 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 3.038 -8.677 0.020\u001b[0m\n", "\u001b[36m\t\t\t-> G 3.030 -8.677 -0.219\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 1.9104\u001b[0m\n", "\u001b[36mNode 46 \t-> B -0.006 0.026 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 3.194 -9.114 0.021\u001b[0m\n", "\u001b[36m\t\t\t-> G 3.186 -9.114 -0.230\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 0.8450\u001b[0m\n", "\u001b[36mNode 47 \t-> B -0.005 0.022 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 3.350 -9.552 0.023\u001b[0m\n", "\u001b[36m\t\t\t-> G 3.341 -9.552 -0.241\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 1.4695\u001b[0m\n", "\u001b[36mNode 48 \t-> B -0.005 0.022 -0.000\u001b[0m\n", "\u001b[36m\t\t\t-> A 3.508 -9.988 0.024\u001b[0m\n", "\u001b[36m\t\t\t-> G 3.499 -9.988 -0.252\u001b[0m\n", "\u001b[36m\tNode mass:\u001b[0m\n", "\u001b[36m\t\tMatrix: 0.3674\u001b[0m\n", " Updated the beam C, modal C and K matrices with the terms from the\n", "gravity linearisation\u001b[0m\n", "\u001b[0m\n", "Aeroelastic system assembled:\u001b[0m\n", "\u001b[34m\tAerodynamic states: 1536\u001b[0m\n", "\u001b[34m\tStructural states: 594\u001b[0m\n", "\u001b[34m\tTotal states: 2130\u001b[0m\n", "\u001b[34m\tInputs: 893\u001b[0m\n", "\u001b[34m\tOutputs: 891\u001b[0m\n", "\u001b[34mFinal system is:\u001b[0m\n", "\u001b[34mState-space system\u001b[0m\n", "\u001b[34mStates: 2130\u001b[0m\n", "\u001b[34mInputs: 893\u001b[0m\n", "\u001b[34mOutputs: 891\u001b[0m\n", "\u001b[34m\u001b[0m\n", "\u001b[36mGenerating an instance of AsymptoticStability\u001b[0m\n", "Variable reference_velocity has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: 1.0\u001b[0m\n", "Variable display_root_locus has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "Variable velocity_analysis has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: []\u001b[0m\n", "Variable iterative_eigvals has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: False\u001b[0m\n", "Variable modes_to_plot has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: []\u001b[0m\n", "Variable postprocessors has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: []\u001b[0m\n", "Variable postprocessors_settings has no assigned value in the settings file.\u001b[0m\n", "\u001b[34m will default to the value: {}\u001b[0m\n", "Dynamical System Eigenvalues\u001b[0m\n", "Calculating eigenvalues using direct method\u001b[0m\n", "\u001b[0m\n", "\u001b[0m\n", "\u001b[0m\n", "|==============|==============|==============|==============|==============|==============|==============|\u001b[0m\n", "| mode | eval_real | eval_imag | freq_n (Hz) | freq_d (Hz) | damping | period (s) |\u001b[0m\n", "|==============|==============|==============|==============|==============|==============|==============|\u001b[0m\n", "| 0 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 1 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 2 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 3 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 4 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 5 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 6 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 7 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 8 | 0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 9 | -0.000000 | 0.000000 | 0.000000 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 10 | -0.000884 | -0.321712 | 0.051202 | 0.051202 | 0.002747 | 19.530439 |\u001b[0m\n", "| 11 | -0.000884 | 0.321712 | 0.051202 | 0.051202 | 0.002747 | 19.530439 |\u001b[0m\n", "| 12 | -0.008470 | -0.391287 | 0.062290 | 0.062275 | 0.021642 | 16.057731 |\u001b[0m\n", "| 13 | -0.008470 | 0.391287 | 0.062290 | 0.062275 | 0.021642 | 16.057731 |\u001b[0m\n", "| 14 | -0.022506 | 0.000000 | 0.003582 | 0.000000 | 1.000000 | inf |\u001b[0m\n", "| 15 | -0.064808 | -53.732514 | 8.551801 | 8.551795 | 0.001206 | 0.116935 |\u001b[0m\n", "| 16 | -0.064808 | 53.732514 | 8.551801 | 8.551795 | 0.001206 | 0.116935 |\u001b[0m\n", "| 17 | -0.101946 | 68.319141 | 10.873341 | 10.873329 | 0.001492 | 0.091968 |\u001b[0m\n", "| 18 | -0.101946 | -68.319141 | 10.873341 | 10.873329 | 0.001492 | 0.091968 |\u001b[0m\n", "| 19 | -0.147587 | 83.265282 | 13.252102 | 13.252081 | 0.001772 | 0.075460 |\u001b[0m\n", "| 20 | -0.147587 | -83.265282 | 13.252102 | 13.252081 | 0.001772 | 0.075460 |\u001b[0m\n", "| 21 | -0.248703 | -109.925782 | 17.495276 | 17.495232 | 0.002262 | 0.057158 |\u001b[0m\n", "| 22 | -0.248703 | 109.925782 | 17.495276 | 17.495232 | 0.002262 | 0.057158 |\u001b[0m\n", "| 23 | -0.293472 | 120.387631 | 19.160344 | 19.160287 | 0.002438 | 0.052191 |\u001b[0m\n", "| 24 | -0.293472 | -120.387631 | 19.160344 | 19.160287 | 0.002438 | 0.052191 |\u001b[0m\n", "| 25 | -0.350319 | 132.903280 | 21.152288 | 21.152214 | 0.002636 | 0.047276 |\u001b[0m\n", "| 26 | -0.350319 | -132.903280 | 21.152288 | 21.152214 | 0.002636 | 0.047276 |\u001b[0m\n", "| 27 | -0.376401 | 138.516954 | 22.045739 | 22.045658 | 0.002717 | 0.045360 |\u001b[0m\n", "| 28 | -0.376401 | -138.516954 | 22.045739 | 22.045658 | 0.002717 | 0.045360 |\u001b[0m\n", "| 29 | -0.494445 | 162.714232 | 25.896894 | 25.896774 | 0.003039 | 0.038615 |\u001b[0m\n", "| 30 | -0.494445 | -162.714232 | 25.896894 | 25.896774 | 0.003039 | 0.038615 |\u001b[0m\n", "| 31 | -0.511651 | -166.238306 | 26.457773 | 26.457648 | 0.003078 | 0.037796 |\u001b[0m\n", "| 32 | -0.511651 | 166.238306 | 26.457773 | 26.457648 | 0.003078 | 0.037796 |\u001b[0m\n", "| 33 | -0.559180 | 175.709816 | 27.965227 | 27.965086 | 0.003182 | 0.035759 |\u001b[0m\n", "| 34 | -0.559180 | -175.709816 | 27.965227 | 27.965086 | 0.003182 | 0.035759 |\u001b[0m\n", "| 35 | -0.569755 | 177.873274 | 28.309556 | 28.309411 | 0.003203 | 0.035324 |\u001b[0m\n", "| 36 | -0.569755 | -177.873274 | 28.309556 | 28.309411 | 0.003203 | 0.035324 |\u001b[0m\n", "| 37 | -0.669914 | -197.999020 | 31.512703 | 31.512523 | 0.003383 | 0.031733 |\u001b[0m\n", "| 38 | -0.669914 | 197.999020 | 31.512703 | 31.512523 | 0.003383 | 0.031733 |\u001b[0m\n", "| 39 | -0.678424 | 199.782714 | 31.796590 | 31.796406 | 0.003396 | 0.031450 |\u001b[0m\n", "| 40 | -0.678424 | -199.782714 | 31.796590 | 31.796406 | 0.003396 | 0.031450 |\u001b[0m\n", "| 41 | -0.715684 | 207.440562 | 33.015387 | 33.015191 | 0.003450 | 0.030289 |\u001b[0m\n", "| 42 | -0.715684 | -207.440562 | 33.015387 | 33.015191 | 0.003450 | 0.030289 |\u001b[0m\n", "| 43 | -0.721193 | -208.623018 | 33.203583 | 33.203385 | 0.003457 | 0.030117 |\u001b[0m\n", "| 44 | -0.721193 | 208.623018 | 33.203583 | 33.203385 | 0.003457 | 0.030117 |\u001b[0m\n", "| 45 | -0.796838 | -224.809928 | 35.779836 | 35.779611 | 0.003544 | 0.027949 |\u001b[0m\n", "| 46 | -0.796838 | 224.809928 | 35.779836 | 35.779611 | 0.003544 | 0.027949 |\u001b[0m\n", "| 47 | -0.801462 | -225.851223 | 35.945565 | 35.945339 | 0.003549 | 0.027820 |\u001b[0m\n", "| 48 | -0.801462 | 225.851223 | 35.945565 | 35.945339 | 0.003549 | 0.027820 |\u001b[0m\n", "| 49 | -0.823218 | 257.880058 | 41.043095 | 41.042886 | 0.003192 | 0.024365 |\u001b[0m\n", "| 50 | -0.823218 | -257.880058 | 41.043095 | 41.042886 | 0.003192 | 0.024365 |\u001b[0m\n", "| 51 | -0.829849 | 232.223377 | 36.959734 | 36.959498 | 0.003573 | 0.027057 |\u001b[0m\n", "| 52 | -0.829849 | -232.223377 | 36.959734 | 36.959498 | 0.003573 | 0.027057 |\u001b[0m\n", "| 53 | -0.833132 | 232.986723 | 37.081226 | 37.080989 | 0.003576 | 0.026968 |\u001b[0m\n", "| 54 | -0.833132 | -232.986723 | 37.081226 | 37.080989 | 0.003576 | 0.026968 |\u001b[0m\n", "| 55 | -0.837692 | 252.830750 | 40.239484 | 40.239264 | 0.003313 | 0.024851 |\u001b[0m\n", "| 56 | -0.837692 | -252.830750 | 40.239484 | 40.239264 | 0.003313 | 0.024851 |\u001b[0m\n", "| 57 | -0.843057 | -274.636584 | 43.709976 | 43.709770 | 0.003070 | 0.022878 |\u001b[0m\n", "| 58 | -0.843057 | 274.636584 | 43.709976 | 43.709770 | 0.003070 | 0.022878 |\u001b[0m\n", "| 59 | -0.855990 | -264.468496 | 42.091689 | 42.091468 | 0.003237 | 0.023758 |\u001b[0m\n", "| 60 | -0.855990 | 264.468496 | 42.091689 | 42.091468 | 0.003237 | 0.023758 |\u001b[0m\n", "| 61 | -0.864725 | 271.184095 | 43.160509 | 43.160289 | 0.003189 | 0.023169 |\u001b[0m\n", "| 62 | -0.864725 | -271.184095 | 43.160509 | 43.160289 | 0.003189 | 0.023169 |\u001b[0m\n", "| 63 | -0.871325 | -283.421756 | 45.108187 | 45.107973 | 0.003074 | 0.022169 |\u001b[0m\n", "| 64 | -0.871325 | 283.421756 | 45.108187 | 45.107973 | 0.003074 | 0.022169 |\u001b[0m\n", "| 65 | -0.878445 | 267.336890 | 42.548217 | 42.547987 | 0.003286 | 0.023503 |\u001b[0m\n", "| 66 | -0.878445 | -267.336890 | 42.548217 | 42.547987 | 0.003286 | 0.023503 |\u001b[0m\n", "| 67 | -0.882869 | -280.833495 | 44.696260 | 44.696039 | 0.003144 | 0.022373 |\u001b[0m\n", "| 68 | -0.882869 | 280.833495 | 44.696260 | 44.696039 | 0.003144 | 0.022373 |\u001b[0m\n", "| 69 | -0.884024 | -245.027542 | 38.997598 | 38.997344 | 0.003608 | 0.025643 |\u001b[0m\n", "| 70 | -0.884024 | 245.027542 | 38.997598 | 38.997344 | 0.003608 | 0.025643 |\u001b[0m\n", "| 71 | -0.886589 | -245.661879 | 39.098557 | 39.098302 | 0.003609 | 0.025577 |\u001b[0m\n", "| 72 | -0.886589 | 245.661879 | 39.098557 | 39.098302 | 0.003609 | 0.025577 |\u001b[0m\n", "| 73 | -0.891210 | -288.915188 | 45.982499 | 45.982280 | 0.003085 | 0.021748 |\u001b[0m\n", "| 74 | -0.891210 | 288.915188 | 45.982499 | 45.982280 | 0.003085 | 0.021748 |\u001b[0m\n", "| 75 | -0.908699 | 251.206723 | 39.981053 | 39.980792 | 0.003617 | 0.025012 |\u001b[0m\n", "| 76 | -0.908699 | -251.206723 | 39.981053 | 39.980792 | 0.003617 | 0.025012 |\u001b[0m\n", "| 77 | -0.910251 | -251.606127 | 40.044621 | 40.044359 | 0.003618 | 0.024972 |\u001b[0m\n", "| 78 | -0.910251 | 251.606127 | 40.044621 | 40.044359 | 0.003618 | 0.024972 |\u001b[0m\n", "| 79 | -0.914184 | -241.156682 | 38.381554 | 38.381278 | 0.003791 | 0.026054 |\u001b[0m\n", "| 80 | -0.914184 | 241.156682 | 38.381554 | 38.381278 | 0.003791 | 0.026054 |\u001b[0m\n", "| 81 | -0.915395 | 290.517030 | 46.237451 | 46.237221 | 0.003151 | 0.021628 |\u001b[0m\n", "| 82 | -0.915395 | -290.517030 | 46.237451 | 46.237221 | 0.003151 | 0.021628 |\u001b[0m\n", "| 83 | -0.933974 | -278.955360 | 44.397373 | 44.397124 | 0.003348 | 0.022524 |\u001b[0m\n", "| 84 | -0.933974 | 278.955360 | 44.397373 | 44.397124 | 0.003348 | 0.022524 |\u001b[0m\n", "| 85 | -0.943144 | -260.320871 | 41.431625 | 41.431353 | 0.003623 | 0.024136 |\u001b[0m\n", "| 86 | -0.943144 | 260.320871 | 41.431625 | 41.431353 | 0.003623 | 0.024136 |\u001b[0m\n", "| 87 | -0.944542 | -260.700481 | 41.492043 | 41.491770 | 0.003623 | 0.024101 |\u001b[0m\n", "| 88 | -0.944542 | 260.700481 | 41.492043 | 41.491770 | 0.003623 | 0.024101 |\u001b[0m\n", "| 89 | -0.953003 | -294.814043 | 46.921357 | 46.921112 | 0.003233 | 0.021312 |\u001b[0m\n", "| 90 | -0.953003 | 294.814043 | 46.921357 | 46.921112 | 0.003233 | 0.021312 |\u001b[0m\n", "| 91 | -0.960626 | -295.652742 | 47.054844 | 47.054595 | 0.003249 | 0.021252 |\u001b[0m\n", "| 92 | -0.960626 | 295.652742 | 47.054844 | 47.054595 | 0.003249 | 0.021252 |\u001b[0m\n", "| 93 | -0.960976 | -265.315963 | 42.226624 | 42.226347 | 0.003622 | 0.023682 |\u001b[0m\n", "| 94 | -0.960976 | 265.315963 | 42.226624 | 42.226347 | 0.003622 | 0.023682 |\u001b[0m\n", "| 95 | -0.961739 | 300.017779 | 47.749558 | 47.749313 | 0.003206 | 0.020943 |\u001b[0m\n", "| 96 | -0.961739 | -300.017779 | 47.749558 | 47.749313 | 0.003206 | 0.020943 |\u001b[0m\n", "| 97 | -0.961940 | -265.596058 | 42.271203 | 42.270925 | 0.003622 | 0.023657 |\u001b[0m\n", "| 98 | -0.961940 | 265.596058 | 42.271203 | 42.270925 | 0.003622 | 0.023657 |\u001b[0m\n", "| 99 | -0.965385 | -266.582863 | 42.428259 | 42.427980 | 0.003621 | 0.023569 |\u001b[0m\n", "\u001b[36mFINISHED - Elapsed time = 16.3462206 seconds\u001b[0m\n", "\u001b[36mFINISHED - CPU process time = 68.0131146 seconds\u001b[0m\n" ] } ], "source": [ "data = sharpy.sharpy_main.main(['', ws.case_route + '/' + ws.case_name + '.sharpy'])" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Post-processing" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Nonlinear Equilibrium\n", "\n", "The files can be opened with Paraview to see the deformation and aerodynamic loading on the flying wing in trim conditions." ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "### Asymptotic Stability" ] }, { "cell_type": "code", "execution_count": 17, "metadata": {}, "outputs": [], "source": [ "eigenvalues_trim = np.loadtxt('./output/horten_u_inf2800_M4N11Msf5/stability/eigenvalues.dat')\n" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Flight Dynamics modes\n", "\n", "The flight dynamics modes can be found close to the origin of the Argand diagram. In particular, the phugoid is the mode that is closest to the imaginary axis. An exercise is left to the user to compare this phugoid predicition with the nonlinear response!" ] }, { "cell_type": "code", "execution_count": 18, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAY4AAAEKCAYAAAAFJbKyAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAeiklEQVR4nO3de5RcZZnv8e+PIAGSkCYgEU1DOJClcg/dXFzHS4JwDKOCzCCDV1BZ8XI41uosHHHhpQfPnEFkpu01XpYIHlGOEwQVoiKCmIQZ5ZaQyB0JSEwkTkDTSDqEW57zx97dXSn6sndSu2pX9++zVq3al7d2PfWmU0+9+937fRURmJmZZbVLswMwM7PW4sRhZma5OHGYmVkuThxmZpaLE4eZmeXixGFmZrns2uwAirbvvvvG7Nmzmx0G/f39TJkypdlhlILrYojrYojrYkgZ6mLlypVPRcQrh9s37hPH7NmzWbFiRbPDYNmyZcybN6/ZYZSC62KI62KI62JIGepC0tqR9vlUlZmZ5eLEYWZmuThxmJlZLk4cZmaWixOHmZnl4sRhZma5OHGYmVkuThxmZpaLE4eZmeXixGFmZrk4cZiZWS5OHGZmlosTh5mZ5eLEYWZmuThxmJlZLk4cZmaWS6kSh6QFkh6WtEbSBaOUO0NSSOpsZHxmZlaixCFpEvA14BTgUOA9kg4dptw04JPAHY2N0GznRMSo62atojSJAzgOWBMRj0XE88Bi4LRhyn0RuATY2sjgzHZGd3c3XV1dg8kiIujq6mLDhg1NjswsvzLNOf4aYF3V+nrg+OoCkuYC7RHxU0nnj3QgSQuBhQAzZ85k2bJl9Y82p82bN5cijjKYiHVx8MEHs3HjRq666ira29tZt24d7e3tTJ48ecLVxUgm4t/FSMpeF2VKHBpm22BbXtIuQA9wzlgHiojLgMsAOjs7o9mTvkM5Jp8vi4lYFwMtjPPPH/q9U6lU6OjomHB1MZKJ+HcxkrLXRZlOVa0H2qvWZwFPVK1PAw4Hlkl6HDgBWOIOcmsFkujp6dluW+26WasoU+K4C5gj6SBJuwFnAUsGdkbE0xGxb0TMjojZwO3AqRGxojnhmmU30OKoVrtu1ipKkzgi4kXgPOAXwIPADyLifkkXSTq1udGZ7biBpNHb20ulUmHbtm1UKhV6e3tZt26dr66yllOmPg4i4gbghpptnx+h7LxGxGS2syTR1tZGpVKhp6dnu9NWu+66K9Jw3Xtm5VWqxGE2XnV3dxMRg0liIHksX768yZGZ5VeaU1Vm411ty8ItDWtVThxmZpaLE4eZmeXixGFmZrk4cZiZWS5OHGZmlosTh1kTeIh1a2VOHGYNNtIQ693d3c0NzCwjJw6zBooI+vr66O3tHRyramA4kr6+Prc8rCX4znGzBqoebqS3t5f29vbBMawGhiMxKzu3OMwabKQh1p00rFWMmTgkzcjwaGtEsGbjwUhDrPs0lbWKLKeqnkgfo/0cmgQcUJeIzMax2iHWOzo6BodYB7c8rDVkSRwPRsTc0QpIWlWneMzGtdoh1pcvXz542qqtrc1Jw1pClsTxhjqVMTNGHmLdScNaxZh9HBGxFUDSuyVNS5c/J+lHko6pLmNm2XiIdWtlea6q+lxEPCPpjcD/AK4EvlFMWGZmVlZ5EsdL6fPbgW9ExPXAbvUPyczMyixP4vijpG8CZwI3SJqc8/VmZjYOZLmP4w1KTsCeCfwCWBARfcAM4FMFx2dmZiWTpcVwNrAS+DawF/AMQERsiIibCozNzMxKaMzLcSPiYwCSXgecAnxH0nRgKXAj8OuIeGmUQ5iZ2TiSuY8iIh6KiJ6IWACcCPwn8G7gjqKCMzOz8tmh0XEj4lnghvRhZmYTyJiJQ9Ki0fZHxL/WLxwzMyu7LC2Oaenza4FjgSXp+juBW4sIyszMyitL5/g/Aki6CTgmIp5J17uBawqNzszMSifPDXwHAM9XrT8PzK5rNGZmVnp5Ose/B9wp6cdAAKcD3y0kKjMzK63MiSMi/knSz4E3pZs+FBGeh8PMbILJeznu79PX7A5Mk/TmiHAHuZnZBJI5cUg6F6gAs4DVwAnAbSQ3A5qZ2QSRp3O8QnI57tqImA/MBZ4sJCozMyutPIlja9VsgJMj4iGSezvqRtICSQ9LWiPpgmH2L5L0gKR7JN0i6cB6vr+ZmY0tT+JYL6kNuA64WdL1wBP1CkTSJOBrJAMpHgq8R9KhNcVWAZ0RcSRwLXBJvd7fzMyyydTHkc7H8cl0Ho5uSUuB6SSj49bLccCaiHgsfc/FwGnAAwMFImJpVfnbgffX8f3NzCyDTIkjIkLSdUBHur68gFheA6yrWl8PHD9K+Y8APy8gDjMzG0Wey3Fvl3RsRNxVUCwaZlsMW1B6P9AJvGWE/QuBhQAzZ85k2bJldQpxx23evLkUcZSB62KI62KI62JI2esiT+KYD3xU0lqgn+SLPtL+hnpYD7RXrc9imD4USScBFwJviYjnhjtQRFwGXAbQ2dkZ8+bNq1OIO27ZsmWUIY4ycF0McV0McV0MKXtd5EkcpxQWReIuYI6kg4A/AmcB760uIGku8E2Sec83FhyPmZkNI8+QI2uLDCQiXpR0HvALYBLw7Yi4X9JFwIqIWAJ8GZgKXJP01/OHiDi1yLjMzGx7WSZyujsijtnZMllExMtmFYyIz1ctn7Sz72FmZjsnS4vj9ZLuGWW/SC7NNTOzCSBL4nhdhjIv7WwgZmbWGrLMAFho34aZmbWWPEOOmJmZOXGYmVk+uROHpCnpgIRmZjYBjZk4JO0i6b2SfiZpI/AQsEHS/ZK+LGlO8WGamVlZZGlxLAUOBj4DvCoi2iNiP5K5x28HLk7HjjIzswkgy+W4J0XEC7UbI+IvwA+BH0p6Rd0jMzOzUhqzxTFc0tiRMmZmNj7kGeQQSe3AYcDhwBHAYRHRWURgZmZWTlk6xz8q6TeS+oDfAeeSDDS4hJrRa83MbPzL0uL4DPD3wFPAxcAeJCPX/qHIwMzMrJyyXFX1joi4IyIejYh3A18FfiKpS5JvIDQzm2CydI7fV7N+I3AcMAP4dUFxmZlZSWXp43jZXOAR8VxEfA44e6QyZmY2PmW6AVDS/5J0QPVGSbsBsyRdSZpAzMxs/MvSOb4A+DDw75L+G7CJpIN8F+AmoCciVhcXopmZlUmW+Ti2Al8Hvp62MvYFtkREX9HBmZlZ+WS+AVDSt4C/BbYAT6TTyd4TEf9WVHBmZlY+ee4cfzPJIIcvSHoNcBRwZDFhmZlZWeVJHLcDewMbI+KPwB+BGwqJyszMSivPDXyXAcslnS/pTZKmFxWUmZmVV57EcRXwA5JWyieA30h6tJCozMystPKcqlofEV+o3iBpcp3jMTOzksvT4lgtqVK9ISKeq3M8ZmZWcnlaHDOBkyR9Grgb+C2wOiKuKSQyMzMrpcyJIyLOhMHTU4eRTOR0PODEYWY2gYyZOCR9HbgXuAe4NyL+StLiuLvg2MzMrISytDhWk9zodxZwuKRn2D6RLC4wPjMzK5ksY1VdVr0uaRZJIjkCeDvgxGFmNoHk6RwHICLWA+vxXeNmZhOSp341M7NcnDjMzCwXJw6zAkTEqOtbtmwZdd2szHY4cUjav95DjkhaIOlhSWskXTDM/smSrk733yFpdj3f36weuru76erqGkwWEUFXVxfd3d0ASGLKlCmDyWLLli1MmTIFSc0K2SyXnWlxfA94SNKl9QhE0iTga8ApwKHAeyQdWlPsI8CmiDgE6AG+VI/3NquXiKCvr4/e3t7B5NHV1UVvby99fX309/cPlp0yZQrbtm1jypQpg9vc8rBWkPuqqgERcZKSn0i1X+476jhgTUQ8BiBpMXAa8EBVmdOA7nT5WuCrkhS15wHMmkQSPT09APT29tLb2wtApVKhp6cHSfT39w8mi1WrVg2+tr+/nz333LPxQZvlpKzfuZK+FBGfHmvbDgcinQEsiIhz0/UPAMdHxHlVZe5Ly6xP1x9NyzxVc6yFwEKAmTNndixe3PxbTTZv3szUqVObHUYpTJS6WLly5eByR0fHdvu2bdvGqlWrmDVrFuvXr2fu3LnsssvE7nKcKH8XWZShLubPn78yIjqH3RkRmR7A3cNsuyfr6zMc/93A5VXrHwD+rabM/cCsqvVHgX1GO25HR0eUwdKlS5sdQmmM97rYtm1bVCqVAAYflUoltm3bFhER/f39g9svvfTSweX+/v4mR95c4/3vIo8y1AWwIkb4Xh3zJ46kj0u6F3idpHuqHr8nGXakXtYD7VXrs4AnRiojaVdgOvCXOsZgtlOiqk+jUqmwbds2KpXKYJ9H9WkqgLlz5w4uV3eYm5VZlj6O7wO/AC4HPlS1/ZmIqOeX9l3AHEkHkcxnfhbw3poyS4CzgduAM4BfpZnRrBQk0dbWtl2fxkCfR1tb23ZJo7+/nzvvvHO7ZOI+DmsFWcaqehp4WlJbRKwtKpCIeFHSeSRJahLw7Yi4X9JFJE2mJcAVwPckrSFpaZxVVDxmO6q7u5uIGLy8diB5DKxHBFu2bBlMEnvuuac7xq2l5Lmq6jZJx0bEXUUFExE3UDMGVkR8vmp5K0lfiFmp1d6TUbtemyScNKyV5Ekc84GPSloL9AMCIiKOLCQyMzMrpTyJ45TCojAzs5aRZ+rYtZL2BuYAu1ftKqzfw8zMyidz4pB0LlAhuUx2NXACydVNJxYTmpmZlVGeW1UrwLHA2oiYD8wFniwkKjMzK608iWNrelUTkiZHxEPAa4sJy8zMyipP5/h6SW3AdcDNkjbx8ju7zcxsnMvTOX56utgtaSnJcB83FhKVmZmV1piJQ9LuwMeAQ4B7gSsiYnnRgZmZWTll6eO4EugkSRqnAP9SaERmZlZqWU5VHRoRRwBIugK4s9iQzMyszLK0OF4YWIiIFwuMxczMWkCWFsdRkv6aLgvYI10fGKtqr8KiMzOz0skyrPqkRgRiZmatYWJPcmxmZrk5cZiZWS5OHGZmlkvmxCHpvHRYdTMzm8DytDheBdwl6QeSFqh2LkwzM5sQMieOiPgsySROVwDnAI9I+j+SDi4oNjMzK6FcfRwREcCf0seLwN7AtZIuKSA2MzMroTwzAH4SOBt4Crgc+FREvCBpF+AR4B+KCdHMzMokU+JI+zOOAv42IrabYzwitkl6RxHBmZlZ+WQ6VZWeoppbmzSq9j9Y16jMzKy08vRx3Cbp2MIiMTOzlpBn6tj5wEclrQX6GRrk8MhCIjMzs1LKkzhOKSwKMzNrGXnmHF+b3jk+B9i9atew/R5mZjY+5bkc91ygAswCVgMnALcBJxYTmpmZlVGezvEKcCywNiLmA3OBJwuJyszMSitP4tgaEVsBJE2OiIeA1xYTlpmZlVWezvH1ktqA64CbJW0CnigmLDMzK6s8neOnp4vdkpYC04EbC4nKzMxKK0+LY1BELK93IGZm1hryXFU1Gfg7YHb16yLiop0NQtIM4Or02I8DZ0bEppoyRwPfAPYCXgL+KSKu3tn3NjOzfPJ0jl8PnEYynHp/1aMeLgBuiYg5wC3peq0twAcj4jBgAfCVtM/FzMwaKM+pqlkRsaCgOE4D5qXLVwLLgE9XF4iI31UtPyFpI/BKoK+gmMzMbBh5Why/kXREQXHMjIgNAOnzfqMVlnQcsBvwaEHxmJnZCJSMmJ6hoPQAyXAjjwHPkXOQQ0m/JJm3vNaFwJUR0VZVdlNE7D3CcfYnaZGcHRG3j1BmIbAQYObMmR2LFy/OEmKhNm/ezNSpU5sdRim4Loa4Loa4LoaUoS7mz5+/MiI6h9uXJ3EcQJosqrdHxB92NkBJDwPzImLDQGKIiJfdXChpL5Kk8c8RcU2WY3d2dsaKFSt2NsSdtmzZMubNm9fsMErBdTHEdTHEdTGkDHUhacTEMWYfh6T/jIg3AvezfdIYSCJ71SHGJSTT0l6cPl8/TBy7AT8Gvps1aZiZWf2N2ceRJg0iYlpE7FX1mBYR9UgakCSMkyU9ApycriOpU9LlaZkzgTcD50hanT6OrtP7m5lZRjt0A2C9RcSfgbcOs30FcG66fBVwVYNDMzOzGnluAFw0zOangZURsbp+IZmZWZnluRy3E/gY8Jr0sZDk3otvSfqH+odmZmZllOdU1T7AMRGxGUDSF4BrSfodVgKX1D88MzMrmzwtjgOA56vWXwAOjIhnSe7rMDOzCSBPi+P7wO2SBi6VfSfw75KmAA/UPTIzMyulPPNxfFHSDcAbSe7h+Fh61RPA+4oIzszMyifv5biPAZOA3YE9Jb05Im6tf1hmZlZWeS7HPReoALOA1cAJwG3AicWEZmZmZZSnc7wCHAusjYj5wFzgyUKiMjOz0sqTOLZGxFZIZgOMiIeAlw1EaGZm41uePo716Yx71wE3S9oEPFFMWGZmVlZ5rqo6PV3slrQUmA7cWEhUZmZWWjs0yGFELK93IGZm1hryXFXVSTJb34HVr8s6A6CZmY0PeVoc/w/4FHAvsK2YcMzMrOzyJI4nI2JJYZGYmVlLyJM4vpDOxncLVYMaRsSP6h6VmZmVVp7E8SHgdcArGDpVFYATh1lOEYGkEdfNyixP4jgqIo4oLBKzCaK7u5u+vj56enqAJGl0dXXR1tZGd3d3c4MzyyDPneO3Szq0sEjMJoCIoK+vj97eXrq6ugDo6uqit7eXvr4+IqLJEZqNLU+L443AOZIeI+njEBC+HNcsO0mDLY3e3l7a29vp7e2lUqnQ09Pj01XWEvK0ON4GHAKcDLwDeHv6bGY5VCePAU4a1krGTBySnpH0V+A+kns47ksf96fPZpbDQJ9Gta6uLp+mspYx5qmqiJjWiEDMJoKBpDFweqqjo4NKpUJvby/gloe1hh0aq8rMdowk2traBvs0li9fPnjaqq2tzUnDWoITh1mDdXd3b3ffxkCfh5OGtYo8neNmVie1ScJJw1qJE4eZmeXixGFmZrk4cZiZWS5OHGZmlosTh5mZ5eLEYdYgtXeG+05xa1VOHGYN0N3dvd2wIgN3kG/YsKHJkZnl58RhVrDaodSrhx158cUX3fKwllOKO8clzQCuBmYDjwNnRsSmEcruBTwI/DgizmtUjGY7qnYo9YFxqSqVCu3t7b75z1pOWVocFwC3RMQckjnNLxil7BeB5Q2JyqxORhpK3awVlSVxnAZcmS5fCbxruEKSOoCZwE0NisusLkYaSt2sFakM51cl9UVEW9X6pojYu6bMLsCvgA8AbwU6RzpVJWkhsBBg5syZHYsXLy4s9qw2b97M1KlTmx1GKUzEuli3bh0bN25kv/32o729fXD9oIMOYsaMGc0OrxQm4t/FSMpQF/Pnz18ZEZ3D7WtYH4ekXwKvGmbXhRkP8QnghohYN9Y54Yi4DLgMoLOzM+bNm5cj0mIsW7aMMsRRBhOxLrq7u+nr62PRokVIGmyBvPrVr55wdTGSifh3MZKy10XDEkdEnDTSPkn/JWn/iNggaX9g4zDF3gC8SdIngKnAbpI2R8Ro/SFmpTDSUOrLl7u7zlpPKa6qApYAZwMXp8/X1xaIiPcNLEs6h+RUlZOGtQwPpW7jRVk6xy8GTpb0CHByuo6kTkmXNzUyMzPbTilaHBHxZ5IO79rtK4Bzh9n+HeA7hQdmZmYvU5YWh5mZtQgnDjMzy8WJw8zMcnHiMDOzXJw4zMwsFycOMzPLxYnDzMxyceIwM7NcnDjMzCwXJw4zM8vFicPMzHJx4jAzs1ycOMzMLBcnDjMzy8WJw8zMcnHiMDOzXBQRzY6hUJKeBNY2Ow5gX+CpZgdREq6LIa6LIa6LIWWoiwMj4pXD7Rj3iaMsJK2IiM5mx1EGroshroshroshZa8Ln6oyM7NcnDjMzCwXJ47GuazZAZSI62KI62KI62JIqevCfRxmZpaLWxxmZpaLE4eZmeXixFEQSTMk3SzpkfR571HK7iXpj5K+2sgYGyVLXUg6WtJtku6XdI+kv29GrEWRtEDSw5LWSLpgmP2TJV2d7r9D0uzGR1m8DPWwSNID6d/ALZIObEacjTJWfVSVO0NSSCrFJbpOHMW5ALglIuYAt6TrI/kisLwhUTVHlrrYAnwwIg4DFgBfkdTWwBgLI2kS8DXgFOBQ4D2SDq0p9hFgU0QcAvQAX2pslMXLWA+rgM6IOBK4FriksVE2Tsb6QNI04JPAHY2NcGROHMU5DbgyXb4SeNdwhSR1ADOBmxoUVzOMWRcR8buIeCRdfgLYCAx712oLOg5YExGPRcTzwGKSOqlWXUfXAm+VpAbG2Ahj1kNELI2ILenq7cCsBsfYSFn+LiD5YXkJsLWRwY3GiaM4MyNiA0D6vF9tAUm7AP8CfKrBsTXamHVRTdJxwG7Aow2IrRFeA6yrWl+fbhu2TES8CDwN7NOQ6BonSz1U+wjw80Ijaq4x60PSXKA9In7ayMDGsmuzA2hlkn4JvGqYXRdmPMQngBsiYl2r/7isQ10MHGd/4HvA2RGxrR6xlcBw/7i118FnKdPqMn9GSe8HOoG3FBpRc41aH+kPyx7gnEYFlJUTx06IiJNG2ifpvyTtHxEb0i/DjcMUewPwJkmfAKYCu0naHBGj9YeUUh3qAkl7AT8DPhsRtxcUajOsB9qr1mcBT4xQZr2kXYHpwF8aE17DZKkHJJ1E8oPjLRHxXINia4ax6mMacDiwLP1h+SpgiaRTI2JFw6Ichk9VFWcJcHa6fDZwfW2BiHhfRBwQEbOB84HvtmLSyGDMupC0G/Bjkjq4poGxNcJdwBxJB6Wf8yySOqlWXUdnAL+K8Xd37pj1kJ6a+SZwakQM+wNjHBm1PiLi6YjYNyJmp98Rt5PUS1OTBjhxFOli4GRJjwAnp+tI6pR0eVMja7wsdXEm8GbgHEmr08fRzQm3vtI+i/OAXwAPAj+IiPslXSTp1LTYFcA+ktYAixj9KryWlLEevkzS+r4m/RuoTbDjRsb6KCUPOWJmZrm4xWFmZrk4cZiZWS5OHGZmlosTh5mZ5eLEYWZmuThxmJlZLk4cVmqSXkqv579P0k92ZsRcSZszvMc1kvbMccy29M7/LGU/KulPkn4r6VFJH8zwmj0kLU9HUkXSEZLWSvp4VZndJN2a3nFe+/rZkp6VtDrrZxohjm5J59ds+6ak/z5CzKslPS9p3515XysnJw4ru2cj4uiIOJxkCI7/WfB7PA98LMuL0tFrZ5CMOZbFkUB3RBwFvAf41wyv+TDwo4h4CSAi7iW5w3gw6aQjq94CjDSHyaMR8bKbKZXYme+A40nuZt5ORDybvt/LhhOx8cGJw1rJbaSjh0p6v6Q701+23xz4RZ7uu07SSiWTQi3M+R7/ARwy0nHSX/APSvo6cDfJHd8Hp3F8eYxjHwE8nC7/niRJDcR8kKTrJa1IP9dr013v4+VDtGwEDqvZdl1adlTDxN8+Un1JulDJJEO/BF5bc5zXA78Ddpf0s7QVdZ/G2QRcNoKI8MOP0j6AzenzJOAakkmeXg/8BHhFuu/rJJNADbxmRvq8B3AfsE/1sUZ5j11JvqQ/PtJxgNnANuCEdN9s4L6Mn2UT8GqSUVH/EfhQuv0VJC2Gg9P1vwH+L8nQ8n8a5jjXAM8BB1ZtmwQ8OUzZ7eKrjX+Uz9kB3AvsCewFrAHOr3rNIpLW0N8B36raPr1q+XFg32b/DflR/4dHx7Wy2yM9Pz8bWAncDHyc5IvtrnTU0D3YfsTdT0o6PV1uB+YAf87wHpC0OK4Y5Th/AtZGztF7JbWTjHZ6A0mr6R6gO939LpIWxA/Tz7NrGse+QF/NcRYAU0hGET4MWAsQES+lfQrTIuKZMcKpjX+4z3kC8ONIJ1UaZsyotwEfIhlX6lJJXwJ+GhH/McZ72zjgxGFl92xEHC1pOvBTkj6OAK6MiM/UFpY0DzgJeENEbJG0DNg9y3vkOE7/DnyOI4FbI+JEJXOu30cyrP5vgKOACyPiiuoXpOV2r1rfnWQmuFNJvrQPJ0lEAyaTbZa4wfjH+JwjzZWxJ9AWyUyNA7NY/g3wz5JuioiLMsRgLcx9HNYSIuJpknmXzwduBc6QtB+ApBmSDkyLTieZu3uLpNeR/HLeEVmP8wxJS2KQpFsk1c5sdwTJfNpExCbg+8Db030bgLcNdFSnV04pLTcpTRgAnyUZdv5xktNIh1e95z4kp6peqNPnvBU4Pb1CahrwzqrXzAeWpu/7amBLRFwFXAock/P9rQU5cVjLiIhVwG9Jfr1/FrhJ0j0kp6/2T4vdCOyabv8iw1z1k1Gm40TEn4Ffpx3DX06//A/h5ZMwDSaO1E9IfqUDfJvk/+KD6SmzT0fEwK/9m4A3pp3lJwNfSbdvlzhIvsyrWx9ZDfs5I+Ju4GpgNfBDklNnA05JXzfwue5M474Q+N87EIO1GA+rblZHkg4HPhwRi+p0vLnAooj4wBjlfgR8JiIertk+m6Tv4fDhXreDMd0NHD9W60bS40BnRDxVr/e2cnCLw6yOIuK+eiWN9HirgKXVlxvXUjJ73HW1SSP1EjB9Z28ArInpmNGSxsANgCRXi42XeeOtilscZmaWi1scZmaWixOHmZnl4sRhZma5OHGYmVkuThxmZpaLE4eZmeXixGFmZrk4cZiZWS7/H/gm2+tWJNTyAAAAAElFTkSuQmCC", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig = plt.figure()\n", "plt.scatter(eigenvalues_trim[:, 0], eigenvalues_trim[:, 1],\n", " marker='x',\n", " color='k')\n", "plt.xlim(-0.5, 0.5)\n", "plt.ylim(-0.5, 0.5)\n", "plt.grid()\n", "plt.xlabel('Real Part, $Re(\\lambda)$ [rad/s]')\n", "plt.ylabel('Imaginary Part, $Im(\\lambda)$ [rad/s]');" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "#### Structural Modes\n", "\n", "Looking further out on the plot, the structural modes appear. There is a curve given the Newmark-$\\beta$ integration scheme and on top of it several modes are damped by the presence of the aerodynamics.\n", "\n", "Try changing `newmark_damp` in the `LinearAssembler` settings to see how this plot changes!" ] }, { "cell_type": "code", "execution_count": 19, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZEAAAEOCAYAAABIESrBAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjMsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy+AADFEAAAgAElEQVR4nO3de5zcdX3v8dc7iwQKIRuaIEiCQRtSuRkIEHzUaqJYghcQWywcqxbriSh0N1npEUo5mdRDj1W6yW4Vj1Go2FZpLAKxUC7yYKEXEiC4hCCg4RJZAQHNhkS5mOzn/PH7zfLbzezuzGZm5/Z+Ph7z2PldZubzZcN+5ntXRGBmZjYek6odgJmZ1S8nETMzGzcnETMzGzcnETMzGzcnETMzGzcnETMzG7eaSSKSZkm6Q9LDkh6S1J6eP1DSbZJ+kv6clp6XpG5JmyVtlHR8dUtgZtZ8aiaJADuBz0bEW4CTgfMlHQlcBNweEXOA29NjgNOAOeljCfDViQ/ZzKy51UwSiYhnIuL+9Pl24GHgUOAM4Or0tquBD6bPzwC+FYl1QKukQyY4bDOzplYzSSRL0mzgOGA98PqIeAaSRAMclN52KPBU5mV96TkzM5sge1U7gOEk7Q9cCyyNiBcljXhrgXMF13CRtISkyYt99tln/mGHHVaOUOvCwMAAkybV5HeFinGZm0O5yvzqq6+yc+fOweO99tqLvffee4/ftxKq9Xv+8Y9//EJEzCh4MSJq5gG8DrgF6MicexQ4JH1+CPBo+vxrwDmF7hvtccQRR0QzueOOO6odwoRzmZtDOco8MDAQbW1tQfIFNIBoa2uLgYGBPQ+wAqr1ewbuixH+ptbMVxclVY4rgYcjojNzaS3w8fT5x4EbMuc/lo7SOhnYFmmzl5nZWCKCpUuX0t3dPeR8d3c3S5cuzX85tTHUTBIBfg/4KPAuSb3p473AF4D3SPoJ8J70GOAm4HFgM/B14DNViNnM6tj69esBaGtrY2BggLa2tiHnbWw10ycSEf9J4X4OgHcXuD+A8ysalJk1LEksXryYBQsWsGrVKiSxatUqAKZNm8Yo/bGWUTNJxMxsouVyOSJiMGHkE4kTSPFqqTnLzGzCDU8YTiClcRIxM7NxcxIxMytg+Ogsj9YqzEnEzGyYXC7HsmXLGBgYAJIEsmzZMpYvX17lyGqPk4iZWUZE0N/fT1dXF/Pnz2dgYIBly5bR1dXF2rVrnUiG8egsM7MMSXR2dnLnnXfS29tLS0sLAPPmzaO3t5d3vvOdQ0Z0NTvXRMzMhpk0aRIbNmwYcq63t5f29nZWrlzpBJLhJGJmNkxE0NHRsdv5zs5OJ5BhnETMzDLynehdXV3MmzdvyLV8H0n+PnMSMTMbQhJTp04d7ANpb2/n0ksvZfr06fT29tLR0THY2Z7L5aodbtU5iZiZDbNixQpOP/102tvb6ezs5MUXX+SFF15g3rx5TJ06lY6ODrq6uti6dWvT10g8OsvMrIAVK1YMjsJauXIl69atY/369fT29gIMrvi7YsWKpq6RuCZiZjaCbCf6ggULdrve3d1Nf39/U9dGnETMzMahu7ubtra2ph/y6yRiZjaK/GitfNKwoZxEzMxGIYnW1taCCcQ7IDqJmJmNKb9eVnd3N+3t7YNb6a5fv55ly5Y1dZ9ITY3OknQV8H7guYg4Oj2XA/4n8Hx6219GxE3ptYuBPwN2AW0RccuEB21mDU8S06ZNG7LsSX4HxNbW1qbuE6mpJAJ8E/gy8K1h51dGxOXZE5KOBM4GjgLeAPxA0hERsWsiAjWz5lJoK91m71SHGmvOioi7gF8WefsZwDUR8UpEPAFsBk6qWHBm1vS8le7uaiqJjOICSRslXSVpWnruUOCpzD196TkzM5sgtdacVchXgc8Dkf78O+ATQKGvAAV7tyQtAZYAzJgxg56enooEWot27NjRVOUFl7lZuMy1oeaTSET8PP9c0teBf0sP+4BZmVtnAk+P8B6rgdUAc+fOjYULF1Yk1lrU09NDM5UXXOZm4TLXhppvzpJ0SObwTGBT+nwtcLakyZIOB+YA90x0fGZmsPvS8M0y7LemaiKSvgMsBKZL6gOWAwslzSNpqnoS+BRARDwkaQ3wI2AncL5HZplZNeRyOfr7+wdHa+Vnube2tjb84ow1lUQi4pwCp68c5f7LgMsqF5GZ2egigv7+frq6ugBYuXLl4KZW7e3tDb8fe00lETOzepOfLwLQ1dU1mEzyExMbXc33iZiZ1bpsIsnLHzf6DoiuiZiZ7aGIYOnSpUPO5Y/z6201arOWk4iZ2R4otFR8d3c33d3dAA2/54iTiJnZHsgvwpjtA8knEGBwocZG5T4RM7M9lMvlhvSBZDX6UvFOImZmZZId2jswMEB7eztdXV0NnUjGbM6SdGAR7zMQEf1liMfMrC4Nb9aSRGdnJ8DgniON2LleTJ/I0+ljtJK3AIeVJSIzszqV3XMkP4u9s7OTSZMmNews9mKasx6OiDdFxOEjPYBfVDpQM7N6kK9x5Gexd3R0DCaQrq4u+vv7G6ppq5iayNvKdI+ZWVMYaxZ7IzVpjVkTiYiXASSdJWlK+vxSSd+TdHz2HjMzS4w0i72REgiUNjrr0ojYLuntwB8AV5NsGGVmZsPkm7CyGnGUVilJJL/M+vuAr0bEDcDe5Q/JzKy+ZftA8sN929rahgz3bZRkUkoS+ZmkrwEfBm6SNLnE15uZNYXhw31XrFgBJEugtLa2Ao2zMOOYSUDS25Q04n0YuAVYnM4JORD4iwrHZ2ZWl7Kz2Pv7+weXQlm+fHlDjdQqZnTWx4GvAD8Gbga2A0TEM8AzlQvNzKy+5TvRsyO18smkUUZqFTM667yIOB7IAdOAb0q6W9LfSHqHpJZKB2lmVs8aeaRW0X0aEfFIRKyMiMXAu4D/BM4C1lcqODOzRtDII7XG1TEeES9FxE0R8ecRcUK5gpF0laTnJG3KnDtQ0m2SfpL+nJael6RuSZslbczPWTEzqyWFRmo10sKMxSzA2DHa9YjoLF84fBP4MvCtzLmLgNsj4guSLkqPPwecBsxJHwtI5qwsKGMsZmZ7rNDCjPmmrfzCjPWsmI71KenPucCJwNr0+APAXeUMJiLukjR72OkzgIXp86uBHpIkcgbwrUjS+DpJrZIOSTv8zcxqRnZhRnitj6TeEwgUkUQiYgWApFuB4yNie3qcA75b0egSr88nhoh4RtJB6flDgacy9/Wl53ZLIpKWAEsAZsyYQU9PT0UDriU7duxoqvKCy9wsXObaUMr2uIcBr2aOXwVmlzWa0hRK4QUbFyNiNbAaYO7cubFw4cIKhlVbenp6aKbygsvcLFzm2lBKEvlH4B5J15H8sT6ToX0XlfLzfDOVpEOA59LzfcCszH0zSfY9MTOzCVLKEN/LgHOBrUA/cG5E/E2lAstYSzLhkfTnDZnzH0tHaZ0MbHN/iJnZxCqlJgLwRPqafYApkt4REWXrXJf0HZJO9OmS+oDlwBeANZL+DPgpydwUgJuA9wKbgV+TJDgzM5tARScRSZ8E2kmajXqBk4G7SSYelkVEnDPCpXcXuDeA88v12WZm1TB83/V624e9lMmG7SRDfLdExCLgOOD5ikRlZtYEcrnckAmH+YmJ9bS6bylJ5OXMLoeTI+IRkrkjZmZWouw+7PlEUo+r+5bSJ9InqRW4HrhN0lY8GsrMbFwaZR/2omoi6X4ibRHRHxE54FLgSuCDFYzNzKyhNcLqvkUlkbQT+/rM8Z0RsTYiXh3lZWZmNopGWN23lD6RdZJOrFgkZmZNpFFW9y2lT2QR8ClJW4BfkSw7EhFxbEUiMzNrYI2yum8pSeS0ikVhZtaEhq/uC0P7ROphzkjRSSQitlQyEDOzZpRPErlcjv7+/sHaSL65q7W1tabnjYzZJyLp/nLcY2ZmhdXznJFiaiJvkbRxlOsCppYpHjOzplPPc0aKSSK/W8Q9u/Y0EDOzZpZPJPkEAvUxZ2TM5qyI2FLEo28igjUza1T1OmeklHkiZmZWAfU8Z6TU/UTMzKzM6nnOSMlJRNJ+JCv6uh/EzKxMsnNG8j/zCaWW54uMmUQkTQLOBj5Csp/IK8BkSc+T7C64OiJ+UtEozcyagKQh80XyCSQ/X2ThwoXVDnE3xfSJ3AG8GbgYODgiZkXEQcDvA+uAL0j6kwrGaGbWFMaaL1KLimnOOiUifjP8ZET8ErgWuFbS68oe2TCSngS2kwwn3hkRJ0g6EPgXYDbwJPDhiNha6VjMzCphrPkid955ZzXDK6iYIb67JZDx3FMmiyJiXkSckB5fBNweEXOA29NjM7O6VW97jJQ0xFfSLEmLJV0o6WpJ91UqsCKdAVydPr8ab5JlZnWu3uaLaKzAJH0K+DhwJDAZuBHYBDwIPBgRP650kGkcTwBbgQC+FhGrJfVHRGvmnq0RMa3Aa5cASwBmzJgxf82aNRMRck3YsWMH+++/f7XDmFAuc3No1DI/9dRTPPfccxx00EHMmjVryPG0adOqUuZFixZtyLQADVFMEnkS+GPgBeALwL7AZyLip2WOc6w43hART0s6CLgN+HNgbTFJJGvu3Lnx6KOPVjja2tHT01OTIzoqyWVuDo1a5rFGZ1WjzJJGTCLFdKy/PyI2pc/PkrQY+L6kbwJdETFQpjhHFRFPpz+fk3QdcBLwc0mHRMQzkg4BnpuIWMzMKmX4HiPZ+SI9PT3VDa6AYjrWNw07vpnkD/iBwH9VKK4hJO0naUr+OfAHJE1qa0ma2kh/3jAR8ZiZVdLwTvRa7VSH4iYbKoa1eUXEK8Clkv5xpHvK7PXAdel/yL2Ab0fEzZLuBdZI+jPgp8BZFYzBzMyGKaY56w5J1wI3ZPtBJO0NzJR0CcmExG9WJkSIiMeBtxY4/wvg3ZX6XDMzG10xSWQx8AngO5LeRDJCal+SprBbgZUR0Vu5EM3MrFaNmUQi4mXgCuCKtPYxHfh1RNTmHHwzM5swRa/iK+nrwIeAXwNPp1vmboyIv69UcGZmVttKWQr+HSQLMP5G0qEkfRTHViYsMzOrB6UkkXXANOC5iPgZ8DOSpeDNzKxJlbJ21mrgznTdrN+XNLVSQZmZWX0oJYn8E7CGpPbyGeC/JT1WkajMzKwulNKc1RcRy7MnJE0uczxmZlZHSqmJ9Epqz55IZ66bmVmTKqUm8nrgFEmfA+4HHgB6I+K7FYnMzMxqXtFJJCI+DINNWEcBxwALACcRM7MmVcwCjFeQbEC1kWQTqhdJaiL3Vzg2MzOrccXURHpJJhWeDRwtaTtDk8o1FYzPzMxqWDFrZ63OHkuaSZJUjgHeBziJmJk1qVI61gGIiD6gD89WNzNreqUM8TUzMxvCScTMzMbNScTMrIbt2rWLHTt2jHhcbeNOIpIOqYVlTyQtlvSopM2SLqp2PGZm5TJ79mz22msvpkyZMpg4Dj74YKZMmYKkKkeX2JOayD8Cj0i6vFzBlEpSC/AV4DTgSOAcSUdWKx4zs3LZtWsXv/rVrwaPp0yZwgMPPMALL7wweK4WaiTjTiIRcQrwJuAfyhdOyU4CNkfE4xHxKslw4zOqGI+ZWVm0tLTw7LPPMn369MFzO3fuHHy+fft29t9//2qENoQiorgbpb+NiM+NdW4iSfojYHFEfDI9/iiwICIuGHbfEmAJwIwZM+avWbNmwmOtlh07dtTEP7SJ5DI3h2Yq84YNGwCYOXMmfX19HHfccUyaNHFd2osWLdoQEScUulbKPJH3AMMTxmkFzk2kQo2Cu2XFdMLkaoC5c+fGwoULKxxW7ejp6aGZygsuc7NohjLv2rWLgw8+eLAJ6/LLL+fCCy8EaqcmUszaWZ8m2YTqzZI2Zi5NAf6rUoEVqQ+YlTmeCTxdpVjMzMpmeAIB2Guv1/5kT5kypSYSSTE1kW8DtwDfAM7NnN8eEb+sSFTFuxeYI+lwkj3fzwb+R3VDMjPbcy0tLey3336DSWT79u3cd999TJ8+ffBctRMIFLd21jZgm6TWiNgyATEVLSJ2SrqAJMm1AFdFxENVDsvMrCyefPJJdu3axUsvvTSYMJ599tkhx9VWSp/I3ZJOjIh7KxbNOETETXgdLzNrUC0tLUMSxvDjaisliSwCPiVpC/Arkk7tiIhjKxKZmZnVvFKSyGkVi8LMzOpSKdvjbpE0DZgD7JO5VFP9JGZmNnGKTiKSPgm0kwyj7QVOBu4G3lWZ0MzMrNaVMuWxHTgR2BIRi4DjgOcrEpWZmdWFUpLIyxHxMoCkyRHxCDC3MmGZmVk9KKVjvU9SK3A9cJukrXh2uJlZUyulY/3M9GlO0h3AVODmikRlZmZ1oZi1s/YBzgN+B3gQuDIi7qx0YGZmVvuK6RO5GjiBJIGcBvxdRSMyM7O6UUxz1pERcQyApCuBeyobkpmZ1YtiaiK/yT+JiJ2j3WhmZs2lmJrIWyW9mD4XsG96nF8764CKRWdmZjWtmKXgWyYiEDMzqz8Tt0mvmZk1HCcRMzMbNycRMzMbt6KTiKQL0qXgzczMgNJqIgcD90paI2mxJFUqqCxJOUk/k9SbPt6buXaxpM2SHpV06kTEY2Zmryk6iUTEX5FsSHUl8KfATyT9jaQ3Vyi2rJURMS993AQg6UjgbOAoYDFwhSSPJDOzuhcRox7XkpL6RCIpybPpYycwDfhXSV+sQGxjOQO4JiJeiYgngM3ASVWIw8ysbHK5HMuWLRtMHBHBsmXLyOVy1Q1sBCo2w0lqAz4OvAB8A7g+In4jaRLwk4ioSI1EUo6k5vMicB/w2YjYKunLwLqI+Kf0viuBf4+Ify3wHkuAJQAzZsyYv2bNmkqEWpN27NjB/vvvX+0wJpTL3BwatcxPPfUUzz33HAcddBCzZs0acjxt2rSqlHnRokUbIuKEghcjYswHyez0K4E3jnD9LcW8zyjv/wNgU4HHGcDrgRaSWtNlwFXpa74C/EnmPa4E/nCszzriiCOimdxxxx3VDmHCuczNoVHLPDAwEO3t7QEMPtrb22NgYKBqZQbuixH+phbVnJW+yXERsWWE6w8X8z6jvP8pEXF0gccNEfHziNgVEQPA13mtyaoPmJV5m5l4kywzq3OSWLly5ZBzK1euZILGMpWslD6RuyWdWLFIRiDpkMzhmSQ1FIC1wNmSJks6nKTT3ysMm1ldi7QPJCvbR1JrSkkii0gSyWOSNkp6UNLGSgWW8cXMZy0ClgFExEPAGuBHJDssnh8RuyYgHjOzisgnkK6uLtrb2xkYGKC9vZ2urq7dEkutKGWP9dMqFsUoIuKjo1y7jKSfxMys7kmitbWV9vb2wSasfNNWa2trlaMrrJQ91rekM9bnAPtkLhXsJzEzs9JEBLlcjohA0uDPfELp6empdoi7KTqJSPok0E7Sgd0LnAzcDbyrMqGZmTWPXC5Hf3//YMLIN221trbW7BwRKK1PpB04EdgSEYuA44DnKxKVmVkTiQj6+/sH+z6yfSP9/f0126kOpfWJvBwRL0tC0uSIeETS3IpFZmbWJLJ9H11dXXR1dQEM6RupVaXURPoktQLXA7dJugHPyzAzK4t6mx+SV8oCjGdGRH9E5IBLSWaIf7BSgZmZNZN6mx+SN65NqSLizohYGxGvljsgM7NmM9b8kFpOJKWMzpoM/CEwO/u6iPjr8odlZtY8xpofUstNWqV0rN8AbAM2AK9UJhwzs+aUnR8CDJkfUstKSSIzI2JxxSIxM2tCwxPH8ONaV0qfyH9LOqZikZiZNZl624CqkFKSyNuB+9P9zCdyAUYzs4ZTzxMMs0ppzlpMsjlVfZTMzKyG1fMEw6wxayKS/jN9+hDwIK/tOvgQr+3tYWZmJarXCYZZYyaRiHh7+nNKRByQeUyJiAMqH6KZWWOq1wmGWeOabGhmZnumnicYZpUy2bCjwOltwIaI6C1fSGZmja+eJxhmldKxfkL6+H56/D7gXuA8Sd+NiC+WOzgzs0ZWrxMMs0ppzvpt4PiI+GxEfJYkocwA3gH86Z4EIeksSQ9JGpB0wrBrF0vanA4tPjVzfnF6brOki/bk883MqmV4wqinBAKlJZHDgOyCi78B3hgRL7Hny6BsAj4E3JU9KelI4GzgKJIhxldIapHUAnyFZN/3I4Fz0nvNzGwCldKc9W1gXbqPCMAHgO9I2g/40Z4EEREPQ8EMfAZwTUS8AjwhaTNwUnptc0Q8nr7umvTePYrDzMxKU3QSiYjPS7qJZOa6gPMi4r708kcqERxwKLAuc9yXngN4atj5BRWKwczMRlBKTQTgcaAF2Af4LUnviIi7xngNAJJ+ABxc4NIlEXFDgfOQJKvhgsLNcCOOh5O0BFgCMGPGDHp6ekYPtoHs2LGjqcoLLnOzcJlrQylDfD8JtAMzgV7gZOBu4F3FvD4iThlHfH3ArMzxTF7bknek84U+ezWwGmDu3LmxcOHCcYRSn3p6emim8oLL3CzqrczZUViFjotRi2UupWO9HTgR2BIRi4DjgOcrEtVr1gJnS5os6XBgDnAPydDiOZIOl7Q3Sef72grHYmY2Lo2wWu9ISkkiL0fEy5DschgRjwBzyxGEpDMl9QFvA26UdAtARDwErCHpML8ZOD8idkXETuAC4BbgYWBNeq+ZWU1plNV6R1JKn0ifpFbgeuA2SVsZpQmpFBFxHXDdCNcuAy4rcP4m4KZyfL6ZWaU0ymq9Iym6JhIRZ0ZEf0TkgEuBK4EPViowM7NG0Qir9Y5kXAswRsSdEbE2Il4d+24zs+bWCKv1jqToJCLpBEnXSbo/3dlwo3c2NDMbWUQM6QNpa2ur29V6R1JKn8g/A39BsjHVQGXCMTNrDLlcjv7+flauXElrayttbW0ArFixoi5X6x1JKUnk+YjwMFozszFkR2RB0v+xdOlSuru7aW9vHzxX7wkESksiyyV9A7idzIKLEfG9skdlZlbHGn1EVlYpHevnAvNIVtP9QPp4fyWCMjOrd408IiurlJrIWyPimIpFYmbWQEYakdVoiaSUmsg679lhZja2Rtk/vRil1ETeDvyppMdJ+kQEREQcW5HIzMzqUH5hxfz+6Z2dnXW7f3oxSkkip5ImjgrFYmZW17LDenO5HAMDA3R0dNDa2koul2u4piwoojlL0nZJL5JsYftg+nMT8FD608ys6RVaaLGjo2PIQouNlkCgiJpIREyZiEDMzOpdswzrzRrX2llmZvaa/H4hQFMM680qdXtcMzPLyDZjFRp1tXTpUlatWtWwicRJxMxsD+RHXkUE3d3dg+fza2V1d3cP3tOIicRJxMxsD0li1apVQ5LIqlWrBq812rDeLCcRM7M9NNrs9EatgeTVRMe6pLMkPSRpQNIJmfOzJb0kqTd9/L/MtfmSHpS0WVK3Gvm3ZGY1a6zZ6Y2uVmoim4APAV8rcO2xiJhX4PxXgSXAOpK91hcD/16xCM3MCsjOTs/XOhp1dnohNZFEIuJhoOj/2JIOAQ6IiLvT42+R7PfuJGJmEy6Xyw2ZTNjIHenD1URz1hgOl/RDSXdK+v303KFAX+aevvScmVlVDE8YzZBAYAJrIpJ+ABxc4NIlEXHDCC97BjgsIn4haT5wvaSjSNbwGm7ENb0kLSFp+mLGjBn09PSUFHs927FjR1OVF1zmZuEy14YJSyIRcco4XvMK6S6KEbFB0mPAESQ1j5mZW2cCT4/yPquB1QBz586NhQsXlhpK3erp6aGZygsuc7NwmWtDTTdnSZohqSV9/iZgDvB4RDwDbJd0cjoq62PASLUZMzOrkJpIIpLOlNQHvA24UdIt6aV3ABslPQD8K3BeRPwyvfZp4BvAZuAx3KluZhU2fFmTRtpcarxqZXTWdcB1Bc5fC1w7wmvuA46ucGhmZsDQvUIkDc4Pye8V0qxqoiZiZlbLIoKtW7cO2Stk6dKlQ/YKaVY1URMxM6tlK1asAJJFFbN7hSxYsKBp5oOMxDURM7NR5Jd6zy6umLdgwYIqRFRbXBMxMxvFSEu9W8I1ETOzcWhra6O7u3uwj6RZOYmYmY0gmxzWr1+/2/W2tramWGRxNG7OMjMrYPny5Wzbto3Ozk46OjpYv3498+bN4/TTT2fbtm10dXXR1tbG8uXLqx1qVTmJmJkNs3z5ctauXUtvby8ABxxwANOnT6e3t5d3vvOddHZ2As2x1PtYnETMzDIigm3bttHb28u8efMGh/MCzJs3j87OTiZNmtT0Q3vznETMzDKym0plEwjAhg0bmDRp0uB95o51M7PdSBpsssrq6Oho6pFYhTiJmJkNMzAwwPz584ecyzdtNfuQ3uGcRMzMMiKCjo6OwT6RXbt20d7ePng8depUN2VluE/EzCxDEq2trbS3tw/pRAeYOnXq4DpalnASMTMbJpfLERGDNY58Z7trILtzc5aZWQHDE4YTSGFOImZmNm5OImbW1Lzl7Z5xEjGzppXL5Vi6dOlg4sjvWNjM292WqiaSiKQvSXpE0kZJ10lqzVy7WNJmSY9KOjVzfnF6brOki6oTuZnVq4jg5ptvpru7ezCRLF26lO7ubm6++WbXSIpUE0kEuA04OiKOBX4MXAwg6UjgbOAoYDFwhaQWSS3AV4DTgCOBc9J7zcyKlt+ZsLu7m0mTJg1uOuUdC4tXE0kkIm6NiJ3p4TpgZvr8DOCaiHglIp4ANgMnpY/NEfF4RLwKXJPea2ZWFEmsWrWKtra2Iefb2tpYtWqVR2MVqRbniXwC+Jf0+aEkSSWvLz0H8NSw8yN+dZC0BFiSHr4iaVN5Qq0L04EXqh3EBHOZm0O5yjwLOCh/0N3d/Vx3d/dTo9xfTdX6Pb9xpAsTlkQk/QA4uMClSyLihvSeS4CdwD/nX1bg/qBwDWrEBsyIWA2sTj/jvog4oYTQ61qzlRdc5mbhMteGCUsiEXHKaNclfRx4P/DueK1Hq4/kW0LeTODp9PlI583MbILURJ+IpMXA54DTI+LXmUtrgbMlTZZ0ODAHuAe4F5gj6XBJe5N0vq+d6LjNzJpdrfSJfBmYDNyWdmati4jzIuIhSWuAH5E0c50fEbsAJF0A3AK0AFdFxENFftbqskdf25qtvOAyNwuXuQbIY6HNzGy8aqI5y8zM6pOTiJmZjVvTJRFJOUk/k9SbPt5b7ZgmiqQLJYWk6QTmFkkAAAahSURBVNWOpdIkfT5dRqdX0q2S3lDtmCpttOWDGpWksyQ9JGlAUk0NfS2nWl7mqemSSGplRMxLHzdVO5iJIGkW8B7gp9WOZYJ8KSKOjYh5wL8B/7vaAU2AgssHNbhNwIeAu6odSKXU+jJPzZpEmtFK4H8xyqTMRhIRL2YO96MJyj3K8kENKyIejohHqx1HhdX0Mk/NmkQuSKv8V0maVu1gKk3S6cDPIuKBascykSRdJukp4CM0R00k6xPAv1c7CCuLQ9l9madDR7h3wtXKPJGyGm2JFeCrwOdJvpl+Hvg7kv/h6toYZf5L4A8mNqLKG2spnYi4BLhE0sXABcDyCQ2wAsa5fFBdK6bMDW6k5Z9qQkMmkbGWWMmT9HWS9vK6N1KZJR0DHA48kE7knAncL+mkiHh2AkMsu2J/z8C3gRtpgCQyzuWD6loJv+dGNdryT1XXdM1Zkg7JHJ5J0jHXsCLiwYg4KCJmR8Rskn+Qx9d7AhmLpDmZw9OBR6oVy0QZZfkgq281vcxTQ9ZExvBFSfNIqoNPAp+qbjhWIV+QNBcYALYA51U5nolQcPmg6oZUWZLOBP4emAHcKKk3Ik4d42V1JSJ27sEyTxXnZU/MzGzcmq45y8zMysdJxMzMxs1JxMzMxs1JxMzMxs1JxMzMxs1JxMzMxs1JxOqGpF3p0u6bJH1/T5Y6l7SjiM/4rqTfKuE9WyV9psh7PyXpWUkPSHpM0seKeM2+ku5MV3VF0jGStkj6dOaevSXdJWm3OWCSZkt6SVJvsWUaIY6cpAuHnfuapN8bIeZeSa82wxYEzchJxOrJS+ny/UcDvwTOr/BnvEqRkxSVzO47ECgqiQDHArmIeCtwDtBZxGs+AXwvInZBshoByezlwQSUrvJ6O/DHI7zHY+ny+LvFL2lP/h4sIFk5eIiIeCn9vJpZpsPKy0nE6tXdpCuZSvoTSfek33i/lv+mnl67XtKGdOOiJSV+xn8AvzPS+6Tf7B+WdAVwP3Al8OY0ji+N8d7HAPklzJ8gSVj5mA+XdIOk+9JyzU0vfQQYvuDgc8BRw85dn947qgLxzxrpv5ekS9JNkX4AzB32Pm8h2b9kH0k3prWrTZJGSmTWSCLCDz/q4gHsSH+2AN8FFgNvAb4PvC69dgXwscxrDkx/7kuyTtpvZ99rlM/Yi+QP9qdHeh9gNsmyKien12YDm4osy1bgDSQrtK4Azk3Pv46kJvHm9Pi9wD8AewPPFnif7wKvAG/MnGsBni9w75D4hsc/SjnnAw8CvwUcAGwGLsy8poOklvSHwNcz56dmnj8JTK/2vyE/yv9oxrWzrH7tm7bnzwY2kOzk92mSP3L3putF7Uvy7TyvLV1fCZKVUOcAvyjiMyCpiVw5yvs8C2yJiN2acUaT7jI5BbiJpDa1Ecillz9IUrO4Ni3PXmkc04H+Ye+zmGTDrRvT12wBiIhdaR/ElIjYPkY4w+MvVM6TgesiXdRR0vDF/04FzgX2By6X9LfAv0XEf4zx2dYAnESsnrwUEfMkTSVZwv98koU0r46I3baClbQQOAV4W0T8WlIPsE8xn1HC+/xqHOU4FrgrIt6lZFO0TcDbgP8G3kqyT8aV2Rek9+2TOd4H+CLJCsXnAkeTJKW8ycDLRcQyGP8Y5Sy4yF468KA1Ip5Oj+eT1J7+r6RbI+Kvi4jB6pj7RKzuRMQ2oA24kGRv7T+SdBCApAMlvTG9dSqwNf2D+Lsk36jHo9j32U5Swxgk6XZJw3ehOwb4YVqWrST7nbwvvfYMcGq+kzsdgaX0vpY0eQD8FfCtiHiSpKnp6Mxn/jZJc9ZvylTOu4Az05FWU4APZF6zCLgj/dw3AL+OiH8CLgeOL/HzrQ45iVhdiogfAg+QfKv/K+BWSRtJmrjye8bcDOyVnv88BUYPFamo94mIXwD/lXYqfylNBL9DMpIsazCJpL5P8u0d4CqS/y8fTpvVPhcR+VrArcDb04729wCr0vNDkgjJH/ZsraRYBcsZEfcD/wL0AteSNK/lnZa+Ll+ue9K4LwH+zzhisDrjpeDNKkTS0cAnIqKjTO93HNARER8d477vARdHxKPDzs8m6as4utDrxhnT/cCCsWo9kp4EToiIF8r12VYbXBMxq5CI2FSuBJK+3w+BO7JDmIdTsvPd9cMTSGoXMHVPJxsOi+n40RJIfrIhyaizgXJ9rtUO10TMzGzcXBMxM7NxcxIxM7NxcxIxM7NxcxIxM7NxcxIxM7NxcxIxM7NxcxIxM7NxcxIxM7Nx+//R0qMkxrXukAAAAABJRU5ErkJggg==", "text/plain": [ "
" ] }, "metadata": { "needs_background": "light" }, "output_type": "display_data" } ], "source": [ "fig = plt.figure()\n", "plt.scatter(eigenvalues_trim[:, 0], eigenvalues_trim[:, 1],\n", " marker='x',\n", " color='k')\n", "plt.xlim(-5, 0.5)\n", "plt.ylim(-200, 200)\n", "plt.grid()\n", "plt.xlabel('Real Part, $Re(\\lambda)$ [rad/s]')\n", "plt.ylabel('Imaginary Part, $Im(\\lambda)$ [rad/s]');" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.7.5" } }, "nbformat": 4, "nbformat_minor": 4 }