%  EM121 Project calculations
%  Jerry Fine
%  May 14, 2011

clc
clear all

%  Parameters

rho_air = 0.0025;        % air density slugs/ft3
rho_al  = 5.2;           % aluminum density slugs/ft3

L_i = 17.643;            % Length blade inches
L_f = L_i/12;            % Length blade ft
t_i = 0.4;               % blade thickness inches
t_f = t_i/12;            % blade thickness ft
A_i = 0.3 * t_i;         % blade c.s in in2
A_f = A_i/144;           % blade c.s. in ft2
V_i = L_i * A_i;         % blade volume in in^3;
V_f = V_i/12^3;          % blade volume in ft^3;

m_b = V_f * rho_al      % blade mass in slugs
m_b_metric = m_b*14.59  % blade mass in kg

y_b = 3.937;             % pin offset in y - inches
R_i = sqrt(y_b^2 + (L_i/2)^2);  % c.g. R inches
R_f = R_i/12;            % R in ft

w = 1300*2*pi/60;        % angular vel. rad/s
C = m_b*R_f*w^2;         % centrifugal force in lbf
C_metric=C*4.448         % centrifugal force in N

CD = 1.0;                % drag coefficient
A = 1/6*CD*rho_air*t_f*w^2*L_f^3  % aero force in lb
A_metric = A*4.448       % aero force in N
dA = 3/4 * L_i;          % aero center inches

phi_r = atan(y_b/(L_i/2));   % angle phi in radians


% input of the coordinates of blade holes

dxb = 0.787;             % blade hole spacing inches
xi = 2.363;              % first hole location
xbs = xi +[0:5]*dxb;     % holes x-coords in blades inches
ybs = y_b*ones(1,6);     % holes y-coords in blade inches

% input of coordinates of disk holes

s = 0.394;               % spacing unit inches
xs = 2.361;              % starting x values in inches
ys = -0.241;             % starting y value in inches
xds = [ xs
        xs
        xs
        xs
        xs + s
        xs + s
        xs + s
        xs + 2*s
        xs + 2*s
        xs + 2*s
        xs + 3*s
        xs + 3*s
        xs + 4*s ];
 yds = [ ys
         ys + s
         ys + 2*s
         ys + 3*s
         ys + 0.5*s
         ys + 1.5*s
         ys + 2.5*s
         ys + s
         ys + 2*s
         ys + 3*s
         ys + 1.5*s
         ys + 2.5*s
         ys + 2*s ];
plot(xds,yds,'go',xbs,ybs,'ro')
grid on
axis([-1 10 -2  5])
hold on
plot(xds(9),yds(9),'bx')
hold off
     
% loop through the possibilities compute a table
% rows correspond to holes on hub
% cols correspond to holes on blade

fprintf('hole     b1      b2       b3        b4       b5      b6\n')

for ny = 1:13
    y2 =yds(ny,1);
    x2 =xds(ny,1);
    fprintf('d%2i',ny)
    for nx = 1:6
        y1 = ybs(1,nx);
        x1 = xbs(1,nx);
        dx = x1 - x2;
        dy = y1 - y2;
        LL = sqrt( dx^2+dy^2);      % Length of link 
        theta = atan2(dy,dx);
        F(ny,nx) = (L_i/2*C*sin(phi_r)+dA*A)/(x1*sin(theta));
        %F(ny,nx) = (L_i/2*C*sin(phi_r))/(x1*sin(theta));  % neglect air
        Am(ny,nx) = LL*F(ny,nx);   % L * F product
        fprintf('   %6.2f',F(ny,nx))  % change to plot Am
    end
    fprintf('\n')
end
        
% Comments

%  It appears to me that the "best" combination
%  which takes into account Length (weight) and
%  force is the extreme right hole on the hub
%  with the 5th hole out on the blade.  Force with
%  that combination is 70.5 lb.  Force with farthest
%  out combination (6th hole out on blade)
%  is lower, 68.2 lb.


        