diff --git a/projectile_path/.gitignore b/projectile_path/.gitignore
index cab3285..29d4d74 100644
--- a/projectile_path/.gitignore
+++ b/projectile_path/.gitignore
@@ -130,4 +130,8 @@ dmypy.json
# Other folders not for github
_not_gh_tests
-_not_gh_contracts
\ No newline at end of file
+_not_gh_contracts
+
+# Added by cargo
+
+/target
diff --git a/projectile_path/tests/prototype_projectile_plot.ipynb b/projectile_path/tests/prototype_projectile_plot.ipynb
index 9c02c67..3d50563 100644
--- a/projectile_path/tests/prototype_projectile_plot.ipynb
+++ b/projectile_path/tests/prototype_projectile_plot.ipynb
@@ -10,7 +10,7 @@
},
{
"cell_type": "code",
- "execution_count": 8,
+ "execution_count": 2,
"id": "7d7134bf",
"metadata": {},
"outputs": [],
@@ -29,7 +29,7 @@
},
{
"cell_type": "code",
- "execution_count": 9,
+ "execution_count": 3,
"id": "ac10d3a4",
"metadata": {},
"outputs": [],
@@ -52,7 +52,7 @@
},
{
"cell_type": "code",
- "execution_count": 10,
+ "execution_count": 4,
"id": "2494abc5",
"metadata": {},
"outputs": [],
@@ -94,7 +94,7 @@
},
{
"cell_type": "code",
- "execution_count": 11,
+ "execution_count": 5,
"id": "92ea27fe",
"metadata": {},
"outputs": [],
@@ -172,7 +172,7 @@
},
{
"cell_type": "code",
- "execution_count": 12,
+ "execution_count": 6,
"id": "c29aebb5",
"metadata": {},
"outputs": [],
@@ -239,7 +239,7 @@
},
{
"cell_type": "code",
- "execution_count": 13,
+ "execution_count": 7,
"id": "df8f26d3",
"metadata": {},
"outputs": [
@@ -308,7 +308,7 @@
},
{
"cell_type": "code",
- "execution_count": 14,
+ "execution_count": 9,
"id": "f9a52462",
"metadata": {
"scrolled": false
@@ -344,6 +344,14 @@
"\n",
"plt.show()"
]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "id": "899ae8fc",
+ "metadata": {},
+ "outputs": [],
+ "source": []
}
],
"metadata": {
diff --git a/projectile_path/tests/prototype_projectile_plot_rust.ipynb b/projectile_path/tests/prototype_projectile_plot_rust.ipynb
new file mode 100644
index 0000000..9e8dba2
--- /dev/null
+++ b/projectile_path/tests/prototype_projectile_plot_rust.ipynb
@@ -0,0 +1,977 @@
+{
+ "cells": [
+ {
+ "cell_type": "markdown",
+ "id": "da964f3c",
+ "metadata": {},
+ "source": [
+ "## Inputs"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 2,
+ "id": "e760fbb0",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "//\n",
+ "// Cairo Inputs\n",
+ "//\n",
+ "\n",
+ "// Number of points to plot\n",
+ "const NUM_PTS: i32 = 20;\n",
+ "\n",
+ "// Launch angle in degrees: -179 <= THETA_0_DEG <= +180\n",
+ "// Must be integer\n",
+ "const THETA_0_DEG: i32 = 105;\n",
+ "\n",
+ "// Launch velocity magnitude: V_0 = ~100 is enough to reach top of plot area on vertical shot\n",
+ "const V_0: f64 = 100.0;"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "027a1621",
+ "metadata": {},
+ "source": [
+ "## Constants"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 3,
+ "id": "30d127c8",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "//\n",
+ "// Math constants\n",
+ "//\n",
+ "\n",
+ "// Number of terms in cosine Taylor series approximation\n",
+ "const N: i32 = 5;\n",
+ "\n",
+ "// pi approximation to match that used in Cairo contract\n",
+ "const PI_APPROX: f64 = 3.141592654;\n",
+ "\n",
+ "\n",
+ "//\n",
+ "// Physical parameters\n",
+ "//\n",
+ "\n",
+ "// Gravitational acceleration in m/s^2\n",
+ "const G: f64 = 9.8;\n",
+ "\n",
+ "// Initial position\n",
+ "const X_0: f64 = 0.0;\n",
+ "const Y_0: f64 = 0.0;\n",
+ "\n",
+ "\n",
+ "//\n",
+ "// Plot parameters\n",
+ "//\n",
+ "\n",
+ "// Min and max values for axes\n",
+ "const X_MAX: f64 = 1000.0;\n",
+ "const X_MIN: f64 = -X_MAX;\n",
+ "const Y_MAX: f64 = 500.0;\n",
+ "const Y_MIN: f64 = -Y_MAX;"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "51910e29",
+ "metadata": {},
+ "source": [
+ "## Math functions"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 4,
+ "id": "1da6fbb0",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "//\n",
+ "// Math functions\n",
+ "//\n",
+ "\n",
+ "// Factorial function\n",
+ "fn factorial(n: i32) -> f64 {\n",
+ " let mut fact = 1.0;\n",
+ " for i in 1..n + 1 {\n",
+ " fact *= i as f64;\n",
+ " }\n",
+ " fact\n",
+ "}\n",
+ "\n",
+ "\n",
+ "// Cosine from Taylor series approximation\n",
+ "// Assume -pi <= theta <= +pi, so no need to shift theta\n",
+ "// Accuracy is improved if instead -pi/2 <= theta <= +pi/2\n",
+ "fn cosine_n_terms(theta: f64, n: i32) -> f64 {\n",
+ " // n = number of terms (not order)\n",
+ " // 2(n-1) = order\n",
+ " // cosine(theta) ~= ((-1)^n)*(theta^(2n))/(2n)!\n",
+ " // ~= 1 - theta^2/2! + theta^4/4! - theta^6/6! + ...\n",
+ " let mut cos_nth = 0.0;\n",
+ " for i in 0..n {\n",
+ " let neg_one: f64 = -1.0;\n",
+ " let power_neg_one = neg_one.powi(i);\n",
+ " let power_theta = theta.powi(2 * i);\n",
+ " let fact = factorial(2 * i);\n",
+ " cos_nth += power_neg_one * power_theta / fact;\n",
+ " }\n",
+ " cos_nth\n",
+ "}\n",
+ "\n",
+ "\n",
+ "// Cosine approximation\n",
+ "// Taylor series approximation is more accurate if -pi/2 <= theta_0 <= +pi/2. So:\n",
+ "// If theta_0 is in 2nd/3rd quadrant:\n",
+ "// (1) move angle to 1st/4th quadrant for cosine approximation, and\n",
+ "// (2) force negative sign for cosine(theta_0)\n",
+ "// If +90 or -90 degrees, use exact value of cos_theta_0\n",
+ "// (Use theta_0_deg for comparisons because calculated theta_0 in radians is slightly rounded)\n",
+ "// Then call cosine_n_terms\n",
+ "fn cosine_approx(theta_0: f64, theta_0_deg: i32, n: i32) -> f64 {\n",
+ " match theta_0_deg {\n",
+ " -179..=-91 => -cosine_n_terms(-PI_APPROX - theta_0, n),\n",
+ " -90 => 0.0,\n",
+ " -89..=89 => cosine_n_terms(theta_0, n),\n",
+ " 90 => 0.0,\n",
+ " 91..=180 => -cosine_n_terms(PI_APPROX - theta_0, n),\n",
+ " _ => {\n",
+ " println!(\"Error: theta_0_deg = {theta_0_deg} outside allowed range -179 to +180\");\n",
+ " 0.0\n",
+ " }\n",
+ " }\n",
+ "}\n",
+ "\n",
+ "\n",
+ "// Sine approximation: need to force correct signs\n",
+ "fn sine_approx(theta_0: f64, cos_theta_0: f64) -> f64 {\n",
+ " if theta_0 >= 0.0 {\n",
+ " (1.0 - cos_theta_0.powi(2)).sqrt()\n",
+ " } else {\n",
+ " -((1.0 - cos_theta_0.powi(2)).sqrt())\n",
+ " }\n",
+ "}"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "4ff3fde5",
+ "metadata": {},
+ "source": [
+ "## Physics functions"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 5,
+ "id": "3c230674",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "//\n",
+ "// Physics functions\n",
+ "//\n",
+ "// Time of projectile in plot area\n",
+ "fn time_in_plot(\n",
+ " theta_0_deg: i32,\n",
+ " x_comps: (f64, f64),\n",
+ " y_comps: (f64, f64, f64),\n",
+ " plot_limits: (f64, f64, f64, f64)\n",
+ ") -> f64 {\n",
+ "\n",
+ " let (x_0, v_0x) = x_comps;\n",
+ " let (y_0, v_0y, g) = y_comps;\n",
+ " let (x_min, x_max, y_min, _y_max) = plot_limits;\n",
+ "\n",
+ " // Max time needed for y-direction\n",
+ " let t_max_y = (v_0y + (v_0y.powi(2) - 2.0 * g * (y_min - y_0)).sqrt()) / g;\n",
+ "\n",
+ " // Find max time needed for x_direction, t_max_x\n",
+ " // Then t_max is minimum of t_max_x and t_max_y\n",
+ " let t_max_x: f64;\n",
+ " let t_max: f64;\n",
+ "\n",
+ " if theta_0_deg.abs() == 90 {\n",
+ " // abs(theta_0_deg) = 90, so v_0x = 0, so t_max_x = infinite, so\n",
+ " t_max = t_max_y;\n",
+ " } else {\n",
+ " match theta_0_deg {\n",
+ " // abs(theta_0_deg) < 90, so v_0x > 0, so projectile moves toward x_max\n",
+ " -89..=89 => t_max_x = (x_max - x_0) / v_0x,\n",
+ " // abs(theta_0_deg) > 90, so v_0x < 0, so projectile moves toward x_min\n",
+ " -179..=-91 | 91..=180 => t_max_x = (x_min - x_0) / v_0x,\n",
+ " _ => {\n",
+ " println!(\"Error: theta_0_deg = {theta_0_deg} outside allowed range -179 to +180\");\n",
+ " t_max_x = 0.0;\n",
+ " }\n",
+ " }\n",
+ " \n",
+ " if t_max_x > t_max_y {\n",
+ " t_max = t_max_y;\n",
+ " } else {\n",
+ " t_max = t_max_x;\n",
+ " }\n",
+ " }\n",
+ "\n",
+ " println!(\"t_max = {t_max}\");\n",
+ " t_max\n",
+ "}\n",
+ "\n",
+ "\n",
+ "fn x_value(x_comps: (f64, f64), t: f64) -> f64 {\n",
+ " let (x_0, v_0x) = x_comps;\n",
+ " x_0 + v_0x * t\n",
+ "}\n",
+ "\n",
+ "\n",
+ "fn y_value(y_comps: (f64, f64, f64), t: f64) -> f64 {\n",
+ " let (y_0, v_0y, g) = y_comps;\n",
+ " y_0 + v_0y * t - 0.5 * g * t.powi(2)\n",
+ "}"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "cae78c44",
+ "metadata": {},
+ "source": [
+ "## Projectile calculations"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 6,
+ "id": "656a02c6",
+ "metadata": {
+ "scrolled": true
+ },
+ "outputs": [
+ {
+ "name": "stdout",
+ "output_type": "stream",
+ "text": [
+ "theta_0 = 1.8325957148333336\n",
+ "v_0x = -25.882306286627742\n",
+ "v_0y = 96.59247497235584\n",
+ "t_max = 23.969806345067223\n",
+ "x_s = [0.0, -32.65230891864121, -65.30461783728242, -97.95692675592362, -130.60923567456484, -163.26154459320605, -195.91385351184724, -228.56616243048848, -261.2184713491297, -293.87078026777084, -326.5230891864121, -359.1753981050533, -391.8277070236945, -424.4800159423357, -457.13232486097695, -489.78463377961816, -522.4369426982594, -555.0892516169006, -587.7415605355417, -620.393869454183]\n",
+ "y_s = [0.0, 114.05942534793961, 212.52160460458768, 295.38653776994425, 362.6542248440092, 414.32466582678273, 450.39786071826467, 470.873809518455, 475.7525122273538, 465.0339688449611, 438.71817937127696, 396.8051438063011, 339.2948621500341, 266.187334402475, 177.48256056362447, 73.18054063348268, -46.718725387950826, -182.2152375006758, -333.3089957046923, -500.00000000000045]\n"
+ ]
+ }
+ ],
+ "source": [
+ "// fn main() {\n",
+ "\n",
+ " // Check that THETA_0_DEG is within allowed range\n",
+ " if THETA_0_DEG > 180 || THETA_0_DEG < -179 {\n",
+ " println!(\"Error: THETA_0_DEG = {THETA_0_DEG} outside allowed range -179 to +180\");\n",
+ " }\n",
+ " \n",
+ " // Calculate theta in radians\n",
+ " let theta_0 = THETA_0_DEG as f64 * PI_APPROX / 180.0;\n",
+ " // Cosine approximation\n",
+ " let cos_theta_0 = cosine_approx(theta_0, THETA_0_DEG, N);\n",
+ " // Sine approximation\n",
+ " let sin_theta_0 = sine_approx(theta_0, cos_theta_0);\n",
+ "\n",
+ " // Projectile velocity components in x- and y-directions\n",
+ " let v_0x = V_0 * cos_theta_0;\n",
+ " let v_0y = V_0 * sin_theta_0;\n",
+ "\n",
+ " println!(\"theta_0 = {theta_0}\");\n",
+ " println!(\"v_0x = {v_0x}\");\n",
+ " println!(\"v_0y = {v_0y}\");\n",
+ " \n",
+ " // Tuples\n",
+ " let x_comps: (f64, f64) = (X_0, v_0x);\n",
+ " let y_comps: (f64, f64, f64) = (Y_0, v_0y, G);\n",
+ " let plot_limits: (f64, f64, f64, f64) = (X_MIN, X_MAX, Y_MIN, Y_MAX);\n",
+ "\n",
+ " // Calculate time in plot\n",
+ " let t_max = time_in_plot(THETA_0_DEG, x_comps, y_comps, plot_limits);\n",
+ " let delta_t = t_max / ((NUM_PTS as f64) - 1.0);\n",
+ "\n",
+ "// // Fill position vectors\n",
+ "// let mut x_s = Vec::new();\n",
+ "// let mut y_s = Vec::new();\n",
+ "\n",
+ "// for i in 0..NUM_PTS {\n",
+ "// let t = i as f64 * delta_t;\n",
+ "// let x = x_value(x_comps, t);\n",
+ "// let y = y_value(y_comps, t);\n",
+ "// x_s.push(x);\n",
+ "// y_s.push(y);\n",
+ "// }\n",
+ "\n",
+ "// Fill position arrays\n",
+ " let mut x_s: [f64; NUM_PTS as usize] = [0.0; NUM_PTS as usize];\n",
+ " let mut y_s: [f64; NUM_PTS as usize] = [0.0; NUM_PTS as usize];\n",
+ "\n",
+ " for i in 0..NUM_PTS as usize {\n",
+ " let t = i as f64 * delta_t;\n",
+ " x_s[i] = x_value(x_comps, t);\n",
+ " y_s[i] = y_value(y_comps, t);\n",
+ " }\n",
+ "\n",
+ " // Need to use {:?}, Debug print, for printing Vec or Array\n",
+ " // because {}, Display print, works only for single values\n",
+ " println!(\"x_s = {x_s:?}\");\n",
+ " println!(\"y_s = {y_s:?}\");\n",
+ "\n",
+ "// }"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "id": "278b8e18",
+ "metadata": {},
+ "source": [
+ "## Plot with Plotters"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 7,
+ "id": "deb3690f",
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ ":dep plotters = { version = \"^0.3.4\", default_features = false, features = [\"evcxr\", \"all_series\", \"all_elements\"] }\n",
+ "extern crate plotters;\n",
+ "use plotters::prelude::*;"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": 8,
+ "id": "ab3f029b",
+ "metadata": {},
+ "outputs": [
+ {
+ "data": {
+ "text/html": [
+ "\n",
+ "
\n",
+ "
\n",
+ " \n",
+ "