Brandon Hoshaw ENGR 322 11/12/13
Computer Assignment 2
Four Bar Using Global Shell % FourBarUsingGlobalShell.m
% Program to plot the coupler curve of the four-bar linkage % in problem 15.68 Beer et al. % Author: Jenny Holte, University of St Thomas % Date: Nov 10, 2013 % % Additional code by: Brandon Hoshaw % Date: Nov 12, 2013 % % Using fsolve and global variables (instead of newtonrap) % theta2_0= atan2(200,-100); % starting orientation of the input link (link 2) omega2 = -10; % constant angular velocity of the input link (link 2) global r2; global r1x; global r1y; global r3; global r4; global theta2; r1x = 400; r1y =120; r2 = sqrt(200^2+100^2); r3 = sqrt(300^2+100^2); r4 = 420; DG =400; theta34i = [150*pi/180; 90*pi/180];
% Make a starting guess at the solution
N=361; % Full 360-degree rotation at 1 degree increments. % Position and velocity Analysis %preallocate array space to save computing time xF = zeros(1,N); yF = zeros(1,N); vF = zeros(1,N); t = zeros(1,N); for I=1:N DeltaTheta2 = (I-1)*2*pi/(N-1); theta2 = theta2_0 - DeltaTheta2;
% due to clockwise rotation decrease % theta2 by Deltatheta2 at each step t(I) = DeltaTheta2/abs(omega2); % time lapsed will be used later % from omega2 = Dtheta2/Dtime
if theta2 > 2*pi theta2 = theta2 - 2*pi; % keep theta 2 between 0 and 2pi end; % Generate numercal solution for theta 3 and theta 4 % Function fourbar evaluates the non-linear functions for the % specifica value of theta2 (the loop equations) theta34 = fsolve(@LoopEqns,theta34i); % Call solver theta34i = theta34; % inital guess for next theta 3 and 4 % compute the position of point G xG(I) = r2*cos(theta2)+DG*cos(theta34(1)-(pi/2)); yG(I) = r2*sin(theta2)+DG*sin(theta34(1)-(pi/2)); % compute the velocity of G v1 = r2*sin(theta2)*omega2; % right side of velocity equation v2 = -r2*cos(theta2)*omega2; % specific to problem 15.68 v =[v1;v2]; J =jacob2x2(theta34); % Define the Jacobian matrix omega34 = inv(J)*v; % Solution to the linear system of eqns vxG = -r2*sin(theta2)*omega2-DG*sin(theta34(1)-(pi/2))*omega34(1); vyG= r2*cos(theta2)*omega2+DG*cos(theta34(1)-(pi/2))*omega34(1); vG(I) = sqrt(vxG^2 + vyG^2);
end plot(xG, yG) axis([-1000 1000 -1000 1000]) title('Coupler Curve(Position of G)') xlabel('x (mm)') ylabel('y (mm)') hgsave('CouplerCurve')
plot(t, vG) title('Velocity of Point G') xlabel('time t(s)') ylabel('Vel F (mm/s)') hgsave('CouplerVelocity')
Jacobian Matrix function [J] = jacob2x2(theta34) % Evaluates the Jacobian of a 2x2 % system of non-linear equations global r3; global r4; J(1,1) = -r3*sin(theta34(1)); % partial deriv df1/dx1 J(1,2) = r4*sin(theta34(2)); % partial deriv df1/dx2 J(2,1) = r3*cos(theta34(1)); % partial deriv df2/dx1 J(2,2) = -r4*cos(theta34(2)); % partial deriv df2/dx2 % end function
Loop Equations % LoopEqns.m
% Functions from Four bar linkage position solution function F = LoopEqns(theta34) global theta2; global r2; global r1x; global r1y; global r3; global r4; F1= r3*cos(theta34(1))-r4*cos(theta34(2))+r2*cos(theta2)+r1x; F2= r3*sin(theta34(1))-r4*sin(theta34(2))+r2*sin(theta2)+r1y; F = [F1;F2];