Planetary gearbox dimension generator final

README.md

@@ -0,0 +1,15 @@
+# PlanetaryGearbox
+
+To run the planetary gearbox optimization file, the required files to download are:
+
+constraint.m
+evalFitness.m
+explainX.m
+fval.m
+
+GA550.m
+goptions.m
+runAllOptimizations.m
+runOneOptimization.m
+
+When all the files are in the same folder in MATLAB, open runAllOptimizations, set the desired minimum constraints for the sun and ring gears, and run the script.

constraint.m

@@ -1,57 +1,84 @@
-function val = constraint(x, bearing_d, sun_extrusion_d)
-%% Planetary Gearbox GA Constraint Function
-% Elliot Stockwell 1/12/21
+function val = constraint(x, bearing_diameter, planet_shaft_diameter, sun_shaft_diameter, maxX)
+%% Gear Ratio Constraint Function
 %
-% x : array of input variables (planet_d1, planet_d2, sun_d1, dp)
-% bearing_d : inner diameter of the bearing
+% X is formatted
+% [ 
+%   front_ring_teeth / 3
+%   front_planet_s2_teeth
+%   front_sun_teeth / 3
+%   rear_planet_s2_teeth
+%   rear_sun_teeth / 3
+%   diametral_pitch
+% ]
 
-pm = 1E2; % penalty multiplier
+%% GA Constants
+PENALTY = 1e2;
+val = 0;
 
-planet_d1 = x(1);
-planet_d2 = x(2);
-sun_d1 = x(3);
-dp = 12 + 4 * x(4);
-ring_d2 = planet_d1 + planet_d2 + sun_d1;
+%% Interpret x
+front_ring_teeth = 3 * x(1);
+front_planet_s2_teeth = x(2);
+front_sun_teeth = 3 * x(3);
+rear_planet_s2_teeth = x(4);
+rear_sun_teeth = 3 * x(5);
+dp = 12 + 4 * x(6);
 
-g = zeros(1, 8);
 
-mod_tolerance = .01;
+%% Calculate derived variables
+rear_ring_teeth = front_ring_teeth;
 
-%% Make sure teeth are integers
-planet_1_mod = mod(dp * planet_d1, 1);
-planet_2_mod = mod(dp * planet_d2, 1);
-sun_mod = mod(dp * sun_d1, 3);
-ring_mod = mod(dp * ring_d2, 3);
+front_planet_s1_teeth = front_ring_teeth - ...
+    front_planet_s2_teeth - front_sun_teeth;
 
-if sun_mod < 3 - mod_tolerance && sun_mod > mod_tolerance
-    g(1) = pm;
-end
+rear_planet_s1_teeth = rear_ring_teeth - ...
+    rear_planet_s2_teeth - rear_sun_teeth;
 
-if ring_mod < 3 - mod_tolerance && ring_mod > mod_tolerance
-    g(2) = pm;
-end
 
-if planet_1_mod < 1 - mod_tolerance && planet_1_mod > mod_tolerance
-    g(3) = pm;
+%% Calculate diameters
+front_ring_diameter = front_ring_teeth / dp;
+front_planet_s1_diameter = front_planet_s1_teeth / dp;
+front_planet_s2_diameter = front_planet_s2_teeth / dp;
+front_sun_diameter = front_sun_teeth / dp;
+
+rear_ring_diameter = rear_ring_teeth / dp;
+rear_planet_s1_diameter = rear_planet_s1_teeth / dp;
+rear_planet_s2_diameter = rear_planet_s2_teeth / dp;
+rear_sun_diameter = rear_sun_teeth / dp;
+
+
+%% Calculate side constraints
+% Variables
+comp_planet_diameter = (front_sun_diameter + front_planet_s1_diameter) / (rear_sun_diameter + rear_planet_s1_diameter) - 1;
+
+% Don't exceed bearing bounds
+if 2 * front_planet_s1_diameter + front_sun_diameter > bearing_diameter
+    val = val + PENALTY;
 end
 
-if planet_2_mod < 1 - mod_tolerance && planet_2_mod > mod_tolerance
-    g(4) = pm;
+if 2 * rear_planet_s1_diameter + rear_sun_diameter > bearing_diameter
+    val = val + PENALTY;
 end
 
-%% Side constraints
-if dp > 24
-    g(5) = pm;
+if comp_planet_diameter < -0.01 || comp_planet_diameter > 0.01
+    val = val + PENALTY;
 end
 
-g(6) = max(0, (sun_d1 + 2 * planet_d1) / bearing_d - 1) * pm;
+% Make sure every gear is larger than the shaft that runs through it
+check_planet_arr = [...
+    front_planet_s1_diameter, ...
+    front_planet_s2_diameter, ...
+    rear_planet_s1_diameter, ...
+    rear_planet_s2_diameter];
 
-g(7) = max(0, (planet_d2 / (0.45 * ring_d2)) - 1) * pm;
+check_sun_arr = [...
+    front_sun_diameter, ...
+    rear_sun_diameter];
 
-if sun_d1 < sun_extrusion_d
-    g(8) = pm;
-end
-%% Sum constriants
-val = sum(g);
+% Count occurances where the constriant is violated
+val = val + PENALTY * sum(check_planet_arr < planet_shaft_diameter);
+
+val = val + PENALTY * sum(check_sun_arr < sun_shaft_diameter);
 
+% Enforce upper bounds on x
+val = val + PENALTY * sum(x > maxX);
 

evalFitness.m

@@ -0,0 +1,9 @@
+function fitnessVal = evalFitness(x, bearing_diameter, ...
+    planet_shaft_diameter, sun_shaft_diameter, maxX, GR_front_target, ...
+    GR_rear_target)
+
+phi = fval(x, GR_front_target, GR_rear_target);
+penalty = constraint(x, bearing_diameter, planet_shaft_diameter, ...
+    sun_shaft_diameter, maxX);
+
+fitnessVal = phi + penalty;

explainX.m

@@ -0,0 +1,66 @@
+function [front_diameters, front_teeth, ...
+    rear_diameters, rear_teeth, dp, frontGR, rearGR] = explainX(x)
+
+%% This function explains the design represented by x
+
+% X is formatted
+% [ 
+%   front_ring_teeth / 3
+%   front_planet_s2_teeth
+%   front_sun_teeth / 3
+%   rear_planet_s2_teeth
+%   rear_sun_teeth / 3
+%   diametral_pitch
+% ]
+
+%% Interpret x
+front_ring_teeth = 3 * x(1);
+front_planet_s2_teeth = x(2);
+front_sun_teeth = 3 * x(3);
+rear_planet_s2_teeth = x(4);
+rear_sun_teeth = 3 * x(5);
+dp = 12 + 4 * x(6);
+
+
+%% Calculate derived variables
+rear_ring_teeth = front_ring_teeth;
+
+front_planet_s1_teeth = front_ring_teeth - ...
+    front_planet_s2_teeth - front_sun_teeth;
+
+rear_planet_s1_teeth = rear_ring_teeth - ...
+    rear_planet_s2_teeth - rear_sun_teeth;
+
+
+%% Calculate diameters
+front_ring_diameter = front_ring_teeth / dp;
+front_planet_s1_diameter = front_planet_s1_teeth / dp;
+front_planet_s2_diameter = front_planet_s2_teeth / dp;
+front_sun_diameter = front_sun_teeth / dp;
+
+rear_ring_diameter = rear_ring_teeth / dp;
+rear_planet_s1_diameter = rear_planet_s1_teeth / dp;
+rear_planet_s2_diameter = rear_planet_s2_teeth / dp;
+rear_sun_diameter = rear_sun_teeth / dp;
+
+
+%% Prepare outputs
+front_diameters = [front_planet_s1_diameter, front_planet_s2_diameter, ...
+    front_sun_diameter, front_ring_diameter];
+
+rear_diameters = [rear_planet_s1_diameter, rear_planet_s2_diameter, ...
+    rear_sun_diameter, rear_ring_diameter];
+
+front_teeth = [front_planet_s1_teeth, front_planet_s2_teeth, ...
+    front_sun_teeth, front_ring_teeth];
+
+rear_teeth = [rear_planet_s1_teeth, rear_planet_s2_teeth, ...
+    rear_sun_teeth, rear_ring_teeth];
+
+frontGR = (front_planet_s1_teeth * front_ring_teeth) / ...
+    (front_sun_teeth * front_planet_s2_teeth) + 1;
+
+rearGR = (rear_planet_s1_teeth * rear_ring_teeth) / ...
+    (rear_sun_teeth * rear_planet_s2_teeth) + 1;
+
+

fval.m

@@ -1,13 +1,45 @@
-function val = fval(x, GR)
-%% Planetary Gear Ratio Target Function
-% Elliot Stockwell 1/12/21
+
+function val = fval(x, GR_front_target, GR_rear_target)
+%% Gear Ratio Objective Function
 %
-% x : array of input variables (planet_d1, planet_d2, sun_d1, dp)
-% GR : target gear ratio
+% X is formatted
+% [ 
+%   front_ring_teeth / 3
+%   front_planet_s2_teeth
+%   front_sun_teeth / 3
+%   rear_planet_s2_teeth
+%   rear_sun_teeth / 3
+%   diametral_pitch
+% ]
+
+
+%% Interpret x
+front_ring_teeth = 3 * x(1);
+front_planet_s2_teeth = x(2);
+front_sun_teeth = 3 * x(3);
+rear_planet_s2_teeth = x(4);
+rear_sun_teeth = 3 * x(5);
+dp = 12 + 4 * x(6);
+
+
+%% Calculate derived variables
+rear_ring_teeth = front_ring_teeth;
+
+front_planet_s1_teeth = front_ring_teeth - ...
+    front_planet_s2_teeth - front_sun_teeth;
+
+rear_planet_s1_teeth = rear_ring_teeth - ...
+    rear_planet_s2_teeth - rear_sun_teeth;
+
+%% Calculate gear ratios
+front_gr = (front_planet_s1_teeth * front_ring_teeth) / ...
+    (front_sun_teeth * front_planet_s2_teeth) + 1;
+
+rear_gr = (rear_planet_s1_teeth * rear_ring_teeth) / ...
+    (rear_sun_teeth * rear_planet_s2_teeth) + 1;
+
 
-planet_d1 = x(1);
-planet_d2 = x(2);
-sun_d1 = x(3);
-ring_d2 = planet_d1 + planet_d2 + sun_d1;
+%% Calculate objective
+val = (front_gr - GR_front_target) ^ 2;
+val = val + (rear_gr - GR_rear_target) ^ 2;
 
-val = (((planet_d1 * ring_d2 / sun_d1 / planet_d2) + 1) - GR)^2;

runAllOptimizations.m

@@ -0,0 +1,120 @@
+%% Master running file for GA methods
+
+%% Define constants
+bearing_diameter = 4.50;        % in
+sun_shaft_diameter = 1.23;      % in
+planet_shaft_diameter = 0.85;   % in
+
+GR_front_target = 6.5;      % dimensionless
+GR_rear_target = 8.00;      % dimensionless
+
+numGenerations = 200;
+num_runs = 500;
+
+%% Define bounds
+%
+% X is formatted
+% [ 
+%   front_ring_teeth / 3
+%   front_planet_s2_teeth
+%   front_sun_teeth / 3
+%   rear_planet_s2_teeth
+%   rear_sun_teeth / 3
+%   diametral_pitch
+% ]
+
+% Remember some of these are / 3
+min_x = [...
+    14, ...
+    14, ...
+    5, ...
+    14, ...
+    5, ...
+    1];
+
+max_x = [...
+    35, ...
+    48, ...
+    35, ...
+    48, ...
+    35, ...
+    4];
+
+
+%% Perform last calculations
+vlb = min_x;
+
+vdiff = max_x - min_x;
+bits = ceil(log(vdiff) / log(2));
+
+vub = vlb + 2 .^ bits - 1;
+maxX = max_x;
+
+
+%% Run GA
+data = zeros(num_runs, 20);
+
+fprintf('\nBeginning GA\n');
+
+parfor i=1:num_runs
+    [FOPT, X] = runOneOptimization(vlb, vub, bits, bearing_diameter, ...
+        planet_shaft_diameter, sun_shaft_diameter, maxX, ...
+        GR_front_target, GR_rear_target, numGenerations)
+
+    [front_diameters, front_teeth, rear_diameters, rear_teeth, ...
+        dp, frontGR, rearGR] = explainX(X);
+
+    data(i, :) = [FOPT, front_diameters, front_teeth, ...
+        rear_diameters, rear_teeth, dp, frontGR, rearGR];
+end
+
+
+%% Save Info
+filename = 'ga_results.csv';
+
+header = {'Badness', ...
+    'Front Planet S1 Diameter (in)', 'Front Planet S2 Diameter (in)', ...
+    'Front Sun Diameter (in)', 'Front Ring Diameter (in)', ...
+    'Front Planet S1 Teeth', 'Front Planet S2 Teeth', ...
+    'Front Sun Teeth', 'Front Ring Teeth', ...
+    'Rear Planet S1 Diameter (in)', 'Rear Planet S2 Diameter (in)', ...
+    'Rear Sun Diameter (in)', 'Rear Ring Diameter (in)', ...
+    'Rear Planet S1 Teeth', 'Rear Planet S2 Teeth', ...
+    'Rear Sun Teeth', 'Rear Ring Teeth', ...
+    'Diametral Pitch (1/in)', 'Front GR', 'Rear GR'};
+
+commaHeader = [header;repmat({','},1,numel(header))];
+commaHeader = commaHeader(:)';
+textHeader = cell2mat(commaHeader);
+
+fid = fopen(filename, 'w');
+fprintf(fid, '%s\n', textHeader(1:end-1));
+fclose(fid);
+
+data_sorted = sortrows(data);
+dlmwrite(filename, data_sorted, '-append');
+
+
+%% Clean up
+fprintf('\nDone.\n');
+
+%% Create histogram
+
+feasible_region_front = [6.25 6.25 6.75 6.75 6.25];
+feasible_region_rear = [7.875 8.5 8.5 7.875 7.875];
+
+figure(1);
+hold on;grid on;grid minor;
+colormap jet;
+set(gcf, 'Position', [100 100 900 700]);
+sct = scatter(data(:, end-1), data(:, end), 25, data(:, 1));
+trg = plot(GR_front_target, GR_rear_target);
+region = plot(feasible_region_front, feasible_region_rear);
+set(trg, 'Marker', '.', 'MarkerSize', 25, 'Color', 'r', ...
+    'LineStyle', 'None');
+set(region, 'Color', 'k', 'LineStyle', '-.', 'LineWidth', 1);
+title('Distribution of Results');
+xlabel('Front Gear Ratio');ylabel('Rear Gear Ratio');
+legend([sct, trg], 'GA Results', 'Target');
+
+