From 2a8d918dc28e132852c45ba96c93c6abda9121ee Mon Sep 17 00:00:00 2001 From: Elliot Stockwell Date: Tue, 16 Feb 2021 17:04:00 -0500 Subject: [PATCH 1/9] Create runOneOptimization.m --- runOneOptimization.m | 20 ++++++++++++++++++++ 1 file changed, 20 insertions(+) create mode 100644 runOneOptimization.m diff --git a/runOneOptimization.m b/runOneOptimization.m new file mode 100644 index 0000000..7645ebf --- /dev/null +++ b/runOneOptimization.m @@ -0,0 +1,20 @@ +function [FOPT, X] = runOneOptimization(vlb, vub, bits, bearing_diameter, ... + planet_shaft_diameter, sun_shaft_diameter, maxX, GR_front_target, ... + GR_rear_target, numGenerations) + +% Chromosome length +l = sum(bits); + +%% Define GA options +parain = zeros(14, 1); +parain(1) = 0; +parain(11) = 4 * l; +parain(13) = (l + 1) / (2 * l * parain(11)); +parain(14) = numGenerations; + +options = goptions(parain); + +[X, FOPT, ~, ~, ~, ~, ~] = GA550('evalFitness', ... + [], options, vlb, vub, bits, bearing_diameter, ... + planet_shaft_diameter, sun_shaft_diameter, maxX, GR_front_target, ... + GR_rear_target); \ No newline at end of file From 2b91aad2a9a982d85a28bbfc11d813553470aa7c Mon Sep 17 00:00:00 2001 From: Elliot Stockwell Date: Tue, 16 Feb 2021 17:04:03 -0500 Subject: [PATCH 2/9] Create runAllOptimizations.m --- runAllOptimizations.m | 120 ++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 120 insertions(+) create mode 100644 runAllOptimizations.m diff --git a/runAllOptimizations.m b/runAllOptimizations.m new file mode 100644 index 0000000..3f741c1 --- /dev/null +++ b/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'); + + From 011b4f732035e7da9d68ff5773633fa4c4fe56c4 Mon Sep 17 00:00:00 2001 From: Elliot Stockwell Date: Tue, 16 Feb 2021 17:04:05 -0500 Subject: [PATCH 3/9] Create goptions.m --- goptions.m | 80 ++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 80 insertions(+) create mode 100644 goptions.m diff --git a/goptions.m b/goptions.m new file mode 100644 index 0000000..c0f55c6 --- /dev/null +++ b/goptions.m @@ -0,0 +1,80 @@ +function OPTIONS=goptions(parain); +%GOPTIONS Default parameters used by the genetic algorithm GENETIC. +% +% Note that since the original version was written, the Matlab Optimization +% Toolbox now uses "optimset" to set generic optimization parameters, so +% this format is somewhat outdated. +% +% The genetic algorithm parameters used for this implementation are: +% +% OPTIONS(1)-Display flag: 0 = none, 1 = some, 2 = all (Default: 1). +% OPTIONS(2)-Termination bit string affinity value (Default: 0.90; set to zero to turn off) +% OPTIONS(3)-Termination tolerance for fitness (Default: 0; not normally used). +% OPTIONS(4)-Termination number of consecutive generations with same best +% fitness (Default: 0; to use, set number, be sure OPTIONS(2) and OPTIONS(3) = 0). +% OPTIONS(5)-Number of elite individuals (Default: 0; no elitism). +% OPTIONS(6)- +% OPTIONS(7)- +% OPTIONS(8)- +% OPTIONS(9)- +% OPTIONS(10)- +% Genetic Algorithm-specific inputs +% OPTIONS(11)-Population size (fixed) +% OPTIONS(12)-Probability of crossover +% OPTIONS(13)-Probability of mutation +% OPTIONS(14)-Maximum number of generations, always used as safeguard +% (Default: 200). +% +% +% Explanation of defaults: +% The default algorithm displays statistical information for each +% generation by setting OPTIONS(1) = 1. Plots are produced when +% OPTIONS(1) = 2. +% The OPTIONS(2) flag is originally set for termination criterion based +% on X; here it is used if bit string affinity is selected. +% The default fitness function termination tolerance, +% OPTIONS(3), is set to 0, which terminates the optimization when 5 +% consecutive best generation fitness values are the same. A positive +% value terminates the optimization when the normalized difference +% between the previous fitness and current generation fitness is less +% than the tolerance. See the code for details. +% OPTIONS(4) has a default value of 5; this means if the best fitness +% value in the population is unchanged for 5 consecutive generations +% the GA is terminated. +% The default algorithm uses a fixed population size, OPTIONS(11), +% and no generational overlap. The default population size is 30. +% Three genetic operations: selection, crossover, and mutation are +% used for procreation. +% The default selection scheme is tournament selection. +% Crossover occurs with probability Pc=OPTIONS(12). The default +% crossover scheme is uniform crossover with Pc = 0.5. +% Each allele of the offspring mutates independently with probability +% Pm=OPTIONS(13); here the default is 0.01. +% The default number of maximum generations, OPTIONS(14) is 200. +% +% Last modified by Bill Crossley 07/20/04 + +% The following lines have been commented out by Steven Lamberson. +% They have been changed to what is seen below them. (06/30/06). +% This change was made in order to fix the following problems: +% 1 - code changed user supplied options(1)=0 to options(1)=1 +% 2 - code changed user supplied options(2)=0 to options(2)=0.9 + +%if nargin<1; parain = []; end +%sizep=length(parain); +%OPTIONS=zeros(1,14); +%OPTIONS(1:sizep)=parain(1:sizep); +%default_options=[1,0.9,0,0,0,0,0,0,0,0,0,0,0,200]; +%OPTIONS=OPTIONS+(OPTIONS==0).*default_options + +if nargin<1; parain = []; end +sizep=length(parain); +OPTIONS=zeros(1,14)-1; +OPTIONS(1:sizep)=parain(1:sizep); +default_options=[1,0.9,0,0,0,0,0,0,0,0,0,0,0,200]; +for i = 1:length(OPTIONS) + if OPTIONS(i) == -1 + OPTIONS(i) = default_options(i); + end +end + From 66aba830dfbab30532509da6748678fec6583ec2 Mon Sep 17 00:00:00 2001 From: Elliot Stockwell Date: Tue, 16 Feb 2021 17:04:09 -0500 Subject: [PATCH 4/9] Create GA550.m --- GA550.m | 555 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 555 insertions(+) create mode 100644 GA550.m diff --git a/GA550.m b/GA550.m new file mode 100644 index 0000000..edd4e9c --- /dev/null +++ b/GA550.m @@ -0,0 +1,555 @@ +function [xopt,fopt,stats,nfit,fgen,lgen,lfit] = GA550(fun, ... + x0,options,vlb,vub,bits,P1,P2,P3,P4,P5,P6,P7P,P8,P9,P10) +%GA550 minimizes a fitness function using a simple genetic algorithm. +% +% X=GA550('FUN',X0,OPTIONS,VLB,VUB) uses a simple +% genetic algorithm to find a minimum of the fitness function +% FUN. FUN can be a user-defined M-file: FUN.M, or it can be a +% string containing the function itself. The user may define all +% or part of an initial population X0. Any undefined individuals +% will be randomly generated between the lower and upper bounds +% (VLB and VUB). If X0 is an empty matrix, the entire initial +% population will be randomly generated. Use OPTIONS to specify +% flags, tolerances, and input parameters. Type HELP GOPTIONS +% for more information and default values. +% +% X=GA550('FUN',X0,OPTIONS,VLB,VUB,BITS) allows the user to +% define the number of BITS used to code non-binary parameters +% as binary strings. Note: length(BITS) must equal length(VLB) +% and length(VUB). If BITS is not specified, as in the previous +% call, the algorithm assumes that the fitness function is +% operating on a binary population. +% +% X=GA550('FUN',X0,OPTIONS,VLB,VUB,BITS,P1,P2,...) allows up +% to ten arguments, P1,P2,... to be passed directly to FUN. +% F=FUN(X,P1,P2,...). If P1,P2,... are not defined, F=FUN(X). +% +% [X,FOPT,STATS,NFIT,FGEN,LGEN,LFIT]=GA550() +% X - design variables of best ever individual +% FOPT - fitness value of best ever individual +% STATS - [min mean max stopping_criterion] fitness values +% for each generation +% NFIT - number of fitness function evalations +% FGEN - first generation population +% LGEN - last generation population +% LFIT - last generation fitness +% +% The algorithm implemented here is based on the book: Genetic +% Algorithms in Search, Optimization, and Machine Learning, +% David E. Goldberg, Addison-Wiley Publishing Company, Inc., +% 1989. +% +% Originally created on 1/10/93 by Andrew Potvin, Mathworks, Inc. +% Modified on 2/3/96 by Joel Grasmeyer. +% Modified on 11/12/02 by Bill Crossley. +% Modified on 7/20/04 by Bill Crossley. + +% Make best_feas global for stopping criteria (4/13/96) +global best_feas +global gen +global fit_hist +% Load input arguments and check for errors +if nargin<4, + error('No population bounds given.') +elseif (size(vlb,1)~=1) | (size(vub,1)~=1), + % Remark: this will change if algorithm accomodates matrix variables + error('VLB and VUB must be row vectors') +elseif (size(vlb,2)~=size(vub,2)), + error('VLB and VUB must have the same number of columns.') +elseif (size(vub,2)~=size(x0,2)) & (size(x0,1)>0), + error('X0 must all have the same number of columns as VLB and VUB.') +elseif any(vlb>vub), + error('Some lower bounds greater than upper bounds') +else + x0_row = size(x0,1); + for i=1:x0_row, + if any(x0(x0_row,:)vub), + error('Some initial population not within bounds.') + end % if initial pop not within bounds + end % for initial pop +end % if nargin<4 + +if nargin<6, + bits = []; +elseif (size(bits,1)~=1) | (size(bits,2)~=size(vlb,2)), + % Remark: this will change if algorithm accomodates matrix variables + error('BITS must have one row and length(VLB) columns') +elseif any(bits~=round(bits)) | any(bits<1), + error('BITS must be a vector of integers >0') +end % if nargin<6 + +% Form string to call for function evaluation +if ~( any(fun<48) | any(fun>122) | any((fun>90) & (fun<97)) | ... + any((fun>57) & (fun<65)) ), + % Only alphanumeric characters implies that 'fun' is a separate m-file + evalstr = [fun '(x']; + for i=1:nargin-6, + evalstr = [evalstr,',P',int2str(i)]; + end +else + % Non-alphanumeric characters implies that the function is contained + % within the single quotes + evalstr = ['(',fun]; +end + +% Determine all options +% Remark: add another options index for type of termination criterion +if size(options,1)>1, + error('OPTIONS must be a row vector') +else + % Use default options for those that were not passed in + options = goptions(options); +end +PRINTING = options(1); +BSA = options(2); +fit_tol = options(3); +nsame = options(4)-1; +elite = options(5); + +% Since operators are tournament selection and uniform crossover and +% default coding is Gray / binary, set crossover rate to 0.50 and use +% population size and mutation rate based on Williams, E. A., and Crossley, +% W. A., "Empirically-derived population size and mutation rate guidelines +% for a genetic algorithm with uniform crossover," Soft Computing in +% Engineering Design and Manufacturing, 1998. If user has entered values +% for these options, then user input values are used. +if options(11) == 0, + pop_size = sum(bits) * 4; +else + pop_size = options(11); +end +if options(12) == 0, + Pc = 0.5; +else + Pc = options(12); +end +if options(13) == 0, + Pm = (sum(bits) + 1) / (2 * pop_size * sum(bits)); +else + Pm = options(13); +end +max_gen = options(14); +% Ensure valid options: e.q. Pc,Pm,pop_size,max_gen>0, Pc,Pm<1 +if any([Pc Pm pop_size max_gen]<0) | any([Pc Pm]>1), + error('Some Pc,Pm,pop_size,max_gen<0 or Pc,Pm>1') +end + +% Encode fitness (cost) function if necessary +ENCODED = any(any(([vlb; vub; x0]~=0) & ([vlb; vub; x0]~=1))) | .... + ~isempty(bits); +if ENCODED, + [fgen,lchrom] = encode(x0,vlb,vub,bits); +else + fgen = x0; + lchrom = size(vlb,2); +end + +% Display warning if initial population size is odd +if rem(pop_size,2)==1, + disp('Warning: Population size should be even. Adding 1 to population.') + pop_size = pop_size +1; +end + +% Form random initial population if not enough supplied by user +if size(fgen,1)=1, + if ENCODED, + disp('Variable coding as binary chromosomes successful.') + disp('') + fgen = decode(fgen,vlb,vub,bits); + end + disp(' Fitness statistics') + if nsame > 0 + disp('Generation Minimum Mean Maximum isame') + elseif BSA > 0 + disp('Generation Minimum Mean Maximum BSA') + else + disp('Generation Minimum Mean Maximum not used') + end +end + +% Set up main loop +STOP_FLAG = 0; +for generation = 1:max_gen+1, + old_gen = new_gen; + + % Decode binary strings if necessary + if ENCODED, + x_pop = decode(old_gen,vlb,vub,bits); + else + x_pop = old_gen; + end + + % Get fitness of each string in population + for i = 1:pop_size, + x = x_pop(i,:); + fitness(i) = eval([evalstr,')']); + nfit = nfit + 1; + end + + % Store minimum fitness value from previous generation (except for + % initial generation) + if generation > 1, + min_fit_prev = min_fit; + min_gen_prev = min_gen; + min_x_prev = min_x; + end + + % identify worst (maximum) fitness individual in current generation + [max_fit,max_index] = max(fitness); + + % impose elitism - currently only one individual; this replaces worst + % individual of current generation with best of previous generation + if (generation > 1 & elite > 0), + old_gen(max_index,:) = min_gen_prev; + x_pop(max_index,:) = min_x_prev; + fitness(max_index) = min_fit_prev; + end + + % identify best (minimum) fitness individual in current generation and + % store bit string and x values + [min_fit,min_index] = min(fitness); + min_gen = old_gen(min_index,:); + min_x = x_pop(min_index,:); + + % Store best fitness and x values + if min_fit < fopt, + fopt = min_fit; + xopt = min_x; + end + + % Compute values for isame or BSA_pop stopping criteria + if nsame > 0 + if generation > 1 + if min_fit_prev == min_fit + isame = isame + 1; + else + isame = 0; + end + end + elseif BSA > 0 + bitlocavg = mean(old_gen,1); + BSA_pop = 2 * mean(abs(bitlocavg - 0.5)); + end + + + % Calculate generation statistics + if nsame > 0 + stats = [stats; generation-1,min(fitness),mean(fitness), ... + max(fitness), isame]; + elseif BSA > 0 + stats = [stats; generation-1,min(fitness),mean(fitness), ... + max(fitness), BSA_pop]; + else + stats = [stats; generation-1,min(fitness),mean(fitness), ... + max(fitness), 0]; + end + + % Display if necessary + if PRINTING>=1, + disp([sprintf('%5.0f %12.6g %12.6g %12.6g %12.6g', stats(generation,1), ... + stats(generation,2),stats(generation,3), stats(generation,4),... + stats(generation,5))]); + end + + % Check for termination + % The default termination criterion is bit string affinity. Also + % available are fitness tolerance across five generations and number of + % consecutive generations with same best fitness. These can be used + % concurrently. + if fit_tol>0, % if fit_tol > 0, then fitness tolerance criterion used + if generation>5, + % Check for normalized difference in fitness minimums + if stats(generation,1) ~= 0, + if abs(stats(generation-5,1)-stats(generation,1))/ ... + stats(generation,1) < fit_tol + if PRINTING >= 1 + fprintf('\n') + disp('GA converged based on difference in fitness minimums.') + end + lfit = fitness; + if ENCODED, + lgen = x_pop; + else + lgen = old_gen; + end + return + end + else + if abs(stats(generation-5,1)-stats(generation,1)) < fit_tol + if PRINTING >= 1 + fprintf('\n') + disp('GA converged based on difference in fitness minimums.') + end + lfit = fitness; + if ENCODED, + lgen = x_pop; + else + lgen = old_gen; + end + return + end + end + end + elseif nsame > 0, % consecutive minimum fitness value criterion + if isame == nsame + if PRINTING >= 1 + fprintf('\n') + disp('GA stopped based on consecutive minimum fitness values.') + end + lfit = fitness; + if ENCODED, + lgen = x_pop; + else + lgen = old_gen; + end + return + end + elseif BSA > 0, % bit string affinity criterion + if BSA_pop >= BSA, + if PRINTING >=1 + fprintf('\n') + disp('GA stopped based on bit string affinity value.') + end + lfit = fitness; + if ENCODED, + lgen = x_pop; + else + lgen = old_gen; + end + return + end + end + + % Tournament selection + new_gen = tourney(old_gen,fitness); + + % Crossover + new_gen = uniformx(new_gen,Pc); + + % Mutation + new_gen = mutate(new_gen,Pm); + + % Always save last generation. This allows user to cancel and + % restart with x0 = lgen + if ENCODED, + lgen = x_pop; + else + lgen = old_gen; + end + + +end % for max_gen + +% Maximum number of generations reached without termination +lfit = fitness; +if PRINTING>=1, + fprintf('\n') + disp('Maximum number of generations reached without termination') + disp('criterion met. Either increase maximum generations') + disp('or ease termination criterion.') +end + + +% end genetic + +function [gen,lchrom,coarse,nround] = encode(x,vlb,vub,bits) +%ENCODE Converts from variable to binary representation. +% [GEN,LCHROM,COARSE,nround] = ENCODE(X,VLB,VUB,BITS) +% encodes non-binary variables of X to binary. The variables +% in the i'th column of X will be encoded by BITS(i) bits. VLB +% and VUB are the lower and upper bounds on X. GEN is the binary +% representation of these X. LCHROM=SUM(BITS) is the length of +% the binary chromosome. COARSE(i) is the coarseness of the +% i'th variable as determined by the variable ranges and +% BITS(i). ROUND contains the absolute indices of the +% X which where rounded due to finite BIT length. +% +% Copyright (c) 1993 by the MathWorks, Inc. +% Andrew Potvin 1-10-93. + +% Remark: what about handling case where length(bits)~=length(vlb)? + + +lchrom = sum(bits); +coarse = (vub-vlb)./((2.^bits)-1); +[x_row,x_col] = size(x); + +gen = []; +if ~isempty(x), + temp = (x-ones(x_row,1)*vlb)./ ... + (ones(x_row,1)*coarse); + b10 = round(temp); + % Since temp and b10 should contain integers 1e-4 is close enough + nround = find(b10-temp>1e-4); + gen = b10to2(b10,bits); +end + +% end encode + + +function [x,coarse] = decode(gen,vlb,vub,bits) +%DECODE Converts from binary Gray code to variable representation. +% [X,COARSE] = DECODE(GEN,VLB,VUB,BITS) converts the binary +% population GEN to variable representation. Each individual +% of GEN should have SUM(BITS). Each individual binary string +% encodes LENGTH(VLB)=LENGTH(VUB)=LENGTH(BITS) variables. +% COARSE is the coarseness of the binary mapping and is also +% of length LENGTH(VUB). +% +% this *.m file created by combining "decode.m" from the MathWorks, Inc. +% originally created by Andrew Potvin in 1993, with "GDECODE.FOR" written +% by William A. Crossley in 1996. +% +% William A. Crossley, Assoc. Prof. School of Aero. & Astro. +% Purdue University, 2001 +% +% gen is an array [population size , string length], each row is one individual's chromosome +% vlb is a row vector [number of parameters], each entry is the lower bound for a variable +% vub is a row vector [number of parameters], each entry is the upper bound for a variable +% bits is a row vector [number of parameters], each entry is number of bits used for a variable +% + +no_para = length(bits); % extract number of parameters using number of rows in bits vector +npop = size(gen,1); % extract population size using number of rows in gen array +x = zeros(npop, no_para); % sets up x as an array [population size, number of parameters] +coarse = zeros(1,no_para); % sets up coarse as a row vector [number of parameters] + +for J = 1:no_para, % extract the resolution of the parameters + coarse(J) = (vub(J)-vlb(J))/(2^bits(J)-1); % resolution of parameter J +end + +for K = 1:npop, % outer loop through each individual (there may be a more efficient way to operate on the + % gen array) BC 10/10/01 + sbit = 1; % initialize starting bit location for a parameter + ebit = 0; % initialize ending bit location + + for J = 1:no_para, % loop through each parameter in the problem + ebit = bits(J) + ebit; % pick the end bit for parameter J + accum = 0.0; % initialize the running sum for parameter J + ADD = 1; % add / subtract flag for Gray code; add if(ADD), subtract otherwise + for I = sbit:ebit, % loop through each bit in parameter J + pbit = I + 1 - sbit; % pbit determines value to be added or subtracted for Gray code + if (gen(K,I)) % if "1" is at current location + if (ADD) % add if appropriate + accum = accum + (2.0^(bits(J)-pbit+1) - 1.0); + ADD = 0; % next time subtract + else + accum = accum - (2.0^(bits(J)-pbit+1) - 1.0); + ADD = 1; % next time add + end + end + end % end of I loop through each bit + x(K,J) = accum * coarse(J) + vlb(J); % decoded parameter J for individual K + sbit = ebit + 1; % next parameter starting bit location + end % end of J loop through each parameter +end % end of K loop through each individual + +%end gdecode + + +function [new_gen,mutated] = mutate(old_gen,Pm) +%MUTATE Changes a gene of the OLD_GEN with probability Pm. +% [NEW_GEN,MUTATED] = MUTATE(OLD_GEN,Pm) performs random +% mutation on the population OLD_POP. Each gene of each +% individual of the population can mutate independently +% with probability Pm. Genes are assumed possess boolean +% alleles. MUTATED contains the indices of the mutated genes. +% +% Copyright (c) 1993 by the MathWorks, Inc. +% Andrew Potvin 1-10-93. + +mutated = find(rand(size(old_gen)) Date: Tue, 16 Feb 2021 17:04:11 -0500 Subject: [PATCH 5/9] Create fval.m --- fval.m | 43 +++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 43 insertions(+) create mode 100644 fval.m diff --git a/fval.m b/fval.m new file mode 100644 index 0000000..cfbc1e0 --- /dev/null +++ b/fval.m @@ -0,0 +1,43 @@ +function val = fval(x, GR_front_target, GR_rear_target) +%% Gear Ratio Objective Function +% +% 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; + + +%% Calculate objective +val = (front_gr - GR_front_target) ^ 2; +val = val + (rear_gr - GR_rear_target) ^ 2; \ No newline at end of file From 4877edf37ae8fa09c273deac17a6a4d34ab3222c Mon Sep 17 00:00:00 2001 From: Elliot Stockwell Date: Tue, 16 Feb 2021 17:04:14 -0500 Subject: [PATCH 6/9] Create explainX.m --- explainX.m | 66 ++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 66 insertions(+) create mode 100644 explainX.m diff --git a/explainX.m b/explainX.m new file mode 100644 index 0000000..083d296 --- /dev/null +++ b/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; + + From aa571c84a139e6f8090e6bce33809df4df9c60b8 Mon Sep 17 00:00:00 2001 From: Elliot Stockwell Date: Tue, 16 Feb 2021 17:04:16 -0500 Subject: [PATCH 7/9] Create evalFitness.m --- evalFitness.m | 9 +++++++++ 1 file changed, 9 insertions(+) create mode 100644 evalFitness.m diff --git a/evalFitness.m b/evalFitness.m new file mode 100644 index 0000000..99f0519 --- /dev/null +++ b/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; \ No newline at end of file From a97ef5d22d6caccb08278afd0b92aba0636c50bd Mon Sep 17 00:00:00 2001 From: Elliot Stockwell Date: Tue, 16 Feb 2021 17:04:18 -0500 Subject: [PATCH 8/9] Create constraint.m --- constraint.m | 83 ++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 83 insertions(+) create mode 100644 constraint.m diff --git a/constraint.m b/constraint.m new file mode 100644 index 0000000..a40dd3e --- /dev/null +++ b/constraint.m @@ -0,0 +1,83 @@ +function val = constraint(x, bearing_diameter, planet_shaft_diameter, sun_shaft_diameter, maxX) +%% Gear Ratio Constraint Function +% +% 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 +% ] + +%% GA Constants +PENALTY = 1e2; +val = 0; + +%% 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; + + +%% 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 2 * rear_planet_s1_diameter + rear_sun_diameter > bearing_diameter + val = val + PENALTY; +end + +if comp_planet_diameter < -0.01 || comp_planet_diameter > 0.01 + val = val + PENALTY; +end + +% 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]; + +check_sun_arr = [... + front_sun_diameter, ... + rear_sun_diameter]; + +% 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); From 4648bcecb23d4742c8c724b43c0553d4469b6c6d Mon Sep 17 00:00:00 2001 From: Elliot Stockwell <72889654+epstockwell@users.noreply.github.com> Date: Tue, 16 Feb 2021 17:08:38 -0500 Subject: [PATCH 9/9] Create README.md --- README.md | 15 +++++++++++++++ 1 file changed, 15 insertions(+) create mode 100644 README.md diff --git a/README.md b/README.md new file mode 100644 index 0000000..920dcac --- /dev/null +++ b/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.