% Ryan D Smith % Program utilized Vector Loop Equation and graphical methods to % approximate angles three and length 4 for a crank-slider type problem. % Closed loop position solution per Hall's Text from Purdue Univ. Mech. Eng. % Velocity and acceleration eqn's from taking derivatives of postion eqn. % VLE == R2 + R3 - R4 = 0 % Assumes crank pivot (crankshaft CL) is at same y value as r4 (slider CL). % CHECK UNITS! Default is angle in degrees (for I/O). Lenghts in inches. % 1 2 3 4 5 6 7 %2345678901234567890123456789012345678901234567890123456789012345678901234 clear clc clf % DECLARATIONS % Initialize Dynamic Inputs rpm = input('Enter crankshaft RPM: '); alpha = input('Enter crankshaft angular acceleration rad/sec^2: '); % Initialize 4 bar link lengths r2 = input('Enter crank throw length in inches: '); % input link r3 = input('Enter con-rod length in inches: '); % con-rod % Enter angles in DEGREES a2 = 0; % Input Angle a3 = 0; % Con-rod Angle % First order Kinematic Variables a3p = 0; % Theta 3 prime: first order Kc r4p = 0; % R4 prime: first order Kc % Second Order Kinematic Variables a3pp = 0; % Theta 3 double prime: second order Kc r4pp = 0; % Length 4 double prime: second order Kc % Misc. Variables index = 0; % counter for incrimenting trials plot_count = 0; % counter for the plot loop junk = 10; % general junk/misc. multi-use variable; initialized arb. value deg2rad = pi / 180; % Converts degrees to radians rad2deg = 180 / pi; % Converts radians to degrees % Convert angles to radians a2 = a2 * deg2rad; rpm = rpm * (pi / 30); % Converts input RPM to rad/sec % Velocity, Accel, Position v_max = -100000000; v_min = 1000000000; a_max = -1000000000; a_min = 100000000; X_avg = 0; % Average piston velocity X_cnt = 0; % Counter for average velocity % Header for results table fprintf('\nTDC = zero degrees crank angle\n') fprintf('Piston displacement (X) measured from TDC = 0\n\n') fprintf(' Crank Angle Position\t Velocity\t Accel.\t\t\n') fprintf(' (deg)\t\t (in.)\t (MPH)\t\t\t(g''s)\n') % Counting loop to incriment 360 input positions junk = 10; % Initialized for counting to 10 for table display for index = 0 : 1 : 360 % Solving for Position unknowns a2 = index * deg2rad; % Input angle converted to radians % Calculating piston position r4 = r2 * cos(a2) + sqrt((r2 ^ 2) * (cos(a2)) ^ 2 - (r2 ^ 2) + (r3 ^ 2)); % Calculating con-rod angle wrt x-axis a3 = asin((-1 * (r2 * sin(a2))) / r3); % Solving for first order Kc % Calculating con-rod angluar velocity a3p = (-1 * r2 * rpm * cos(a2)) / (r3 * cos(a3)); % Calculating piston velocity r4p = (-1 * r3 * sin(a3)) - (r2 * rpm * sin(a2)); % Solving for second order Kc % Calculating con-rod angular acceleration a3pp = ((-1 * alpha * cos(a2)) + (r2 * (rpm ^ 2) * sin(a2)) + (r3 * (a3p ^ 2) * sin(a3))) / (r3 * cos(a3)); % Calculating piston acceleration r4pp = -1 * ((r3 * sin(a3)) + (r2 * alpha * sin(a2)) + (r2 * (rpm ^ 2) * cos(a2)) + (r3 * (a3p ^ 2) * cos(a3))); % Converting results for display or next use. Note: variable names % correspond to same variable; _disp added to denote display variable. % First occurance here may represent initialization. a2_disp = index; % Uses index instead of converting a2g back to degrees (same thing) a3_disp = a3 * rad2deg; % Convert final angle 3 to degrees for display r4_disp = r4; % Final piston position for % Putting results in degrees into matrix for later use a2_rad(index + 1) = index; % Radian matrix of a2 a2_mat(index + 1) = a2_disp; % Use index + 1 b/c matrix pos. 0 undefined. Unravel later; a3_rad(index + 1) = a3; % Radian matrix of a3 a3_mat(index + 1) = a3_disp; r4_rad(index + 1) = r4; r4_mat(index + 1) = r4_disp; a3p_mat(index + 1) = a3p; % rad/sec r4p_mat(index + 1) = r4p; a3pp_mat(index + 1) = a3pp; % rad/sec/sec r4pp_mat(index + 1) = r4pp; % Calculating Pos, Vel, Accel X = r4 - (r2 + r3); % inches from TDC X_dot = [(r4p) / (12)] * (3600 / 5280); % mi/hr X_2dot = ((r4pp) / 12) / 32.2; % units = g's X_mat(index + 1) = X; X_dot_mat(index + 1) = X_dot; X_2dot_mat(index + 1) = X_2dot; if(X_dot < v_min) v_min = X_dot; v_min_a2 = index; elseif (X_dot > v_max) v_max = X_dot; v_max_a2 = index; end if(X_2dot < a_min) a_min = X_2dot; a_min_a2 = index; elseif(X_2dot > a_max) a_max = X_2dot; a_max_a2 = index; end if(X_dot >= 0) X_avg = X_avg + X_dot; X_cnt = X_cnt + 1; end % Displaying results for every 10 deg. input angle from 0 to 360 deg. if junk == 10 fprintf('\t%3.0f\t\t', a2_disp) if(X < 0) fprintf(' %3.3f\t\t', X) else fprintf(' %3.3f\t\t', X) end if(X_dot < 0) fprintf('%3.3f\t\t', X_dot) else fprintf(' %3.3f\t\t', X_dot) end if(X_2dot < 0) fprintf(' %3.3f\n', X_2dot) else fprintf(' %3.3f\n', X_2dot) end junk = 1; % Reset junk to start counting for next block of 10 else junk = junk + 1; % Incrementing the counter to display text every 10 deg. end % End If end % End For fprintf('\nMinimum Velocity (MPH): %3.3f @ %3.3f crank angle\n', v_min, v_min_a2) fprintf('Maximum Velocity (MPH): %3.3f @ %3.3f crank angle\n', v_max, v_max_a2) fprintf('Average Velocity (MPH): %3.3f\n', X_avg/X_cnt) fprintf('Minimum Accel. (g''s): %3.3f @ %3.3f crank angle\n', a_min, a_min_a2) fprintf('Maximum Accel. (g''s): %3.3f @ %3.3f crank angle\n', a_max, a_max_a2) % Plotting results for 0 - 360 deg. a2 subplot(2,3,1), plot(a2_mat, X_mat) title('Pos. (below TDC) Vs. Crank Angle (TDC = 0)') ylabel('Pos (inches)'), xlabel('Crank Angle (deg.)'),grid subplot(2,3,3), plot(a2_mat, X_dot_mat) title('Velocity Vs. Crank Angle (TDC = 0)') ylabel('Velocity (mi/hr)'), xlabel('Crank Angle (deg.)'),grid subplot(2,3,5), plot(a2_mat, X_2dot_mat) title('Accel. Vs. Crank Angle (TDC = 0)') ylabel('Accel. (g''s)'), xlabel('Crank Angle (deg.)'),grid