Source code for scenarioDragSensitivity

#
#  ISC License
#
#  Copyright (c) 2016, Autonomous Vehicle Systems Lab, University of Colorado at Boulder
#
#  Permission to use, copy, modify, and/or distribute this software for any
#  purpose with or without fee is hereby granted, provided that the above
#  copyright notice and this permission notice appear in all copies.
#
#  THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES
#  WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF
#  MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR
#  ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES
#  WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR PROFITS, WHETHER IN AN
#  ACTION OF CONTRACT, NEGLIGENCE OR OTHER TORTIOUS ACTION, ARISING OUT OF
#  OR IN CONNECTION WITH THE USE OR PERFORMANCE OF THIS SOFTWARE.
#

r"""
Overview
--------

This script executes a parametric analysis of the control law examined in :ref:`scenarioDragRendezvous`,
considering the performance of that control
law across both increasing initial displacements and variations in atmospheric density. 

This script is found in the folder ``src/examples`` and executed by using::

      python3 scenarioDragSensitivity.py

The simulation layout is identical to that used in :ref:`scenarioDragRendezvous`. The simulator used in
that scenario is run using a grid
of true anomaly offsets and atmospheric densities using Python's multiprocessing library, and a
set of surface plots reflecting the controls' terminal performance
and trajectories are produced.

Illustration of Simulation Results
----------------------------------

In this scenario, the differential drag scenario used in :ref:`scenarioDragRendezvous` is examined
across a range of initial along-track orbit offsets and atmospheric densities.
The resulting Hill-frame trajectories corresponding to every fifth simulation run are shown in the following image.

.. image:: /_images/Scenarios/scenarioDragSensitivity_hillTrajectories.svg
   :align: center

To visualize the sensitivity of terminal position and velocity errors to both increasing baseline and
variations in density, the following surface plots - which show the
scale of terminal errors as a function of atmospheric density and maneuver baseline - are shown below:

.. image:: /_images/Scenarios/scenarioDragSensitivity_positionSensitivity.svg
   :align: center

.. image:: /_images/Scenarios/scenarioDragSensitivity_velocitySensitivity.svg
   :align: center

"""

import multiprocessing as mp
import os
import pickle

import numpy as np
from matplotlib import pyplot as plt

from scenarioDragRendezvous import drag_simulator

fileName = os.path.basename(os.path.splitext(__file__)[0])

[docs] def sim_wrapper(arg): """ Because multiprocessing.map only considers single-input functions as targets for maps, this function maps args and kwargs from the dicts provided by multiprocessing to the inputs demanded by drag_simulator. """ args, kwargs = arg result = drag_simulator(*args, **kwargs) return result
def drag_sensitivity_analysis(ctrlType, nuOffsetNum, densityNum, rerunSims=False): alt_offsets= [0] # np.arange(-100,100,10) nu_offsets = np.arange(0.001, 1, (1-0.001)/nuOffsetNum) density_multipliers = np.logspace(-1, 0.4, num=densityNum) pool = mp.Pool(mp.cpu_count()) X,Y,Z = np.meshgrid(alt_offsets, nu_offsets, density_multipliers) positions = np.vstack([X.ravel(), Y.ravel(), Z.ravel()]).T kwargs = {'ctrlType': ctrlType} arg = [(position, kwargs) for position in positions] if rerunSims: sim_results = pool.map(sim_wrapper, arg) with open(ctrlType+"_sweep_results.pickle", "wb") as output_file: pickle.dump(sim_results, output_file,-1) else: with open(ctrlType+"_sweep_results.pickle", "rb") as fp: sim_results = pickle.load(fp) pool.close()
[docs] def results_to_ranges_and_plot(results_list): """ Converts a results dict from scenarioDragRendezvous to a set of initial and final distance and speed errors, and returns a plot of all of the Hill-frame trajectories. """ fig = plt.figure() pos_errs = np.empty(len(results_list)) vel_errs = np.empty(len(results_list)) init_relpos = np.empty(len(results_list)) init_relvel = np.empty(len(results_list)) dens_list = np.empty(len(results_list)) density_colorvals = np.logspace(-1,0.4,num=20) # normalizer = cmap.norma for ind, result in enumerate(results_list): hill_pos = result['relState.r_DC_H'] hill_vel = result['relState.v_DC_H'] init_relpos[ind] = np.linalg.norm(hill_pos[1,:]) init_relvel[ind] = np.linalg.norm(hill_vel[1,:]) pos_errs[ind] = np.linalg.norm(hill_pos[-1,:]) vel_errs[ind] = np.linalg.norm(hill_vel[-1,:]) dens_list[ind] = result['dens_mult'] if ind%5==0: plt.plot(hill_pos[:,0],hill_pos[:,1] ) plt.grid() plt.xlabel('Hill X (m)') plt.ylabel('Hill Y (m)') return init_relpos, init_relvel, pos_errs, vel_errs, dens_list, fig
def comparison_sweep(savePickle): with open("lqr_sweep_results.pickle", "rb") as fp: lqr_sim_results = pickle.load(fp) if not savePickle: os.remove("lqr_sweep_results.pickle") figlist = {} lqr_init_range, lqr_init_vel, lqr_err_range, lqr_err_vel, lqr_dens, lqr_fig = results_to_ranges_and_plot(lqr_sim_results) figlist[fileName+'_hillTrajectories'] = lqr_fig unique_ranges = np.unique(lqr_init_range.round(decimals=2)) x_shape = len(unique_ranges) unique_dens = np.unique(lqr_dens) y_shape = len(unique_dens) import matplotlib.colors as colors fig = plt.figure() X,Y = np.meshgrid(unique_ranges, unique_dens) Z = lqr_err_range.reshape([x_shape, y_shape]) Z_vel = lqr_err_vel.reshape([x_shape, y_shape]) # Position error contours contour = plt.contourf(X.T,Y.T,Z, 30,norm=colors.LogNorm(Z.min(), Z.max())) plt.ylabel('Density Multiplier') plt.xlabel('Initial Displacement (m)') cbar = fig.colorbar(contour) cbar.set_label("Terminal Positioning Error (m)") figlist[fileName+'_positionSensitivity'] = fig # Velocity error contours fig2 = plt.figure() contour2 = plt.contourf(X.T,Y.T,Z_vel, 30,norm=colors.LogNorm(Z.min(), Z.max())) plt.ylabel('Density Multiplier') plt.xlabel('Initial Displacement (m)') cbar2 = fig2.colorbar(contour2) cbar2.set_label("Terminal Velocity Error (m/s)") figlist[fileName+'_velocitySensitivity'] = fig2 return figlist def run(show_plots, nuOffsetNum, densityNum, rerunSims=False, savePickle=False): drag_sensitivity_analysis('lqr', nuOffsetNum, densityNum, rerunSims=rerunSims) plots = comparison_sweep(savePickle) if show_plots: plt.show() return plots if __name__=='__main__': nuOffsetNum = 15 densityNum = 20 run(True, nuOffsetNum, densityNum, rerunSims=True, savePickle=False)