[1]:

import fmdtools.sim.propagate as propagate
import fmdtools.analyze as an
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import pandas as pd
import numpy as np

Urban Drone Demo

Model Overview

The drone model is defined in drone_mdl_urban.py, along with some visualization functions.

[2]:

from drone_mdl_urban import Drone

[3]:

mdl = Drone()

This is the model structure:

[4]:

mg = an.graph.FunctionArchitectureGraph(mdl)
fig, ax = mg.draw()
../../_images/examples_multirotor_Urban_Drone_Demo_5_0.png

We can also view the grid environment using its show methods:

[5]:

collections={"all_occupied": {"color": "red"}, "all_allowed": {"color": "blue"}, "start": {"color": "blue"}, "end": {"color": "blue"}}
fig, ax = mdl.flows['environment'].c.show("height", collections=collections)
../../_images/examples_multirotor_Urban_Drone_Demo_7_0.png

Which shows the Start, End, and allowed/unsafe locations in the 1000x1000-m grid. In this display, line thickness corresponds to building height, and hatching corresponds to whether or not the space is occupied. We can also display this using show.coord3d:

[6]:

fig, ax = mdl.flows['environment'].c.show_z("height", collections=collections)
../../_images/examples_multirotor_Urban_Drone_Demo_9_0.png

Nominal Simulation

Below we show how this drone performs in the nominal scenario.

[7]:

results_nom, hist_nom =propagate.nominal(mdl)
fig, axs = hist_nom.plot_line("fxns.store_ee.s.soc",
                            'flows.ee_1.s.rate',
                            'flows.dofs.s.z',
                            'fxns.plan_path.m.mode',
                            'fxns.plan_path.s.pt',
                            'flows.des_traj.s.power')
../../_images/examples_multirotor_Urban_Drone_Demo_11_0.png

As shown, the flight ends fairly quickly (in 10 minutes), with the drone successively proceeding through points in the flight plan.

We can also view this flightpath in 3-d space using History.plot_trajectories:

[8]:

fig, ax = hist_nom.plot_trajectories('dofs.s.x', 'dofs.s.y', 'dofs.s.z')
../../_images/examples_multirotor_Urban_Drone_Demo_13_0.png

Trajectory plots can be overlaid on top of environment plots. In this case, we defined the method plot_env_with_traj and plot_env_with_traj_z for this case.

[9]:

from drone_mdl_urban import plot_env_with_traj, plot_env_with_traj_z

[10]:

plot_env_with_traj_z(hist_nom, mdl)

[10]:

(<Figure size 400x400 with 1 Axes>, <Axes3D: title={'center': 'trajectory'}>)
../../_images/examples_multirotor_Urban_Drone_Demo_16_1.png 
[11]:

plot_env_with_traj(hist_nom, mdl)

[11]:

(<Figure size 640x480 with 2 Axes>,
 <Axes: title={'center': 'trajectory'}, xlabel='x', ylabel='y'>)
../../_images/examples_multirotor_Urban_Drone_Demo_17_1.png

As shown, this is a rather simple straight-line path. If we wanted a more complex scenario, we could make the path more complex by adding multiple destinations or planning the path based on allowed flight/landing locations.

As it is, we may also want to adjust the timestep/speed to get more resolution, since the drone only has a few discrete timesteps in the air.

The results for the simulation are:

[12]:

results_nom

[12]:

endclass:
--rate:                              1.0
--cost:                              0.0
--expected_cost:                     0.0
--repcost:                             0
--unsafe_flight_time:                  0
--body_strikes:                      0.0
--head_strikes:                      0.0
--property_restrictions:               0
--safecost:                          0.0
--landcost:                            0
--p_safety:                          0.0
--severities: {'hazardous': 0.0, 'minor': 1.0}

Resilience model

A number of different faults have been implemented in the system.

For example, here we inject a mechanical fault in the left-rear rotor during flight:

[13]:

mdl.fxns['affect_dof'].m.faultmodes

[13]:

{'lf_short': Fault(prob=1.0000000000000002e-06, cost=200, phases={}, units='sim'),
 'lf_openc': Fault(prob=1.0000000000000002e-06, cost=200, phases={}, units='sim'),
 'lf_ctlup': Fault(prob=2.0000000000000003e-06, cost=500, phases={}, units='sim'),
 'lf_ctldn': Fault(prob=2.0000000000000003e-06, cost=500, phases={}, units='sim'),
 'lf_ctlbreak': Fault(prob=2.0000000000000003e-06, cost=1000, phases={}, units='sim'),
 'lf_mechbreak': Fault(prob=1.0000000000000002e-06, cost=500, phases={}, units='sim'),
 'lf_mechfriction': Fault(prob=5.000000000000001e-07, cost=500, phases={}, units='sim'),
 'lf_propwarp': Fault(prob=1.0000000000000001e-07, cost=200, phases={}, units='sim'),
 'lf_propstuck': Fault(prob=2.0000000000000002e-07, cost=200, phases={}, units='sim'),
 'lf_propbreak': Fault(prob=3.0000000000000004e-07, cost=200, phases={}, units='sim'),
 'lr_short': Fault(prob=1.0000000000000002e-06, cost=200, phases={}, units='sim'),
 'lr_openc': Fault(prob=1.0000000000000002e-06, cost=200, phases={}, units='sim'),
 'lr_ctlup': Fault(prob=2.0000000000000003e-06, cost=500, phases={}, units='sim'),
 'lr_ctldn': Fault(prob=2.0000000000000003e-06, cost=500, phases={}, units='sim'),
 'lr_ctlbreak': Fault(prob=2.0000000000000003e-06, cost=1000, phases={}, units='sim'),
 'lr_mechbreak': Fault(prob=1.0000000000000002e-06, cost=500, phases={}, units='sim'),
 'lr_mechfriction': Fault(prob=5.000000000000001e-07, cost=500, phases={}, units='sim'),
 'lr_propwarp': Fault(prob=1.0000000000000001e-07, cost=200, phases={}, units='sim'),
 'lr_propstuck': Fault(prob=2.0000000000000002e-07, cost=200, phases={}, units='sim'),
 'lr_propbreak': Fault(prob=3.0000000000000004e-07, cost=200, phases={}, units='sim'),
 'rf_short': Fault(prob=1.0000000000000002e-06, cost=200, phases={}, units='sim'),
 'rf_openc': Fault(prob=1.0000000000000002e-06, cost=200, phases={}, units='sim'),
 'rf_ctlup': Fault(prob=2.0000000000000003e-06, cost=500, phases={}, units='sim'),
 'rf_ctldn': Fault(prob=2.0000000000000003e-06, cost=500, phases={}, units='sim'),
 'rf_ctlbreak': Fault(prob=2.0000000000000003e-06, cost=1000, phases={}, units='sim'),
 'rf_mechbreak': Fault(prob=1.0000000000000002e-06, cost=500, phases={}, units='sim'),
 'rf_mechfriction': Fault(prob=5.000000000000001e-07, cost=500, phases={}, units='sim'),
 'rf_propwarp': Fault(prob=1.0000000000000001e-07, cost=200, phases={}, units='sim'),
 'rf_propstuck': Fault(prob=2.0000000000000002e-07, cost=200, phases={}, units='sim'),
 'rf_propbreak': Fault(prob=3.0000000000000004e-07, cost=200, phases={}, units='sim'),
 'rr_short': Fault(prob=1.0000000000000002e-06, cost=200, phases={}, units='sim'),
 'rr_openc': Fault(prob=1.0000000000000002e-06, cost=200, phases={}, units='sim'),
 'rr_ctlup': Fault(prob=2.0000000000000003e-06, cost=500, phases={}, units='sim'),
 'rr_ctldn': Fault(prob=2.0000000000000003e-06, cost=500, phases={}, units='sim'),
 'rr_ctlbreak': Fault(prob=2.0000000000000003e-06, cost=1000, phases={}, units='sim'),
 'rr_mechbreak': Fault(prob=1.0000000000000002e-06, cost=500, phases={}, units='sim'),
 'rr_mechfriction': Fault(prob=5.000000000000001e-07, cost=500, phases={}, units='sim'),
 'rr_propwarp': Fault(prob=1.0000000000000001e-07, cost=200, phases={}, units='sim'),
 'rr_propstuck': Fault(prob=2.0000000000000002e-07, cost=200, phases={}, units='sim'),
 'rr_propbreak': Fault(prob=3.0000000000000004e-07, cost=200, phases={}, units='sim')}

[14]:

results_fault, hist_fault =propagate.one_fault(mdl, "affect_dof", "lr_mechbreak", time=2.0)

[15]:

fig, axs = hist_fault.plot_line("fxns.store_ee.s.soc",
                                'flows.ee_1.s.rate',
                                'flows.dofs.s.z',
                                'fxns.plan_path.m.mode',
                                'fxns.plan_path.s.pt',
                                'flows.des_traj.s.power')
../../_images/examples_multirotor_Urban_Drone_Demo_23_0.png 
[16]:

plot_env_with_traj_z(hist_fault, mdl)

[16]:

(<Figure size 400x400 with 1 Axes>, <Axes3D: title={'center': 'trajectory'}>)
../../_images/examples_multirotor_Urban_Drone_Demo_24_1.png 
[17]:

plot_env_with_traj(hist_fault, mdl)

[17]:

(<Figure size 640x480 with 2 Axes>,
 <Axes: title={'center': 'trajectory'}, xlabel='x', ylabel='y'>)
../../_images/examples_multirotor_Urban_Drone_Demo_25_1.png 
[18]:

results_fault

[18]:

endclass.rate:    1.0000000000000002e-06
endclass.cost:                   28600.0
endclass.expected_cost: 2860.0000000000005
endclass.repcost:                    500
endclass.unsafe_flight_time:         281
endclass.body_strikes:               0.0
endclass.head_strikes:               0.0
endclass.property_restrictions:        0
endclass.safecost:               28100.0
endclass.landcost:                     0
endclass.p_safety:                   0.0
endclass.severities: {'hazardous': 0.0, 'minor': 1.0000000000000002e-06}

Here we inject a large list of faults in the system and evaluate their relative consequences in terms of metrics calculated in find_classification:

[19]:

from fmdtools.analyze.phases import PhaseMap, from_hist
phasemaps = from_hist(hist_nom, fxn_modephases=[])
phasemaps

[19]:

{'ctl_dof': PhaseMap({'nominal': [0.0, 30.0]}, {}),
 'plan_path': PhaseMap({'taxi': [0.0, 1.0], 'move': [1.1, 1.6], 'land': [1.7, 2.4], 'taxi1': [2.5, 30.0]}, {})}

[20]:

from drone_mdl_urban import make_move_quad
move_quad=make_move_quad(hist_nom, phasemaps['plan_path'].phases['move'])
move_quad

[20]:

{'samp': 'quadrature',
 'quad': {'nodes': [-0.20000000000000018, 1.0], 'weights': [1.0, 0.92]}}

[21]:

from fmdtools.sim.sample import FaultDomain, FaultSample


fd = FaultDomain(mdl)
fd.add_all()

fs = FaultSample(fd, phasemap = phasemaps['plan_path'])
fs.add_fault_phases("move", method = "quad",
                    args=(move_quad['quad']['nodes'], move_quad['quad']['weights']))

fs

[21]:

FaultSample of scenarios:
 - manage_health_lostfunction_t2
 - manage_health_lostfunction_t2
 - store_ee_nocharge_t2
 - store_ee_nocharge_t2
 - store_ee_lowcharge_t2
 - store_ee_lowcharge_t2
 - store_ee_s1p1_short_t2
 - store_ee_s1p1_short_t2
 - store_ee_s1p1_degr_t2
 - store_ee_s1p1_degr_t2
 - ... (118 total)

[22]:

endresults, hists = propagate.fault_sample(mdl, fs, staged=False)

SCENARIOS COMPLETE: 100%|██████████| 118/118 [01:29<00:00,  1.32it/s]

[23]:

plot_env_with_traj_z(hists , mdl)

[23]:

(<Figure size 400x400 with 1 Axes>, <Axes3D: title={'center': 'trajectory'}>)
../../_images/examples_multirotor_Urban_Drone_Demo_32_1.png 
[24]:

statsfmea = an.tabulate.FMEA(endresults, fs,
                            weight_metrics=['rate'],
                            avg_metrics=['unsafe_flight_time', 'cost', 'repcost',
                                         'landcost', 'body_strikes',
                                         'head_strikes', 'property_restrictions'])
fmeatab = statsfmea.as_table(sort_by="cost")
fmeatab

[24]:
rate unsafe_flight_time cost repcost landcost body_strikes head_strikes property_restrictions
plan_path vision_lack_of_detection 5.989583e-03 281.0 3.377386e+06 0.0 0.0 0.0 0.0 0.0
dist_ee degr 2.395833e-06 281.0 1.303438e+04 1000.0 0.0 0.0 0.0 0.0
short 1.437500e-06 281.0 1.269100e+04 1500.0 0.0 0.0 0.0 0.0
break 9.583333e-07 281.0 1.240733e+04 1500.0 0.0 0.0 0.0 0.0
affect_dof lf_ctlbreak 9.583333e-07 281.0 1.219775e+04 1000.0 0.0 0.0 0.0 0.0
lr_ctlbreak 9.583333e-07 281.0 1.219775e+04 1000.0 0.0 0.0 0.0 0.0
rf_ctlbreak 9.583333e-07 281.0 1.219775e+04 1000.0 0.0 0.0 0.0 0.0
rr_ctlbreak 9.583333e-07 281.0 1.219775e+04 1000.0 0.0 0.0 0.0 0.0
lr_ctlup 9.583333e-07 281.0 1.198817e+04 500.0 0.0 0.0 0.0 0.0
lr_ctldn 9.583333e-07 281.0 1.198817e+04 500.0 0.0 0.0 0.0 0.0
rf_ctldn 9.583333e-07 281.0 1.198817e+04 500.0 0.0 0.0 0.0 0.0
rf_ctlup 9.583333e-07 281.0 1.198817e+04 500.0 0.0 0.0 0.0 0.0
rr_ctlup 9.583333e-07 281.0 1.198817e+04 500.0 0.0 0.0 0.0 0.0
lf_ctldn 9.583333e-07 281.0 1.198817e+04 500.0 0.0 0.0 0.0 0.0
rr_ctldn 9.583333e-07 281.0 1.198817e+04 500.0 0.0 0.0 0.0 0.0
lf_ctlup 9.583333e-07 281.0 1.198817e+04 500.0 0.0 0.0 0.0 0.0
hold_payload deform 2.875000e-09 281.0 1.184170e+04 1500.0 0.0 0.0 0.0 0.0
break 7.187500e-10 281.0 1.184043e+04 1500.0 0.0 0.0 0.0 0.0
affect_dof lr_mechbreak 4.791667e-07 281.0 1.171408e+04 500.0 0.0 0.0 0.0 0.0
lf_mechbreak 4.791667e-07 281.0 1.171408e+04 500.0 0.0 0.0 0.0 0.0
rr_mechbreak 4.791667e-07 281.0 1.171408e+04 500.0 0.0 0.0 0.0 0.0
rf_mechbreak 4.791667e-07 281.0 1.171408e+04 500.0 0.0 0.0 0.0 0.0
plan_path degloc 2.875000e-08 281.0 1.165673e+04 1000.0 0.0 0.0 0.0 0.0
ctl_dof degctl 2.875000e-08 281.0 1.165673e+04 1000.0 0.0 0.0 0.0 0.0
noctl 7.187500e-09 281.0 1.164418e+04 1000.0 0.0 0.0 0.0 0.0
plan_path noloc 7.187500e-09 281.0 1.164418e+04 1000.0 0.0 0.0 0.0 0.0
manage_health lostfunction 1.796875e-10 281.0 1.164010e+04 1000.0 0.0 0.0 0.0 0.0
affect_dof lr_short 4.791667e-07 281.0 1.159121e+04 200.0 0.0 0.0 0.0 0.0
lf_openc 4.791667e-07 281.0 1.159121e+04 200.0 0.0 0.0 0.0 0.0
rr_openc 4.791667e-07 281.0 1.159121e+04 200.0 0.0 0.0 0.0 0.0
rf_short 4.791667e-07 281.0 1.159121e+04 200.0 0.0 0.0 0.0 0.0
rf_openc 4.791667e-07 281.0 1.159121e+04 200.0 0.0 0.0 0.0 0.0
lf_short 4.791667e-07 281.0 1.159121e+04 200.0 0.0 0.0 0.0 0.0
lr_openc 4.791667e-07 281.0 1.159121e+04 200.0 0.0 0.0 0.0 0.0
rr_short 4.791667e-07 281.0 1.159121e+04 200.0 0.0 0.0 0.0 0.0
rf_mechfriction 2.395833e-07 281.0 1.157704e+04 500.0 0.0 0.0 0.0 0.0
lf_mechfriction 2.395833e-07 281.0 1.157704e+04 500.0 0.0 0.0 0.0 0.0
rr_mechfriction 2.395833e-07 281.0 1.157704e+04 500.0 0.0 0.0 0.0 0.0
lr_mechfriction 2.395833e-07 281.0 1.157704e+04 500.0 0.0 0.0 0.0 0.0
lr_propbreak 1.437500e-07 281.0 1.140136e+04 200.0 0.0 0.0 0.0 0.0
lf_propbreak 1.437500e-07 281.0 1.140136e+04 200.0 0.0 0.0 0.0 0.0
rr_propbreak 1.437500e-07 281.0 1.140136e+04 200.0 0.0 0.0 0.0 0.0
rf_propbreak 1.437500e-07 281.0 1.140136e+04 200.0 0.0 0.0 0.0 0.0
store_ee lowcharge 1.677083e-07 281.0 1.137459e+04 100.0 0.0 0.0 0.0 0.0
affect_dof lf_propstuck 9.583333e-08 281.0 1.137424e+04 200.0 0.0 0.0 0.0 0.0
lr_propstuck 9.583333e-08 281.0 1.137424e+04 200.0 0.0 0.0 0.0 0.0
rf_propstuck 9.583333e-08 281.0 1.137424e+04 200.0 0.0 0.0 0.0 0.0
rr_propstuck 9.583333e-08 281.0 1.137424e+04 200.0 0.0 0.0 0.0 0.0
store_ee s1p1_nocharge 1.437500e-07 281.0 1.136108e+04 100.0 0.0 0.0 0.0 0.0
affect_dof rr_propwarp 4.791667e-08 281.0 1.134712e+04 200.0 0.0 0.0 0.0 0.0
rf_propwarp 4.791667e-08 281.0 1.134712e+04 200.0 0.0 0.0 0.0 0.0
lr_propwarp 4.791667e-08 281.0 1.134712e+04 200.0 0.0 0.0 0.0 0.0
lf_propwarp 4.791667e-08 281.0 1.134712e+04 200.0 0.0 0.0 0.0 0.0
store_ee s1p1_lowcharge 9.583333e-08 281.0 1.133405e+04 100.0 0.0 0.0 0.0 0.0
s1p1_short 7.187500e-08 281.0 1.132054e+04 100.0 0.0 0.0 0.0 0.0
s1p1_break 7.187500e-08 281.0 1.132054e+04 100.0 0.0 0.0 0.0 0.0
s1p1_degr 7.187500e-08 281.0 1.132054e+04 100.0 0.0 0.0 0.0 0.0
nocharge 4.791667e-08 281.0 1.130702e+04 100.0 0.0 0.0 0.0 0.0
plan_path vision_undesired_detection 5.989583e-03 0.0 0.000000e+00 0.0 0.0 0.0 0.0 0.0 
[ ]: