adipu/guadaloop_lev_control
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

useOfAi

Browse Source
This commit is contained in:
pulipakaa24
2025-12-12 13:42:41 -06:00
parent e5f479e266
commit a20f560cd7
8 changed files with 1735 additions and 465 deletions
Download Patch File Download Diff File Expand all files Collapse all files

174
RL Testing/Function Fitting.ipynb
View File

@@ -10,7 +10,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 30, "execution_count": 2,
"id": "e4ff9106", "id": "e4ff9106",
"metadata": {}, "metadata": {},
"outputs": [], "outputs": [],
@@ -24,7 +24,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 31, "execution_count": 4,
"id": "32946653", "id": "32946653",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -47,60 +47,6 @@
"memory usage: 209.7 KB\n", "memory usage: 209.7 KB\n",
"\n", "\n",
"After adding mirrored data:\n", "After adding mirrored data:\n",
"Total rows: 8281\n",
"Roll angle range: -4.0° to 4.0°\n",
"Unique roll angles: [np.float64(-4.0), np.float64(-3.0), np.float64(-2.0), np.float64(-1.5), np.float64(-1.0), np.float64(-0.5), np.float64(0.0), np.float64(0.5), np.float64(1.0), np.float64(1.5), np.float64(2.0), np.float64(3.0), np.float64(4.0)]\n"
]
}
],
"source": [
"magDf = pd.read_csv(\"Ansys Results 12-9.csv\")\n",
"magDf.info()\n",
"\n",
"# Condition: Drop rows where 'GapHeight' is less than 5\n",
"condition = magDf['GapHeight [mm]'] < 5\n",
"\n",
"# Drop the rows that satisfy the condition (i.e., keep rows where the condition is False)\n",
"magDf = magDf[~condition]\n",
"\n",
"# Create mirrored data with negative roll angles\n",
"# For non-zero roll angles, create symmetric data by:\n",
"# - Negating roll angle\n",
"# - Swapping currL and currR\n",
"# - Negating torque (force remains the same due to symmetry)\n",
"\n",
"# Filter for non-zero roll angles\n",
"non_zero_roll = magDf[magDf['rollDeg [deg]'] != 0]\n",
"\n",
"# Create the mirrored rows\n",
"mirrored_data = non_zero_roll.copy()\n",
"mirrored_data['rollDeg [deg]'] = -non_zero_roll['rollDeg [deg]']\n",
"mirrored_data['currL [A]'], mirrored_data['currR [A]'] = non_zero_roll['currR [A]'].values, non_zero_roll['currL [A]'].values\n",
"mirrored_data['YokeTorque.Torque [mNewtonMeter]'] = -non_zero_roll['YokeTorque.Torque [mNewtonMeter]']\n",
"# Force remains unchanged\n",
"\n",
"# Combine original and mirrored data\n",
"magDf = pd.concat([magDf, mirrored_data], ignore_index=True)\n",
"\n",
"# Sort by rollDeg for better organization\n",
"magDf = magDf.sort_values(['rollDeg [deg]', 'currL [A]', 'currR [A]', 'GapHeight [mm]']).reset_index(drop=True)\n",
"\n",
"print(f\"\\nAfter adding mirrored data:\")\n",
"print(f\"Total rows: {len(magDf)}\")\n",
"print(f\"Roll angle range: {magDf['rollDeg [deg]'].min():.1f}° to {magDf['rollDeg [deg]'].max():.1f}°\")\n",
"print(f\"Unique roll angles: {sorted(magDf['rollDeg [deg]'].unique())}\")"
]
},
{
"cell_type": "code",
"execution_count": 32,
"id": "4fe774dd",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"<class 'pandas.core.frame.DataFrame'>\n", "<class 'pandas.core.frame.DataFrame'>\n",
"RangeIndex: 8281 entries, 0 to 8280\n", "RangeIndex: 8281 entries, 0 to 8280\n",
"Data columns (total 6 columns):\n", "Data columns (total 6 columns):\n",
@@ -118,6 +64,26 @@
} }
], ],
"source": [ "source": [
"magDf = pd.read_csv(\"Ansys Results 12-9.csv\")\n",
"magDf.info()\n",
"\n",
"condition = magDf['GapHeight [mm]'] < 5\n",
"\n",
"magDf = magDf[~condition]\n",
"\n",
"nzRoll = magDf[magDf['rollDeg [deg]'] != 0]\n",
"\n",
"mirrored_data = nzRoll.copy()\n",
"mirrored_data['rollDeg [deg]'] = -nzRoll['rollDeg [deg]']\n",
"mirrored_data['currL [A]'], mirrored_data['currR [A]'] = nzRoll['currR [A]'].values, nzRoll['currL [A]'].values\n",
"mirrored_data['YokeTorque.Torque [mNewtonMeter]'] = -nzRoll['YokeTorque.Torque [mNewtonMeter]']\n",
"\n",
"magDf = pd.concat([magDf, mirrored_data], ignore_index=True)\n",
"\n",
"magDf = magDf.sort_values(['rollDeg [deg]', 'currL [A]', 'currR [A]', 'GapHeight [mm]']).reset_index(drop=True)\n",
"\n",
"print()\n",
"print(\"After adding mirrored data:\")\n",
"magDf.info()" "magDf.info()"
] ]
}, },
@@ -126,12 +92,12 @@
"id": "ecd1f124", "id": "ecd1f124",
"metadata": {}, "metadata": {},
"source": [ "source": [
"### After removing non-model rows, we have 4459 rows, which matches the number of simulations conducted in Ansys." "### After removing non-model rows, we have 4459 rows * ~2 for the mirrored roll angles, which matches the number of simulations conducted in Ansys."
] ]
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 33, "execution_count": null,
"id": "1a76017e", "id": "1a76017e",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -155,7 +121,9 @@
"plt.title(r\"$-15A, -15A, 0^\\circ roll$\")\n", "plt.title(r\"$-15A, -15A, 0^\\circ roll$\")\n",
"plt.xlabel(\"Inverse Gap Height (1/mm)\")\n", "plt.xlabel(\"Inverse Gap Height (1/mm)\")\n",
"plt.ylabel(\"Force on magnetic yoke (N)\")\n", "plt.ylabel(\"Force on magnetic yoke (N)\")\n",
"plt.show()" "plt.show()\n",
"\n",
"# Experimented here with ways to represent gap-height-to-force relationship."
] ]
}, },
{ {
@@ -307,7 +275,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 35, "execution_count": null,
"id": "bd16a120", "id": "bd16a120",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -386,6 +354,7 @@
"from sklearn.preprocessing import PolynomialFeatures\n", "from sklearn.preprocessing import PolynomialFeatures\n",
"from sklearn.linear_model import LinearRegression\n", "from sklearn.linear_model import LinearRegression\n",
"\n", "\n",
"# The following was written by AI - see [4]\n",
"# 1. Load your Ansys CSV - using invGap for modeling\n", "# 1. Load your Ansys CSV - using invGap for modeling\n",
"X = magDf[['currL [A]', 'currR [A]', 'rollDeg [deg]', 'invGap']]\n", "X = magDf[['currL [A]', 'currR [A]', 'rollDeg [deg]', 'invGap']]\n",
"y = magDf[['YokeForce.Force_z [newton]', 'YokeTorque.Torque [mNewtonMeter]']]\n", "y = magDf[['YokeForce.Force_z [newton]', 'YokeTorque.Torque [mNewtonMeter]']]\n",
@@ -413,7 +382,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 36, "execution_count": null,
"id": "dc79d9fe", "id": "dc79d9fe",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -437,21 +406,17 @@
} }
], ],
"source": [ "source": [
"# Get predictions from the model\n",
"y_pred = model.predict(X_poly)\n", "y_pred = model.predict(X_poly)\n",
"\n", "\n",
"# Extract force and torque predictions\n",
"force_pred = y_pred[:, 0]\n", "force_pred = y_pred[:, 0]\n",
"torque_pred = y_pred[:, 1]\n", "torque_pred = y_pred[:, 1]\n",
"\n", "\n",
"# Extract actual values\n",
"force_actual = y.iloc[:, 0].values\n", "force_actual = y.iloc[:, 0].values\n",
"torque_actual = y.iloc[:, 1].values\n", "torque_actual = y.iloc[:, 1].values\n",
"\n", "\n",
"# Create comparison plots\n",
"fig, axes = plt.subplots(1, 2, figsize=(14, 5))\n", "fig, axes = plt.subplots(1, 2, figsize=(14, 5))\n",
"\n", "\n",
"# Plot 1: Force comparison\n", "# Force comparison\n",
"axes[0].scatter(force_actual, force_pred, alpha=0.5, s=10)\n", "axes[0].scatter(force_actual, force_pred, alpha=0.5, s=10)\n",
"axes[0].plot([force_actual.min(), force_actual.max()], \n", "axes[0].plot([force_actual.min(), force_actual.max()], \n",
" [force_actual.min(), force_actual.max()], \n", " [force_actual.min(), force_actual.max()], \n",
@@ -462,7 +427,7 @@
"axes[0].legend()\n", "axes[0].legend()\n",
"axes[0].grid(True, alpha=0.3)\n", "axes[0].grid(True, alpha=0.3)\n",
"\n", "\n",
"# Plot 2: Torque comparison\n", "# Torque comparison\n",
"axes[1].scatter(torque_actual, torque_pred, alpha=0.5, s=10)\n", "axes[1].scatter(torque_actual, torque_pred, alpha=0.5, s=10)\n",
"axes[1].plot([torque_actual.min(), torque_actual.max()], \n", "axes[1].plot([torque_actual.min(), torque_actual.max()], \n",
" [torque_actual.min(), torque_actual.max()], \n", " [torque_actual.min(), torque_actual.max()], \n",
@@ -476,7 +441,6 @@
"plt.tight_layout()\n", "plt.tight_layout()\n",
"plt.show()\n", "plt.show()\n",
"\n", "\n",
"# Calculate R² scores\n",
"from sklearn.metrics import r2_score\n", "from sklearn.metrics import r2_score\n",
"print(f\"Force R² Score: {r2_score(force_actual, force_pred):.6f}\")\n", "print(f\"Force R² Score: {r2_score(force_actual, force_pred):.6f}\")\n",
"print(f\"Torque R² Score: {r2_score(torque_actual, torque_pred):.6f}\")" "print(f\"Torque R² Score: {r2_score(torque_actual, torque_pred):.6f}\")"
@@ -492,7 +456,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 37, "execution_count": null,
"id": "837a814f", "id": "837a814f",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -508,7 +472,7 @@
} }
], ],
"source": [ "source": [
"# For the first subset (same as you plotted earlier)\n", "# For the first subset\n",
"subset_indices = range(0, 13)\n", "subset_indices = range(0, 13)\n",
"x_subset = magDf.iloc[subset_indices][\"GapHeight [mm]\"]\n", "x_subset = magDf.iloc[subset_indices][\"GapHeight [mm]\"]\n",
"y_actual_subset = magDf.iloc[subset_indices][\"YokeForce.Force_z [newton]\"]\n", "y_actual_subset = magDf.iloc[subset_indices][\"YokeForce.Force_z [newton]\"]\n",
@@ -548,7 +512,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 38, "execution_count": null,
"id": "fd12d56b", "id": "fd12d56b",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -590,7 +554,6 @@
"from sklearn.model_selection import cross_val_score\n", "from sklearn.model_selection import cross_val_score\n",
"from sklearn.pipeline import Pipeline\n", "from sklearn.pipeline import Pipeline\n",
"\n", "\n",
"# Split data for proper validation\n",
"from sklearn.model_selection import train_test_split\n", "from sklearn.model_selection import train_test_split\n",
"X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)\n", "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)\n",
"\n", "\n",
@@ -602,26 +565,21 @@
"test_scores_torque = []\n", "test_scores_torque = []\n",
"\n", "\n",
"for degree in degrees:\n", "for degree in degrees:\n",
" # Create pipeline\n",
" pipe = Pipeline([\n", " pipe = Pipeline([\n",
" ('poly', PolynomialFeatures(degree=degree)),\n", " ('poly', PolynomialFeatures(degree=degree)),\n",
" ('model', LinearRegression())\n", " ('model', LinearRegression())\n",
" ])\n", " ])\n",
" \n", " \n",
" # Fit on training data\n",
" pipe.fit(X_train, y_train)\n", " pipe.fit(X_train, y_train)\n",
" \n", " \n",
" # Score on both sets\n",
" train_score = pipe.score(X_train, y_train)\n", " train_score = pipe.score(X_train, y_train)\n",
" test_score = pipe.score(X_test, y_test)\n", " test_score = pipe.score(X_test, y_test)\n",
" \n", " \n",
" # Store scores (using overall R² for multi-output)\n",
" train_scores_force.append(train_score)\n", " train_scores_force.append(train_score)\n",
" test_scores_force.append(test_score)\n", " test_scores_force.append(test_score)\n",
" \n", " \n",
" print(f\"Degree {degree}: Train R² = {train_score:.6f}, Test R² = {test_score:.6f}\")\n", " print(f\"Degree {degree}: Train R² = {train_score:.6f}, Test R² = {test_score:.6f}\")\n",
"\n", "\n",
"# Plot the results\n",
"plt.figure(figsize=(10, 6))\n", "plt.figure(figsize=(10, 6))\n",
"plt.plot(degrees, train_scores_force, 'o-', label='Training Score', linewidth=2, markersize=8)\n", "plt.plot(degrees, train_scores_force, 'o-', label='Training Score', linewidth=2, markersize=8)\n",
"plt.plot(degrees, test_scores_force, 's-', label='Test Score', linewidth=2, markersize=8)\n", "plt.plot(degrees, test_scores_force, 's-', label='Test Score', linewidth=2, markersize=8)\n",
@@ -633,7 +591,6 @@
"plt.xticks(degrees)\n", "plt.xticks(degrees)\n",
"plt.show()\n", "plt.show()\n",
"\n", "\n",
"# Find best degree\n",
"best_degree = degrees[np.argmax(test_scores_force)]\n", "best_degree = degrees[np.argmax(test_scores_force)]\n",
"print(f\"\\nBest polynomial degree: {best_degree}\")\n", "print(f\"\\nBest polynomial degree: {best_degree}\")\n",
"print(f\"Best test R² score: {max(test_scores_force):.6f}\")" "print(f\"Best test R² score: {max(test_scores_force):.6f}\")"
@@ -661,7 +618,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 40, "execution_count": null,
"id": "d94aa4c1", "id": "d94aa4c1",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -694,6 +651,7 @@
} }
], ],
"source": [ "source": [
"# The following was written by AI - see [5]\n",
"# Train best degree polynomial on full dataset\n", "# Train best degree polynomial on full dataset\n",
"poly_best = PolynomialFeatures(degree=best_degree)\n", "poly_best = PolynomialFeatures(degree=best_degree)\n",
"X_poly_best = poly_best.fit_transform(X)\n", "X_poly_best = poly_best.fit_transform(X)\n",
@@ -794,6 +752,7 @@
} }
], ],
"source": [ "source": [
"# The following was written by AI - see [6]\n",
"# Select 4 different conditions to visualize\n", "# Select 4 different conditions to visualize\n",
"# Each condition has 13 gap height measurements (one complete parameter set)\n", "# Each condition has 13 gap height measurements (one complete parameter set)\n",
"condition_indices = [\n", "condition_indices = [\n",
@@ -874,7 +833,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 42, "execution_count": null,
"id": "35f4438e", "id": "35f4438e",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -896,6 +855,7 @@
} }
], ],
"source": [ "source": [
"# The following was written by AI - see [7]\n",
"# Calculate R² score for each 13-row segment (each parameter set)\n", "# Calculate R² score for each 13-row segment (each parameter set)\n",
"num_segments = len(magDf) // 13\n", "num_segments = len(magDf) // 13\n",
"segment_scores = []\n", "segment_scores = []\n",
@@ -953,7 +913,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 43, "execution_count": null,
"id": "6ceaf2b2", "id": "6ceaf2b2",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -977,6 +937,7 @@
} }
], ],
"source": [ "source": [
"# The following was written by AI - see [8]\n",
"# Get worst segment data\n", "# Get worst segment data\n",
"worst_data = magDf.iloc[worst_segment['start']:worst_segment['end']]\n", "worst_data = magDf.iloc[worst_segment['start']:worst_segment['end']]\n",
"worst_currL = worst_segment['currL']\n", "worst_currL = worst_segment['currL']\n",
@@ -1059,7 +1020,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 44, "execution_count": null,
"id": "c120048d", "id": "c120048d",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -1085,14 +1046,12 @@
} }
], ],
"source": [ "source": [
"# Calculate overall torque prediction accuracy\n",
"y_pred_full = model_best.predict(X_poly_best)\n", "y_pred_full = model_best.predict(X_poly_best)\n",
"torque_actual_full = y.iloc[:, 1].values\n", "torque_actual_full = y.iloc[:, 1].values\n",
"torque_pred_full = y_pred_full[:, 1]\n", "torque_pred_full = y_pred_full[:, 1]\n",
"\n", "\n",
"torque_r2_overall = r2_score(torque_actual_full, torque_pred_full)\n", "torque_r2_overall = r2_score(torque_actual_full, torque_pred_full)\n",
"\n", "\n",
"# Create scatter plot for torque\n",
"fig, ax = plt.subplots(1, 1, figsize=(10, 8))\n", "fig, ax = plt.subplots(1, 1, figsize=(10, 8))\n",
"\n", "\n",
"ax.scatter(torque_actual_full, torque_pred_full, alpha=0.4, s=15, color='#1f77b4')\n", "ax.scatter(torque_actual_full, torque_pred_full, alpha=0.4, s=15, color='#1f77b4')\n",
@@ -1127,7 +1086,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 45, "execution_count": null,
"id": "70b9b4e1", "id": "70b9b4e1",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -1151,6 +1110,7 @@
} }
], ],
"source": [ "source": [
"# The following was written by AI - see [9]\n",
"# Select 4 different gap heights to visualize\n", "# Select 4 different gap heights to visualize\n",
"# Using currL = -15A, currR = -15A configuration\n", "# Using currL = -15A, currR = -15A configuration\n",
"gap_heights_to_plot = [8, 10, 15, 21] # mm\n", "gap_heights_to_plot = [8, 10, 15, 21] # mm\n",
@@ -1227,7 +1187,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 46, "execution_count": null,
"id": "3b8d7f6a", "id": "3b8d7f6a",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -1253,6 +1213,7 @@
} }
], ],
"source": [ "source": [
"# The following was written by AI - see [10]\n",
"# Select 4 different gap heights to visualize\n", "# Select 4 different gap heights to visualize\n",
"# Using roll = 0.5 degrees (small but non-zero torque)\n", "# Using roll = 0.5 degrees (small but non-zero torque)\n",
"gap_heights_to_plot = [8, 10, 15, 21] # mm\n", "gap_heights_to_plot = [8, 10, 15, 21] # mm\n",
@@ -1346,7 +1307,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 47, "execution_count": null,
"id": "902343b2", "id": "902343b2",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -1379,6 +1340,7 @@
} }
], ],
"source": [ "source": [
"# Following generated by AI - see [11]\n",
"# Calculate R² score for torque in each 13-row segment\n", "# Calculate R² score for torque in each 13-row segment\n",
"torque_segment_scores = []\n", "torque_segment_scores = []\n",
"torque_segment_info = []\n", "torque_segment_info = []\n",
@@ -1450,7 +1412,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 48, "execution_count": null,
"id": "9d5580f3", "id": "9d5580f3",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -1473,6 +1435,7 @@
} }
], ],
"source": [ "source": [
"# Following generated by AI - see [11]\n",
"# Get worst torque segment data\n", "# Get worst torque segment data\n",
"worst_torque_data = magDf.iloc[worst_torque_segment['start']:worst_torque_segment['end']]\n", "worst_torque_data = magDf.iloc[worst_torque_segment['start']:worst_torque_segment['end']]\n",
"worst_torque_currL = worst_torque_segment['currL']\n", "worst_torque_currL = worst_torque_segment['currL']\n",
@@ -1551,7 +1514,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 49, "execution_count": null,
"id": "778f58df", "id": "778f58df",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -1590,6 +1553,7 @@
} }
], ],
"source": [ "source": [
"# The following was generated by AI - see [12]\n",
"print(\"=\" * 60)\n", "print(\"=\" * 60)\n",
"print(\"MODEL EXPORT SUMMARY\")\n", "print(\"MODEL EXPORT SUMMARY\")\n",
"print(\"=\" * 60)\n", "print(\"=\" * 60)\n",
@@ -1626,7 +1590,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 50, "execution_count": null,
"id": "edb239e1", "id": "edb239e1",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -1652,6 +1616,7 @@
} }
], ],
"source": [ "source": [
"# The following was generated by AI - see [12]\n",
"# Generate standalone predictor class\n", "# Generate standalone predictor class\n",
"predictor_code = f'''\"\"\"\n", "predictor_code = f'''\"\"\"\n",
"Magnetic Levitation Force and Torque Predictor\n", "Magnetic Levitation Force and Torque Predictor\n",
@@ -1829,7 +1794,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 51, "execution_count": null,
"id": "4a1ddacd", "id": "4a1ddacd",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -1864,6 +1829,7 @@
} }
], ],
"source": [ "source": [
"# The following was generated by AI - see [12]\n",
"# Import the standalone predictor (reload to get latest version)\n", "# Import the standalone predictor (reload to get latest version)\n",
"import importlib\n", "import importlib\n",
"import maglev_predictor\n", "import maglev_predictor\n",
@@ -1927,15 +1893,12 @@
"source": [ "source": [
"### Export Sklearn Model (Lossless & Faster)\n", "### Export Sklearn Model (Lossless & Faster)\n",
"\n", "\n",
"Save the actual sklearn model using joblib for lossless serialization. This is faster than the standalone predictor because:\n", "Save the actual sklearn model using joblib for lossless serialization. This is also faster than the standalone predictor"
"- Uses optimized sklearn implementations\n",
"- No manual feature generation overhead\n",
"- Smaller memory footprint for predictions"
] ]
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 52, "execution_count": null,
"id": "5fd33470", "id": "5fd33470",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -1974,6 +1937,7 @@
} }
], ],
"source": [ "source": [
"# The following was generated by AI - see [12]\n",
"import joblib\n", "import joblib\n",
"import pickle\n", "import pickle\n",
"\n", "\n",
@@ -2030,7 +1994,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 53, "execution_count": null,
"id": "889df820", "id": "889df820",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -2061,6 +2025,7 @@
} }
], ],
"source": [ "source": [
"# The following was generated by AI - see [12]\n",
"import time\n", "import time\n",
"\n", "\n",
"# Load the sklearn model\n", "# Load the sklearn model\n",
@@ -2141,7 +2106,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 54, "execution_count": null,
"id": "badbc379", "id": "badbc379",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -2164,6 +2129,7 @@
} }
], ],
"source": [ "source": [
"# The following was generated by AI - see [13]\n",
"# Find equilibrium gap height for 5.8 kg pod using polynomial root finding\n", "# Find equilibrium gap height for 5.8 kg pod using polynomial root finding\n",
"import numpy as np\n", "import numpy as np\n",
"from maglev_predictor import MaglevPredictor\n", "from maglev_predictor import MaglevPredictor\n",
@@ -2220,19 +2186,11 @@
" print(f\" Gap: {gap_val:.6f} mm → Force: {force:.6f} N\")\n", " print(f\" Gap: {gap_val:.6f} mm → Force: {force:.6f} N\")\n",
"print()" "print()"
] ]
},
{
"cell_type": "code",
"execution_count": null,
"id": "3ec4bb72",
"metadata": {},
"outputs": [],
"source": []
} }
], ],
"metadata": { "metadata": {
"kernelspec": { "kernelspec": {
"display_name": "LevSim", "display_name": "RLenv",
"language": "python", "language": "python",
"name": "python3" "name": "python3"
}, },
@@ -2246,7 +2204,7 @@
"name": "python", "name": "python",
"nbconvert_exporter": "python", "nbconvert_exporter": "python",
"pygments_lexer": "ipython3", "pygments_lexer": "ipython3",
"version": "3.10.19" "version": "3.10.9"
} }
}, },
"nbformat": 4, "nbformat": 4,

1658
RL Testing/Use of AI.txt
View File

File diff suppressed because it is too large Load Diff

111
RL Testing/find_equilibrium.py
View File

@@ -1,111 +0,0 @@
"""
Find equilibrium gap height for magnetic levitation system.
Given:
- Pod mass: 5.8 kg
- Required force: 5.8 * 9.81 = 56.898 N
- All currents: 0 A
- Roll angle: 0 degrees
Find: Gap height (mm) that produces this force
"""
import numpy as np
from maglev_predictor import MaglevPredictor
# Initialize predictor
predictor = MaglevPredictor()
# Target force
target_force = 5.8 * 9.81 # 56.898 N
print("=" * 70)
print("EQUILIBRIUM GAP HEIGHT FINDER")
print("=" * 70)
print(f"Target Force: {target_force:.3f} N (for 5.8 kg pod)")
print(f"Parameters: currL = 0 A, currR = 0 A, roll = 0°")
print()
# Define objective function
def force_error(gap_height):
"""Calculate difference between predicted force and target force."""
force, _ = predictor.predict(currL=0, currR=0, roll=0, gap_height=gap_height)
return force - target_force
# First, let's scan across gap heights to understand the behavior
print("Scanning gap heights from 5 to 30 mm...")
print("-" * 70)
gap_range = np.linspace(5, 30, 26)
forces = []
for gap in gap_range:
force, torque = predictor.predict(currL=0, currR=0, roll=0, gap_height=gap)
forces.append(force)
if abs(force - target_force) < 5: # Close to target
print(f"Gap: {gap:5.1f} mm → Force: {force:7.3f} N ← CLOSE!")
else:
print(f"Gap: {gap:5.1f} mm → Force: {force:7.3f} N")
forces = np.array(forces)
print()
print("-" * 70)
# Check if target force is within the range
min_force = forces.min()
max_force = forces.max()
if target_force < min_force or target_force > max_force:
print(f"⚠ WARNING: Target force {target_force:.3f} N is outside the range!")
print(f" Force range: {min_force:.3f} N to {max_force:.3f} N")
print(f" Cannot achieve equilibrium with 0A currents.")
else:
# Find the gap height using root finding
# Use brentq for robust bracketing (requires sign change)
# Find bracketing interval
idx_above = np.where(forces > target_force)[0]
idx_below = np.where(forces < target_force)[0]
if len(idx_above) > 0 and len(idx_below) > 0:
# Find the transition point
gap_low = gap_range[idx_below[-1]]
gap_high = gap_range[idx_above[0]]
print(f"✓ Target force is achievable!")
print(f" Bracketing interval: {gap_low:.1f} mm to {gap_high:.1f} mm")
print()
# Use simple bisection method for accurate root finding
tol = 1e-6
while (gap_high - gap_low) > tol:
gap_mid = (gap_low + gap_high) / 2
force_mid, _ = predictor.predict(currL=0, currR=0, roll=0, gap_height=gap_mid)
if force_mid > target_force:
gap_low = gap_mid
else:
gap_high = gap_mid
equilibrium_gap = (gap_low + gap_high) / 2
# Verify the result
final_force, final_torque = predictor.predict(
currL=0, currR=0, roll=0, gap_height=equilibrium_gap
)
print("=" * 70)
print("EQUILIBRIUM FOUND!")
print("=" * 70)
print(f"Equilibrium Gap Height: {equilibrium_gap:.6f} mm")
print(f"Predicted Force: {final_force:.6f} N")
print(f"Target Force: {target_force:.6f} N")
print(f"Error: {abs(final_force - target_force):.9f} N")
print(f"Torque at equilibrium: {final_torque:.6f} mN·m")
print()
print(f"✓ Pod will levitate at {equilibrium_gap:.3f} mm gap height")
print(f" with no current applied (permanent magnets only)")
print("=" * 70)
else:
print("⚠ Could not find bracketing interval for bisection.")
print(" Target force may not be achievable in the scanned range.")

11
RL Testing/lev_PPO.ipynb
View File

@@ -293,7 +293,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 7, "execution_count": null,
"id": "7ee8fb34", "id": "7ee8fb34",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -333,6 +333,7 @@
} }
], ],
"source": [ "source": [
"# The following was generated by AI - see [19]\n",
"environ = LevPodEnv(use_gui=False, initial_gap_mm=14, max_steps=500) # Start below target\n", "environ = LevPodEnv(use_gui=False, initial_gap_mm=14, max_steps=500) # Start below target\n",
"model = ActorCriticNetwork(environ.observation_space.shape[0], environ.action_space.shape[0]).to(device)\n", "model = ActorCriticNetwork(environ.observation_space.shape[0], environ.action_space.shape[0]).to(device)\n",
"train_data, reward, gap_error = rollout(model, environ)\n", "train_data, reward, gap_error = rollout(model, environ)\n",
@@ -380,7 +381,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 8, "execution_count": null,
"id": "e6f27ed4", "id": "e6f27ed4",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -415,6 +416,7 @@
} }
], ],
"source": [ "source": [
"# The following was generated by AI - see [20]\n",
"# Visualize the new reward function\n", "# Visualize the new reward function\n",
"gap_errors_mm = np.linspace(0, 20, 100)\n", "gap_errors_mm = np.linspace(0, 20, 100)\n",
"gap_rewards = np.exp(-0.5 * (gap_errors_mm / 3.0)**2)\n", "gap_rewards = np.exp(-0.5 * (gap_errors_mm / 3.0)**2)\n",
@@ -443,7 +445,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 9, "execution_count": null,
"id": "fb554183", "id": "fb554183",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [
@@ -457,6 +459,7 @@
} }
], ],
"source": [ "source": [
"# The following was generated by AI - see [21]\n",
"import os\n", "import os\n",
"from datetime import datetime\n", "from datetime import datetime\n",
"\n", "\n",
@@ -701,7 +704,7 @@
}, },
{ {
"cell_type": "code", "cell_type": "code",
"execution_count": 12, "execution_count": null,
"id": "3678193c", "id": "3678193c",
"metadata": {}, "metadata": {},
"outputs": [ "outputs": [

3
RL Testing/lev_pod_env.py
View File

@@ -155,6 +155,7 @@ class LevPodEnv(gym.Env):
return obs, {} return obs, {}
# The following was generated by AI - see [14]
def step(self, action): def step(self, action):
# Check if PyBullet connection is still active (GUI might be closed) # Check if PyBullet connection is still active (GUI might be closed)
try: try:
@@ -374,6 +375,7 @@ class LevPodEnv(gym.Env):
return obs, reward, terminated, truncated, info return obs, reward, terminated, truncated, info
# The following was generated by AI - see [15]
def _get_obs(self, initial_reset=False): def _get_obs(self, initial_reset=False):
""" """
Returns observation: [gaps(4), velocities(4)] Returns observation: [gaps(4), velocities(4)]
@@ -433,6 +435,7 @@ class LevPodEnv(gym.Env):
return obs return obs
# The following was generated by AI - see [16]
def _create_modified_urdf(self): def _create_modified_urdf(self):
""" """
Create a modified URDF with bolt positions adjusted based on initial gap height. Create a modified URDF with bolt positions adjusted based on initial gap height.

1
RL Testing/mag_lev_coil.py
View File

@@ -1,3 +1,4 @@
# The following was generated by AI - see [17]
class MagLevCoil: class MagLevCoil:
def __init__(self, r_resistance, l_inductance, source_voltage, maxCurrent): def __init__(self, r_resistance, l_inductance, source_voltage, maxCurrent):
self.R = r_resistance self.R = r_resistance

1
RL Testing/test_env.py
View File

@@ -1,3 +1,4 @@
# The following was generated by AI - see [18]
""" """
Test script for LevPodEnv Test script for LevPodEnv
Runs a simple episode with constant actions to verify the environment works Runs a simple episode with constant actions to verify the environment works

241
RL Testing/visualize_urdf.py
View File

@@ -1,241 +0,0 @@
"""
URDF Structure Visualizer for Lev Pod using PyBullet
Loads and displays the pod.urdf file in PyBullet's GUI
"""
import pybullet as p
import pybullet_data
import time
import os
# Initialize PyBullet in GUI mode
physicsClient = p.connect(p.GUI)
# Set up the simulation environment
p.setAdditionalSearchPath(pybullet_data.getDataPath())
p.setGravity(0, 0, 0)
# Configure camera view - looking at inverted maglev system
p.resetDebugVisualizerCamera(
cameraDistance=0.5,
cameraYaw=45,
cameraPitch=1, # Look up at the hanging pod
cameraTargetPosition=[0, 0, 0]
)
# Create the maglev track with collision physics (ABOVE, like a monorail)
# The track BOTTOM surface is at Z=0, pod hangs below
track_collision = p.createCollisionShape(
shapeType=p.GEOM_BOX,
halfExtents=[1.0, 0.2, 0.010] # 2m long × 0.4m wide × 2cm thick
)
track_visual = p.createVisualShape(
shapeType=p.GEOM_BOX,
halfExtents=[1.0, 0.2, 0.010],
rgbaColor=[0.3, 0.3, 0.3, 0.8] # Gray, semi-transparent
)
trackId = p.createMultiBody(
baseMass=0, # Static object
baseCollisionShapeIndex=track_collision,
baseVisualShapeIndex=track_visual,
basePosition=[0, 0, 0.010] # Track bottom at Z=0, center at Z=10mm
)
# Set track surface properties (steel)
p.changeDynamics(trackId, -1,
lateralFriction=0.3, # Steel-on-steel
restitution=0.1) # Minimal bounce
# Load the lev pod URDF
urdf_path = "pod.xml"
if not os.path.exists(urdf_path):
print(f"Error: Could not find {urdf_path}")
print(f"Current directory: {os.getcwd()}")
p.disconnect()
exit(1)
# INVERTED SYSTEM: Pod hangs BELOW track
# Track bottom is at Z=0
# Gap height = distance from track bottom DOWN to magnetic yoke (top of pod body)
# URDF collision spheres at +25mm from center = bolts that contact track from below
# URDF box top at +25mm from center = magnetic yoke
#
# For 10mm gap:
# - Yoke (top of pod at +25mm from center) should be at Z = 0 - 10mm = -10mm
# - Therefore: pod center = -10mm - 25mm = -35mm
# - Bolts (at +25mm from center) end up at: -35mm + 25mm = -10mm (touching track at Z=0)
start_pos = [0, 0, -0.10085] # Pod center 100.85mm BELOW track → yoke at 10mm gap
start_orientation = p.getQuaternionFromEuler([0, 0, 0])
podId = p.loadURDF(urdf_path, start_pos, start_orientation)
print("\n" + "=" * 60)
print("PyBullet URDF Visualizer")
print("=" * 60)
print(f"Loaded: {urdf_path}")
print(f"Position: {start_pos}")
print("\nControls:")
print(" • Mouse: Rotate view (left drag), Pan (right drag), Zoom (scroll)")
print(" • Ctrl+Mouse: Apply forces to the pod")
print(" • Press ESC or close window to exit")
print("=" * 60)
# Get and display URDF information
num_joints = p.getNumJoints(podId)
print(f"\nNumber of joints: {num_joints}")
# Get base info
base_mass, base_lateral_friction = p.getDynamicsInfo(podId, -1)[:2]
print(f"\nBase link:")
print(f" Mass: {base_mass} kg")
print(f" Lateral friction: {base_lateral_friction}")
# Get collision shape info
collision_shapes = p.getCollisionShapeData(podId, -1)
print(f"\nCollision shapes: {len(collision_shapes)}")
for i, shape in enumerate(collision_shapes):
shape_type = shape[2]
dimensions = shape[3]
local_pos = shape[5]
shape_names = {
p.GEOM_BOX: "Box",
p.GEOM_SPHERE: "Sphere",
p.GEOM_CAPSULE: "Capsule",
p.GEOM_CYLINDER: "Cylinder",
p.GEOM_MESH: "Mesh"
}
print(f" Shape {i}: {shape_names.get(shape_type, 'Unknown')}")
print(f" Dimensions: {dimensions}")
print(f" Position: {local_pos}")
# Enable visualization of collision shapes
p.configureDebugVisualizer(p.COV_ENABLE_RENDERING, 1)
p.configureDebugVisualizer(p.COV_ENABLE_GUI, 1)
p.configureDebugVisualizer(p.COV_ENABLE_WIREFRAME, 0)
# Add coordinate frame visualization at the pod's origin
axis_length = 0.1
p.addUserDebugLine([0, 0, 0], [axis_length, 0, 0], [1, 0, 0], lineWidth=3, parentObjectUniqueId=podId, parentLinkIndex=-1)
p.addUserDebugLine([0, 0, 0], [0, axis_length, 0], [0, 1, 0], lineWidth=3, parentObjectUniqueId=podId, parentLinkIndex=-1)
p.addUserDebugLine([0, 0, 0], [0, 0, axis_length], [0, 0, 1], lineWidth=3, parentObjectUniqueId=podId, parentLinkIndex=-1)
print("\n" + "=" * 60)
print("Visualization is running. Interact with the viewer...")
print("Close the PyBullet window to exit.")
print("=" * 60 + "\n")
# Store collision shape local positions and types for tracking
collision_local_positions = []
collision_types = []
for shape in collision_shapes:
collision_local_positions.append(shape[5]) # Local position (x, y, z)
collision_types.append(shape[2]) # Shape type (BOX, CYLINDER, etc.)
# Identify the 4 yoke top collision boxes and 4 sensor cylinders
# Both are at Z ≈ 0.08585m, but yokes are BOXES and sensors are CYLINDERS
yoke_indices = []
yoke_labels = []
sensor_indices = []
sensor_labels = []
for i, (local_pos, shape_type) in enumerate(zip(collision_local_positions, collision_types)):
# Check if at sensor/yoke height (Z ≈ 0.08585m)
if abs(local_pos[2] - 0.08585) < 0.001: # Within 1mm tolerance
if shape_type == p.GEOM_BOX:
# Yoke collision boxes (BOX shapes)
yoke_indices.append(i)
x_pos = "Front" if local_pos[0] > 0 else "Back"
y_pos = "Right" if local_pos[1] > 0 else "Left"
yoke_labels.append(f"{x_pos} {y_pos}")
elif shape_type == p.GEOM_CYLINDER or shape_type == p.GEOM_MESH:
# Sensor cylinders (may be loaded as MESH by PyBullet)
# Distinguish from yokes by X or Y position patterns
# Yokes have both X and Y non-zero, sensors have one coordinate near zero
if abs(local_pos[0]) < 0.06 or abs(local_pos[1]) < 0.02:
sensor_indices.append(i)
# Label sensors by position
if abs(local_pos[0]) < 0.001: # X ≈ 0 (center sensors)
label = "Center Right" if local_pos[1] > 0 else "Center Left"
else:
label = "Front" if local_pos[0] > 0 else "Back"
sensor_labels.append(label)
print(f"\nIdentified {len(yoke_indices)} yoke collision boxes for gap height tracking")
print(f"Identified {len(sensor_indices)} sensor cylinders for gap height tracking")
# Run the simulation with position tracking
import numpy as np
step_count = 0
print_interval = 240 # Print every second (at 240 Hz)
try:
while p.isConnected():
p.stepSimulation()
# Extract positions periodically
if step_count % print_interval == 0:
# Get pod base position and orientation
pos, orn = p.getBasePositionAndOrientation(podId)
# Convert quaternion to rotation matrix
rot_matrix = p.getMatrixFromQuaternion(orn)
rot_matrix = np.array(rot_matrix).reshape(3, 3)
# Calculate world positions of yoke tops and sensors
print(f"\n--- Time: {step_count/240:.2f}s ---")
print(f"Pod center: [{pos[0]*1000:.1f}, {pos[1]*1000:.1f}, {pos[2]*1000:.1f}] mm")
print("\nYoke Gap Heights:")
yoke_gap_heights = []
for i, yoke_idx in enumerate(yoke_indices):
local_pos = collision_local_positions[yoke_idx]
# Transform local position to world coordinates
local_vec = np.array(local_pos)
world_offset = rot_matrix @ local_vec
world_pos = np.array(pos) + world_offset
# Add 0.005m (5mm) to get top surface of yoke box (half-height of 10mm box)
yoke_top_z = world_pos[2] + 0.005
# Gap height: distance from track bottom (Z=0) down to yoke top
gap_height = 0.0 - yoke_top_z # Negative means below track
gap_height_mm = gap_height * 1000
yoke_gap_heights.append(gap_height_mm)
print(f" {yoke_labels[i]} yoke: pos=[{world_pos[0]*1000:.1f}, {world_pos[1]*1000:.1f}, {yoke_top_z*1000:.1f}] mm | Gap: {gap_height_mm:.2f} mm")
# Calculate average yoke gap height (for Ansys model input)
avg_yoke_gap = np.mean(yoke_gap_heights)
print(f" Average yoke gap: {avg_yoke_gap:.2f} mm")
print("\nSensor Gap Heights:")
sensor_gap_heights = []
for i, sensor_idx in enumerate(sensor_indices):
local_pos = collision_local_positions[sensor_idx]
# Transform local position to world coordinates
local_vec = np.array(local_pos)
world_offset = rot_matrix @ local_vec
world_pos = np.array(pos) + world_offset
# Add 0.005m (5mm) to get top surface of cylinder (half-length of 10mm cylinder)
sensor_top_z = world_pos[2] + 0.005
# Gap height: distance from track bottom (Z=0) down to sensor top
gap_height = 0.0 - sensor_top_z
gap_height_mm = gap_height * 1000
sensor_gap_heights.append(gap_height_mm)
print(f" {sensor_labels[i]} sensor: pos=[{world_pos[0]*1000:.1f}, {world_pos[1]*1000:.1f}, {sensor_top_z*1000:.1f}] mm | Gap: {gap_height_mm:.2f} mm")
# Calculate average sensor gap height
avg_sensor_gap = np.mean(sensor_gap_heights)
print(f" Average sensor gap: {avg_sensor_gap:.2f} mm")
step_count += 1
time.sleep(1./240.) # Run at 240 Hz
except KeyboardInterrupt:
print("\nExiting...")
p.disconnect()
print("PyBullet session closed.")