Force Analysis Padmanabhan Kumar & Tristan Hill

Introduction & Problem Statement

A computer case requires a door to be operated my a motor through a mechanism such that the door rests on the top of the case while it is completely open. The four bar mechanism was chosen for the purpose. The synthesis of the four bar was done using burmester's synthesis and thus the link lengths were obtained. The bearing reactions(forces) are to be computed by the matrix method of force analysis.

Four Position Problem - Burmester's Synthesis

the precision positions should go here

The linkage was synthesized starting with a prescribed motion of the door in the form of four precision positions. At each precision position, the angle of the door is specified. A single fourbar solution consists of two dyads, an input and an output, and each dyad is described by its ground point and corresponding coupler point. From the precession positions, Burmester's synthesis is used to produce the possible ground point and coupler point curves, M and K, respectively. Burmester's curves give n squared possible combinations of input and output dyads, where n is equal to length of either the ground or coupler point curve.

If there is 100 possible ground points, then there are 10000 possible fourbars. Out of all the possible combinations only a small number are viable solutions, leaving most of the combinations with problematic with either branch defect, Grashoff condition, or both. The sets with branch defect pass through two of the four positions in one branch and the other two in the other branch, while the sets showing Grashoff condition do not allow for the full rotation of the input link.

As the input parameters are varied the ratio of viable to non-viable solution sets given by the Burmester curves varies greatly. For this problem the ratio was very low, ranging from approximately 15:10000 to 200:10000 viable to non viable solutions. An exhaustive search through the possible sets is used to locate viable dyad pairs. The distance formula is used to judge if a given solution is good by determining if the coupler point passes through each precision point with the prescribed angle. For this problem a single solution was chosen from the possible viable solutions for further analysis.

Design of Fourbar mechanism for computer case door

The mechanism was synthesized using burmester's synthesis for the door(coupler) to pass through 4 different positions. The position of the door in different position is as shown below.


Position, Velocity and Acceleration for the resulting four bar mechanism

Position Analysis

The loop equation of the four bar is given, there are two unknown angles to be solved

\begin{align} \vec{r_1} + \vec{r_2} + \vec{r_3} + \vec{r_4} = 0 \end{align}

Vectors 1 and 2 are combined into the know vector

\begin{align} \vec{r_k} = - (\vec{r_1} + \vec{r_2}) \end{align}

The angles and magnitudes are extracted from the vectors

\begin{align} \vec{r_k} \Rightarrow |\vec{r_k}| = r_k , \angle \vec{r_k} = \theta_k \end{align}
\begin{align} \vec{r_3} \Rightarrow |\vec{r_3}| = r_3 , \angle \vec{r_3} = \theta_3 \end{align}
\begin{align} \vec{r_4} \Rightarrow |\vec{r_4}| = r_4 , \angle \vec{r_4} = \theta_4 \end{align}

The coefficients A, B, and C are define to solve for the unknown angles

\begin{align} A = r_3^2-r_k^2-r_4^2-2r_k\cos(\theta_k)r_k\\ \end{align}
\begin{align} B=4r_k\sin(\theta_k)r_k\\ \end{align}
\begin{align} C=r_3^2-r_k^2-r_4^2+2r_k\cos(\theta_k)r_4\\ \end{align}
\begin{align} m = \left[ \frac {-B \pm \sqrt {B^2 - 4AC}}{2A} \right] \end{align}
\begin{align} \theta_4=\arctan(m) \end{align}
\begin{align} \theta_3=\arctan2\left[(\frac{r_k\sin(\theta_k)-r_4\sin(\theta_4)}{r_3}),(\frac{r_k\cos(\theta_k)-r_4\cos(\theta_4)}{r_3})\right] \end{align}

Recompose the vectors from their components

\begin{align} r_3,\theta_3 \Rightarrow \vec{r_3} \end{align}
\begin{align} r_4,\theta_4 \Rightarrow \vec{r_4} \end{align}

Position Analysis - an alternative method

This is the method used in Design of Machinery, Fourth Edition by Norton
The loop equation of the four bar is given by

\begin{align} \vec{r_1} + \vec{r_2} + \vec{r_3} + \vec{r_4} = 0 \end{align}


\begin{split} & \vec{r_1} - Ground Vector && \vec{r_2} - Vector for Link2 && \vec{r_3} - Vector for Link3 &&\vec{r_4} - Vector for Link4 \end{split}

For a fourbar, the link 2 is considered as the input i.e. the crank. Depending on the variation of the crank angle, the coupler angle and the follower angle can be found from the loop equation adn hence their respective positions.

\begin{align} \theta_4 = 2 \arctan \left[ \frac {-B \pm \sqrt {B^2 - 4AC}}{2A} \right] \end{align}
\begin{align} \theta_3 = 2 \arctan \left[ \frac {-E \pm \sqrt {E^2 - 4DF}}{2D} \right] \end{align}


\begin{align} A = \cos\theta_2 - k1-k2\cos\theta_2 + k3 \end{align}
\begin{align} B = -2\sin\theta_2 \end{align}
\begin{align} C = k1 -(k2+1)\cos\theta_2+k3 \end{align}
\begin{align} D = \cos\theta_2-k1-k4\cos\theta_2+k5 \end{align}
\begin{align} E = -2\sin\theta_2 \end{align}
\begin{align} F = -k1+(k4-1)\cos\theta_2+k5 \end{align}
\begin{align} k1= \frac {r_4}{r_1} \end{align}
\begin{align} k2= \frac{r_4}{r_3} \end{align}
\begin{align} k3= \frac{\left(r_1^2-r_2^2+r_3^2+r_4^2)\right}{2r_1r_3} \end{align}
\begin{align} k4= \frac{r_4}{r_2} \end{align}
\begin{align} k5= \frac{\left(r_3^2-r_4^2+r_1^2+r_2^2)\right}{2r_1r_2} \end{align}

Velocity Analysis

The Link1 is the ground link which is static. Therefore velocity of Link1 is zero

For Link2

\begin{align} \vec{r_2} = r_2e^{i\theta_2} \end{align}
\begin{align} v_2 = r_2 \dot{\theta_2} i e^{i\theta_2} \end{align}
\begin{align} v_{2x} = -r_2 \dot{\theta_2} \sin{\theta_2} \end{align}
\begin{align} v_{2y} = r_2 \dot{\theta_2} \cos{\theta_2} \end{align}
\begin{align} a_{2x} = -r_2 \ddot{\theta_2} \sin\theta_2 - r_2 \dot{\theta_2}^2 \cos\theta_2 \end{align}
\begin{align} a_{2y} = r_2 \ddot{\theta_2} \cos\theta_2 - r_2 \dot{\theta_2}^2 \sin\theta_2 \end{align}

For Link 3

\begin{align} \vec{r_3} = r_3e^{i\theta_3} \end{align}
\begin{align} v_3 = r_3 \dot{\theta_3} i e^{i\theta_3} \end{align}
\begin{align} v_{3x} = -r_3 \dot{\theta_3} \sin{\theta_3} \end{align}
\begin{align} v_{3y} = r_3 \dot{\theta_3} \cos{\theta_3} \end{align}
\begin{align} a_{3x} = -r_3 \ddot{\theta_3} \sin\theta_3 - r_3 \dot{\theta_3}^2 \cos\theta_3 \end{align}
\begin{align} a_{3y} = r_3 \ddot{\theta_3} \cos\theta_3 - r_3 \dot{\theta_3}^2 \sin\theta_3 \end{align}

For Link4

\begin{align} \vec{r_4} = r_4e^{i\theta_4} \end{align}
\begin{align} v_4 = r_4 \dot{\theta_4} i e^{i\theta_4} \end{align}
\begin{align} v_{4x} = -r_4 \dot{\theta_4} \sin\theta_4 \end{align}
\begin{align} v_{4y} = r_4 \dot{\theta_4} \cos{\theta_4} \end{align}
\begin{align} a_{4x} = -r_4 \ddot{\theta_4} \sin\theta_4 - r_4 \dot{\theta_4}^2 \cos\theta_4 \end{align}
\begin{align} a_{4y} = r_4 \ddot{\theta_4} \cos\theta_4 - r_4 \dot{\theta_4}^2 \sin\theta_4 \end{align}

For finding velocities

\begin{align} L_1 = \vec{r_1} + \vec{r_2} + \vec{r_3} + \vec{r_4} = 0 \end{align}
\begin{align} \frac {d L_1}{dt} = r_2 \dot{\theta_2} i e^{i\theta_2} + r_3 \dot{\theta_3} i e^{i\theta_3} + r_4 \dot{\theta_4} i e^{i\theta_4} = 0 \end{align}
\begin{align} \frac {d L_{1x}}{dt} = -r_2 \dot{\theta_2} \sin\theta_2 - r_3 \dot{\theta_3} \sin{\theta_3} - r_4 \dot{\theta_4} \sin{\theta_4} = 0 \end{align}
\begin{align} \frac {d L_{1y}}{dt} = r_2 \dot{\theta_2} \cos\theta_2 + r_3 \dot{\theta_3} \cos{\theta_3} + r_4 \dot{\theta_4} \cos{\theta_4} = 0 \end{align}

writing Eqns (36) and (37) in matrix form, we get

\begin{align} \begin{bmatrix} \--r_3\sin\theta_3& \--r_4\sin\theta_4\\ \-r_3\cos\theta_3& \-r_4\cos\theta_4 \end{bmatrix} \begin{bmatrix} \-\dot{\theta_3}\\ \-\dot{\theta_4} \end{bmatrix} = \begin{bmatrix} \-r_2\sin\theta_2\\ \--r_2\cos\theta_2 \end{bmatrix} \dot{\theta_2} \end{align}
\begin{align} \begin{bmatrix} \-\dot{\theta_3}\\ \-\dot{\theta_4} \end{bmatrix} = \begin{bmatrix} \--r_3\sin\theta_3& \--r_4\sin\theta_4\\ \-r_3\cos\theta_3& \-r_4\cos\theta_4 \end{bmatrix}^{-1} \begin{bmatrix} \-r_2\sin\theta_2\\ \--r_2\cos\theta_2 \end{bmatrix} \dot{\theta_2} \end{align}

Thus the velocities are found.

Acceleration analysis

For finding accelerations

we differentiate the loop equation twice with respect to time.

\begin{align} \frac {d^2 L_{1x}}{dt} = -r_2 \ddot{\theta_2} \sin\theta_2 - r_2 \dot{\theta_2}^2\cos\theta_2 - r_3 \ddot{\theta_3} \sin{\theta_3} - r_3 \dot{\theta_3}^2\cos\theta_3 - r_4 \ddot{\theta_4} \sin{\theta_4} - r_4 \dot{\theta_4}^2\cos\theta_4= 0 \end{align}
\begin{align} \frac {d^2 L_{1y}}{dt} = r_2 \ddot{\theta_2} \cos\theta_2 - r_2 \dot{\theta_2}^2\sin\theta_2 + r_3 \ddot{\theta_3} \cos{\theta_3} - r_3 \dot{\theta_3}^2\sin\theta_3 + r_4 \ddot{\theta_4} \cos{\theta_4} - r_4 \dot{\theta_4}^2\sin\theta_4= 0 \end{align}

writing Eqns (40) & (41) in matrix form we get

\begin{align} \begin{bmatrix} \--r_3\sin\theta_3& \--r_4\sin\theta_4\\ \-r_3\cos\theta_3& \-r_4\cos\theta_4 \end{bmatrix} \begin{bmatrix} \-\ddot{\theta_3}\\ \-\ddot{\theta_4} \end{bmatrix} = \begin{bmatrix} \-r_2\sin\theta_2\\ \--r_2\cos\theta_2 \end{bmatrix} \ddot{\theta_2} + \begin{bmatrix} \-r_2\cos\theta_2\\ \-r_2\sin\theta_2 \end{bmatrix} \dot{\theta_2}^2 + \begin{bmatrix} \-r_3\cos\theta_3\\ \-r_3\sin\theta_3 \end{bmatrix} \dot{\theta_3}^2 + \begin{bmatrix} \-r_4\cos\theta_4\\ \-r_4\sin\theta_4 \end{bmatrix} \dot{\theta_4}^2 \end{align}
\begin{align} \begin{bmatrix} \-\ddot{\theta_3}\\ \-\ddot{\theta_4} \end{bmatrix} = \begin{bmatrix} \--r_3\sin\theta_3& \--r_4\sin\theta_4\\ \-r_3\cos\theta_3& \-r_4\cos\theta_4 \end{bmatrix}^{-1} \left( \begin{bmatrix} \-r_2\sin\theta_2\\ \--r_2\cos\theta_2 \end{bmatrix} \ddot{\theta_2} + \begin{bmatrix} \-r_2\cos\theta_2\\ \-r_2\sin\theta_2 \end{bmatrix} \dot{\theta_2}^2 + \begin{bmatrix} \-r_3\cos\theta_3\\ \-r_3\sin\theta_3 \end{bmatrix} \dot{\theta_3}^2 + \begin{bmatrix} \-r_4\cos\theta_4\\ \-r_4\sin\theta_4 \end{bmatrix} \dot{\theta_4}^2\right) \end{align}

Force Analysis

After finding the respective velocities and accelerations, the modelling of the components was done using Pro-E to obtain the masses and the mass moment of inertia of the respective components. These are required to compute the forces in the joints whilst using the matrix method of force analysis. The component properties are shown below.

\begin{split} &m_2=71.83 g &&I_{g2}=1.1454*10^{-4} kgm^{2} &&m_3=2 kg &&I_{g3}=716.41*10^{-4} kgm^{2} &&m_4=92.69 g &&I_{g4}=2.3878*10^{-4} kgm^{2} \end{split}

Matrix method

'C' is defined by the matrix shown below.

\begin{align} \setcounter{MaxMatrixCols}{20} C = \begin{bmatrix} \-1&0&1&0&0&0&0&0&0&0&0&0&0\\ \-0& 1& 0& 1& 0& 0& 0& 0& 0& 0& 0& 0& 0\\ \--r_{21y}&r_{21x}&-r_{22y}&r_{22x}& 0& 0& 0& 0& 0& 0& 0& 0& 1\\ \-0& 0& 0& 0& 1& 0& 1& 0& 0& 0& 0& 0& 0\\ \-0& 0& 0& 0& 0& 1& 0& 1& 0& 0& 0& 0& 0\\ \-0& 1& 0& 1& 0& -r_{32y}& r_{32x}& -r_{33y}& r_{33x}& 0& 0& 0& 0\\ \-0& 0& 0& 0& 0& 0& 0& 0& 1& 0& 1& 0& 0\\ \-0& 0& 0& 0& 0& 0& 0& 0& 0& 1& 0& 1& 0\\ \-0& 0& 0& 0& 0& 0& 0& 0& -r_{43y}& r_{43x}& -r_{44y}& r_{44x}& 0\\ \-0& 0& 1& 0& 1& 0& 0& 0& 0& 0& 0& 0& 0\\ \-0& 0& 0& 1& 0& 1& 0& 0& 0& 0& 0& 0& 0\\ \-0& 0& 0& 0& 1& 0& 1& 0& 0& 0& 0& 0& 0\\ \-0& 0& 0& 0& 0& 1& 0& 1& 0& 0& 0& 0& 0\\ \end{bmatrix} \end{align}

and 'm' is defined by the matrix

\begin{align} m = \begin{bmatrix} \-m_2a_{g2x}\\ \-m_2a_{g2y} + m_2g\\ \-I_{g2}\alpha_2\\ \-m_3a_{g3x}\\ \-m_3a_{g3y}+m_3g\\ \-I_{g3}\alpha_3\\ \-m_4a_{g4x}\\ \-m_4a_{g4y}+m_4g\\ \-I_{g4}\alpha_4\\ \-0\\ \-0\\ \-0\\ \-0\\ \end{bmatrix} \end{align}

'F' is the force matrix that is defined by

\begin{align} F = \begin{bmatrix} \-F_{12x}\\ \-F_{12y}\\ \-F_{32x}\\ \-F_{32y}\\ \-F_{23x}\\ \-F_{23y}\\ \-F_{43x}\\ \-F_{43y}\\ \-F_{34x}\\ \-F_{34y}\\ \-F_{14x}\\ \-F_{14y}\\ \-T_{in}\\ \end{bmatrix} \end{align}

According to the matrix method,

\begin{equation} C * F = m \end{equation}

'F' is the force matrix and is given by

\begin{equation} F = C^{-1} * m \end{equation}

Virtual Work

Input Torque is computed through method of Virtual Work for comparison purposes. For this calculation, only the mass and inertia of the door (link 3) is considered.
The equation of virtual work becomes:

\begin{align} T_{in} \cdot \omega_2-(m_3A_{g3}) \cdot V_{g3} -(I_{g3}*\alpha_3) \cdot \omega_3 + (m_3g) \cdot V_{g3} = 0 \end{align}

Describe method here


The force analysis was done to obtain the bearing reactions in the joints for the problem as shown. From the results, it is very evident that the bearing reactions vary with respect to the change of position of the links and their respective weights. The problem is solved such that, if the mechanism were to be controlled by a motor driving the crank (Tin). The maximum values of the force reactions in the joints are as follows.

MATLAB Simulation


  • Black - Linkage
  • Yellow - Center of Gravity Vector
  • Green - CG to Coupler Point Vectors
  • Red - Velocity Vector
  • Blue - Acceleration Vector
  • Red Curve - Coupler Point Curve
  • Blue Curve - Ground Point Curve
  • Stars - Precision Points

CAD based animation of the mechanism

\begin{split} &F_{j1}=46.3542 N &&F_{j2}=45.6825 N &&F_{j3}=27.4609 N &&F_{j4}=26.6758 N &&T_in=3.5385 Nm \end{split}

While running the program, we find the variation in the Tin such that when the coupler (door) is to be lifted the value is positive. When the Link2 has turned through an angle that the weight of the door appends to reduce the torque required to drive the mechanism, the torque required for the input link reduces further and goes to the negative. This means that the weight of the door assists the motor in driving the Link2 by its weight. This also suggests that the mechanism requires the motor to drive through only that angle. The mechanism is self actuated by the weight of the coupler, after Link2 has turned through that angle, to bring the coupler (door) to rest on the top of the box. This also suggests that a drive is required to turn Link2 only through that angle.

Bearing Reactions

The variation of the forces on the joints with respect to the input link angle is shown in the graph below.

bearingrea ctions-2.png

Matrix Method vs. Virtual Work

The results of the Matrix Method and the Method of Work are compared for verification. The input torque versus input link angle, theta2, is graphed side by side for comparison. The two curves are almost identical, and the small amount of error is due to the fact that the mass and inertia of link 2 and link 4 were considered in the matrix method but neglected in the method of Virtual Work.


Solving in MATLAB - The various components of this analysis were solved in MATLAB.

Linkage Synthesis - Burmester

The ground and coupler points that define the linkage are solved in MATLAB using Burmester's synthesis

Position Analysis - Non Linear - solved with trigonometric substitution

The algorithm for the position analysis described above is computed in matlab. The sign sensitive atan2(Y,X) function is utilized.

Velocity Analysis - Linear - solved with matrix inverse

The velocity analysis described above is computed in MATLAB. The matrix inverse is utilized to solve the system of equations.

Acceleration Analysis - Linear - solved with matrix inverse

The acceleration analysis described above is computed in MATLAB. The matrix inverse is utilized to solve the system of equations.

Force Analysis - Linear - solved with matrix inverse

The acceleration analysis described above is computed in MATLAB. The needed input torque and resulting bearing reactions are solved using the Matrix Method. The method of Virtual Work is used for verification of the matrix method.


Code is available on Google docs, here