fbpx

Macauley’s Method for Fast Beam Deflections

by Dr Seán Carroll
|
Published: 28th July 2021
|
Tutorial

1.0 What is Macauley’s Method?

In this tutorial, we will again visit the topic of deflection calculation and introduce Macauley’s Method. We last discussed beam deflection in this post, where we worked out how to calculate deflections from first principles. We’ll again rely on integrating the differential equation of the deflection curve,

(1)   \begin{equation*} EI\frac{\mathrm{d}^2v}{\mathrm{d}x^2} = M(x) \end{equation*}

except this time, we’ll speed the whole process up considerably using a technique called Macauley’s method. It’s also known as the singularity function method. In this tutorial, I’ll assume you’re familiar with the process (and baked in assumptions) of integrating the differential equation of the deflection curve covered previously – if not, pause here and go and read this post first

Recall that the basic process of calculating beam deflection involves evaluating an expression for the bending moment as a function of x, the coordinate along the beam. We can substitute this expression into the equation above and perform a double integration of both sides of the equation to produce a closed-form solution (an equation) for the deflection as a function of x

In theory, this all works beautifully, but practically, there’s a catch – unless our bending moment can be described by a continuous function of x, we will need to break up our bending moment diagram into segments where a single continuous function defines each segment. The double integration required for each segment will yield two unknown constants of integration. This is all perfectly solvable, but the process gets very tedious for anything other than very simple loading. Take for example the beam we analysed in the previous post. 

Simply-supported-beam | DegreeTutors.com
Fig. 1: Simply supported beam subject to two point loads and a uniformly distributed load. (First analyses here).

This beam must be divided into 3 segments, resulting in six constants of integration that must be evaluated using boundary and continuity conditions. The final output of the analysis is three equations describing the beam’s deflection in the three beam regions, A-B, B-C and C-D. 

This is where Macauley’s method comes in. I’ll outline the headlines here – but this is only really going to make sense when we work through an example in the next section – so don’t worry if the following leaves you scratching your head! As the name suggests, Macauley’s method is really just a method or process we can use to evaluate the differential equation of the deflection curve. The method works by introducing a singularity function (hence the alternative name), which allows us to discount or ignore certain terms in the bending moment function M(x)

So, instead of building three different functions for the bending moment in the question above, we would build a single equation for M(x) by cutting the beam between C and D, thereby capturing the moment influence of all loads on the structure. The ‘trick’ that Macauley’s method introduces is that we multiply each term in the bending moment equation by a singularity function (will make more sense below) that we let equal to zero when it has a value less than one. The way we implement this allows us to cleverly eliminate terms from the moment equation when we’re evaluating sections of the beam, where those terms are not relevant.

For example, if we cut the beam between B and C and consider equilibrium of the sub-structure to the left of the cut, the 75\:kN point load would feature in the moment equation but the 50\:kN would not. Macauley’s method simply allows us to construct a single moment equation and eliminate any non-relevant terms. This method results in us performing double integration on one equation (not three) yielding only two unknown constants of integration – greatly reducing our workload. In the next section we’ll flesh out this idea with an example. 

2.0 Macauley’s Method Example

To demonstrate the benefits of Macauley’s method, let’s start by re-analysing the structure pictured in Fig. 1 above. We’ll break the process up into discrete repeatable steps that you can replicate on other beams. 

2.1 Step 1: Sub-structure setup

We need to generate an equation that describes the internal bending moment as a function of x. To do this, we make a cut in the structure at some distance x to the right of the left support. We’ll treat the left-most tip of the beam as the origin of our x-axis. When implementing Macauley’s method, we make sure we cut the structure at a location that allows us to consider all of the actions (applied forces and moments) to the left of the cut. This typically means making a cut some infinitesimally small distance to the left of the right tip of the beam. When taking moments about the cut, notice that all of our lever arm terms contain x. This is critical as these lever arm terms will be treated as our singularity functions.

Macauley's Method 1 | DegreeTutors.com
Fig. 2: Sub-structure for Macauley’s method analysis

2.2 Step 2: Build an expression for M(x) with singularity functions

Taking moments about the cut and recalling the vertical reaction at A was previously determined to be V_A=139.375\:kN, yields the following expression for M(x)

(2)   \begin{equation*}M(x) = 139.375[x] - 75[x-3] - \frac{20}{2}[x]^2 - 50[x-6]\end{equation*}

Pay special attention to the terms in square brackets, these are our singularity functions. Macauley’s method says that when the terms in brackets evaluate to a negative number, we treat the term in brackets as a zero. Think about this for a minute – what is this actually saying? Imagine we make a cut in the structure at x=5\:m to determine the internal bending moment at this location by considering equilibrium of the sub-structure to the left of the cut. Naturally the 50\:kN force at C would not feature in this equation. The ‘rule’ we just implemented caters for this since the singularity function for the 50\:kN force would evaluate to 5-6=-1 and according to Macauley’s method, this would be eliminated from the equation! Therefore after applying Macauley’s method to the moment expression we would have,

(3)   \begin{equation*}M(x) = 139.375[5] - 75[2] - \frac{20}{2}[5]^2 - \underbrace{50[-1]}_{\text{eliminated}}\end{equation*}


(4)   \begin{equation*}M(x) = 139.375[x] - 75[x-3] - \frac{20}{2}[x]^2 \end{equation*}

You can think of what we’ve done so far as building the most complete expression for the internal bending moment and then implementing a rule that allows us to eliminate terms from the expression as required.

2.3 Step 3: Integrate the differential equation of the deflection curve

Now we can substitute our expression for the internal bending moment, complete with singularity functions (which are really just the lever arm expressions) into the differential equation of the deflection curve and integrate the equation. 

(5)   \begin{equation*}EI\frac{\mathrm{d}^2v}{\mathrm{d}x^2} = 139.375[x] - 75[x-3] - \frac{20}{2}[x]^2 - 50[x-6]\end{equation*}

Before integrating, note that we are keeping the singularity functions intact and not expanding them or multiplying through the terms to the left of the brackets. After integrating twice we end up with,

(6)   \begin{equation*}EI\frac{\mathrm{d}v}{\mathrm{d}x} = \frac{139.375}{2}[x]^2 - \frac{75}{2}[x-3]^2 - \frac{20}{6}[x]^3 - \frac{50}{2}[x-6]^2 + C_1\end{equation*}


(7)   \begin{equation*}EI\: v = \frac{139.375}{6}[x]^3 - \frac{75}{6}[x-3]^3 - \frac{20}{24}[x]^4 - \frac{50}{6}[x-6]^3 + C_1x + C_2\end{equation*}

where C_1 and C_2 are the constants of integration. 

2.4 Step 4: Apply boundary conditions and solve for constants of integration

We can now apply our boundary conditions in the usual way to determine the unknown constants C_1 and C_2

Boundary condition 1:

At x=0 the deflection v=0. Note that when x=0, all of our singularity functions either become zero or negative which according to Macauley’s method we treat as zero,

(8)   \begin{equation*}EI\: v = \frac{139.375}{6}[0]^3 - \frac{75}{6}[-3]^3 - \frac{20}{24}[0]^4 - \frac{50}{6}[-6]^3 + C_1(0) + C_2\end{equation*}

Therefore, C_2=0.

Boundary condition 2: 

At x=L, v=0. This boundary condition yields,

(9)   \begin{equation*}EI\: v = \frac{139.375}{6}[L]^3 - \frac{75}{6}[L-3]^3 - \frac{20}{24}[L]^4 - \frac{50}{6}[L-6]^3 + C_1(L) + C_2\end{equation*}

Since L=8, we end up with C_1=-856.4. Therefore our final equation for the deflection becomes,

(10)   \begin{equation*}\boxed{EI\: v = \frac{139.375}{6}[x]^3 - \frac{75}{6}[x-3]^3 - \frac{20}{24}[x]^4 - \frac{50}{6}[x-6]^3 - 856.4x}\end{equation*}

This essentially marks the end of the process and we can now use our equation to evaluate the deflection as required. So for example, let’s consider the deflection at the mid-span location, at x=4. Substituting x=4 into our equation and observing Macauley’s method of eliminating terms where the singularity function evaluates to a negative would yield,

(11)   \begin{equation*}EI\: v = \frac{139.375}{6}[4]^3 - \frac{75}{6}[1]^3 - \frac{20}{24}[4]^4 - \frac{50}{6}[-2]^3 - 856.4(4)\end{equation*}


(12)   \begin{equation*}EI\: v = 1486.7-12.5-213.3-0-3425.6\end{equation*}


(13)   \begin{equation*}\boxed{v = \frac{-2164.7}{EI}}\end{equation*}

We can confirm that this is the same value we derived in our previous analysis of this beam. However, we only needed to integrate a single equation and evaluate two unknown constants of integration. The process was a lot quicker this time around. In a nutshell, this is Macauley’s method – fundamentally we’re still integrating the differential equation of the deflection curve, we’re just being a little smarter about it. 

3.0 Plotting the Deflected Shape

Once we have our equation for the deflection, an excellent next step would be plotting the equation and visualising the deflected shape of the beam. We can write a simple algorithm to achieve this. I’ll use Python inside a Jupyter notebook, but you can use any language you’re comfortable with. If you’re new to Python, take a look at this project that will help you get up and running with little or no experience. 

The first thing I do inside a blank Jupyter notebook is import some packages to give me additional functionality. On line 2 below, I import Numpy to help me work with arrays of numbers, and on lines 3-5, I’m importing Plotly, a graphing library that helps me produce nice graphs.

With these preliminaries out of the way, we can get on to the core of the code. In lines 3-6 below we set up a range of x-coordinates along the length of the beam. On line 10 we initialise a container to hold the value of deflection at each x-coordinate and on line 11 we establish a for loop to cycle through each x-coordinate and calculate the value of deflection at each point.

We’ve divided the expression for deflection into different parts which allows us to use if conditions to test if the current value of x renders the singularity function zero or not. This way we’re implementing Macauley’s method and only taking into account the relevant parts of the deflection equation for each position on the beam. 

After stepping through all of our test conditions and building up our deflection for the current location, we save that value on line 39 and move onto the next x-location with the next iteration of our for loop and repeat the process all over again.

Now that we’ve calculated an array of deflection values for each position along the beam, we can use Plotly to visualise them. The code below is pretty standard plotting code for Plotly. You can use whatever programme or library you like to visualise the deflection; all you’re doing is plotting a range of x versus y values. I include the code below more for completeness than anything else. 

Macauley's Method 2 | DegreeTutors.com
Fig. 3: Beam deflection obtained using Macauley’s method

Comparison of the deflection shown in Fig. 3 above with the plot of deflection from our previous analysis here will confirm that we’ve achieved the same result. 

4.0 A Tougher Example – Macauley’s Method with Partial UDLs and Point Moments

The previous example was sufficient to demonstrate the advantages of Macauley’s method over our previous approach to calculating the deflection. However, before we wrap up our discussion, we should probably work through a slightly more challenging example question. This will give us the chance to address partial uniformly distributed loads (patch loads) and point moments. To use Macauley’s method to deal with these actions, we need to implement specific strategies. So with this in mind, consider the structure below. Our task is to plot the deflected shape for the structure as we did for the previous example. We’ll also determine the magnitude and location of the maximum deflection between A and B. 

Macauley's Method 4 | DegreeTutors.com
Fig. 4: Example beam

4.1 Beam reactions

We can start by calculating the vertical reactions at A and B. Taking the sum of the moments about point A,

(14)   \begin{equation*}\sum M_{cut}=0\:\:(\curvearrowright\:+ve)\end{equation*}


(15)   \begin{equation*}-60+(50\times 4\times 4) + (100\times 8) +(75\times 13) - 10V_B=0\end{equation*}


(16)   \begin{equation*}V_B=251.5\end{equation*}

Evaluating the sum of the forces in the vertical direction yields V_A = 123.5\:kN.

4.2 Macauley’s method and patch loading

The arrangement of loads on our structure causes us two problems that we’ll need to consider if we’re to use Macauley’s method,

  • the fact that the distributed load doesn’t run all the way to the end of the beam on the right hand side causes us a problem which we’ll solve next.
  • the existence of a moment needs to be addressed – because this does not have a lever arm, we don’t naturally have a singularity function that we can apply to it – this one is an easy fix which we’ll see shortly. 

Let’s simplify things for demonstration purposes and imagine our beam only had a distributed load acting on it. We can set up the sub-structure to the left of the cut as shown below in Fig. 5.

Macauley's Method 5 | DegreeTutors.com
Fig. 5: Sub-structure with uniformly distributed load only

Our expression for the bending moment at the cut would be,

(17)   \begin{equation*} M(x) = V_A[x-3] -4w[x-7] + V_B[x-13] \end{equation*}

Notice that our expression has the full UDL, 4W ‘baked’ into it. This means that this equation can only be valid for values of 9\leq x\leq16. This typically wouldn’t be a problem because our singularity functions usually eliminate the influence of a load when the cut is to the left of the load. However, notice that the singularity function only eliminates the influence of the UDL when x<7. So there is a window when 7<x<9 when our equation will give us the wrong bending moment.

To overcome this we need our UDL to always run to the end of the beam, in other words our vertical cut must always cut through whatever UDLs are applied to the beam. An easy way to satisfy this requirement is to extend our UDL to the end of the beam and then apply a UDL in the opposite direction to cancel the influence of the extra UDL we’ve added – this makes more sense when you see it in action, Fig. 6 below. 

Macauley's Method 6 | DegreeTutors.com
Fig. 6: Sub-structure with all example loading shown. Note that the forces in green are the resultants of the UDL loads.

The combined influence of both UDLs in Fig. 6 is the same as the influence of the single UDL shown in Fig. 4.

4.3 Macauley’s method and point moments

We said above that the existence of point moments caused a problem because a moment doesn’t get multiplied by a lever arm in our moment equation and therefore there is no singularity function to eliminate the moment when it doesn’t contribute to M(x). Now in our example question this problem never actually arises because the moment at C will always feature in the equation for M(x), i.e. we never evaluate the internal bending moment to the left of the 60\:kNm, because it’s applied at the end of the beam. 

However, to illustrate the problem and solution for future reference, let’s imagine the same beam with a moment, M_o applied at say x=5m, Fig. 7 below. In this case the internal moment expression would be,

(18)   \begin{equation*}M(x) = V_A[x-3] - M_o + V_B[x-13] \end{equation*}

Macauley's Method 7 | DegreeTutors.com
Fig. 7: Sub-structure with single point moment applied.

Currently we have no way of eliminating the moment when x<5. However if we introduce a ‘sudo’ lever arm, [x-5]^0, we give ourselves a singularity function that enables us to apply Macauley’s method to this moment and eliminate it. Note that [x-5]^0 = 1, so we never actually alter the magnitude of the moment when it should feature in our equation. So, the complete expression for the internal moment, that we can apply Macauley’s method to is,

(19)   \begin{equation*}M(x) = V_A[x-3] - M_o[x-5]^0 + V_B[x-13] \end{equation*}

 

But as we said, we don’t need to implement this ‘fix’ in our example because the moment will always feature in our equation because it’s positioned on the extreme left end of the beam. So, with these two details covered, we can get back to solving our question.

4.4 Build an expression for M(x) with singularity functions and integrate

We have the complete sub-structure with all loading and lever arms shown in Fig. 6 above. Taking moments about the cut we can build our expression for M(x) as follows,

(20)   \begin{equation*}M(x) = -60[x-0]^0 + 123.5[x-3]-\frac{50}{2}[x-5]^2-100[x-11]+\frac{50}{2}[x-9]^2+251.5[x-13]\end{equation*}

We now substitute this expression into our differential equation of the deflection curve,

(21)   \begin{equation*}EI\frac{\mathrm{d}^2v}{\mathrm{d}x^2} = -60[x-0]^0 + 123.5[x-3]-\frac{50}{2}[x-5]^2-100[x-11]+\frac{50}{2}[x-9]^2+251.5[x-13]\end{equation*}

Integrating to obtain an expression for the slope,

(22)   \begin{equation*}EI\frac{\mathrm{d}v}{\mathrm{d}x} = -60[x] + \frac{123.5}{2}[x-3]^2-\frac{50}{6}[x-5]^3-\frac{100}{2}[x-11]^2+\frac{50}{6}[x-9]^3+\frac{251.5}{2}[x-13]^2 + C_1\end{equation*}

Integrating again to obtain our expression for the displacement,

(23)   \begin{equation*}EI v = -\frac{60}{2}[x]^2 + \frac{123.5}{6}[x-3]^3-\frac{50}{24}[x-5]^4-\frac{100}{6}[x-11]^3+\frac{50}{24}[x-9]^4+\frac{251.5}{6}[x-13]^3 + C_1x + C_2\end{equation*}

Boundary condition 1: v(x=3)=0:

Applying the boundary condition for the displacement at A yields,

(24)   \begin{equation*}3C_1 + C_2 = 270\end{equation*}

Boundary condition 2: v(x=13)=0:

The boundary condition for the displacement at B yields,

(25)   \begin{equation*}13C_1 + C_2 = -7380\end{equation*}

Solving the simultaneous boundary condition equations we obtain C_1 = -765 and C_2=2565. Therefore our final expression for the displacement of the beam is obtained as,

(26)   \begin{equation*}\boxed{v = \frac{1}{EI}\left \{  -\frac{60}{2}[x]^2 +\frac{123.5}{6}[x-3]^3 - \frac{50}{24}[x-5]^4 - \frac{100}{6}[x-11]^3 + \frac{50}{24}[x-9]^4 + \frac{251.5}{6}[x-13]^3 - 765[x] + 2565  \right\}}\end{equation*}

4.5 Plotting the deflected shape

We will employ the exact same technique as before to plot the deflected shape. The code below is structured the same as the code we’ve seen above but the contents are obviously altered to represent our new deflection equation. The code to produce the plot is the exact same and is not repeated below. 

Macauley's Method 8 | DegreeTutors.com
Fig. 8. Deflected shape for the beam. Note that in the calculation, the flexural rigidity EI is assumed to be 1 for convenience

It might seem odd that the 75\:kN point load at G is never directly featured in our internal bending moment equation, M(x). However, its influence is taken into account indirectly through the reactions at A and B. You can confirm this by simply increasing the point load to say 100\:kN and repeating the calculation. This will require you to calculate new reactions, re-evaluate the constants of integration and slightly amend the code block above to capture the changes – but it’s not as laborious as it sounds. I’ve plotted a comparison of the two deflections in Fig. 9 below. 

Macauley's Method 9 | DegreeTutors.com
Fig. 9: Comparison of deflected shapes for a 75 kN point load at G (blue) and a 100 kN point load at G (green)

The influence of the larger point load at G is clearly visible and qualitatively it’s exactly what we would expect. 

4.6 Maximum deflection between A and B

The last thing we want to tackle is identifying the location and magnitude of the maximum span deflection between A and B. This is relatively easy to obtain as the maximum deflection occurs where the slope of the deflected shape is zero. Fortunately we have an equation for the deflected shape already – so we just need to let it equal to zero and solve for x…easier said than done! Take a look at the equation again (repeated from above),

(27)   \begin{equation*}EI\frac{\mathrm{d}v}{\mathrm{d}x} = -60[x] + \frac{123.5}{2}[x-3]^2-\frac{50}{6}[x-5]^3-\frac{100}{2}[x-11]^2+\frac{50}{6}[x-9]^3+\frac{251.5}{2}[x-13]^2 -765\end{equation*}

Solving this equation algebraically for x is too difficult. A much more practical approach is to visually identify the region of interest, say between 7< x< 10 meters and use trial and error, i.e. plug values of x into the equation and see which one makes the equation equal to zero. In the dark ages one might have done this by hand – but lucky for us we can put a robot on the job instead – via a Python script!

Our approach will be to use a for loop to cycle along each position on the beam and calculate the slope, in much the same way that we calculated the deflection. At the end of each loop iteration we’ll test if the sign of the slope has changed from negative to positive. When it has, it means that the point of zero slope was crossed between the previous and current iteration. We’ll take the x-position on the current iteration as the point at which the slope was zero and therefore the location of the maximum deflection.

Our script outputs the following text string,

The maximum deflection between A and B occurs at x=7.781 m

We can now substitute this value for x in our deflection equation and obtain the deflection at this location,

(28)   \begin{equation*}EI\:v_{max} = -\frac{60}{2}[7.781]^2+\frac{123.5}{6}[4.781]^3-\frac{50}{24}[2.781]^4-765[7.781]+2565\end{equation*}


(29)   \begin{equation*}\boxed{v_{max} = \frac{-3078.5}{EI}}\end{equation*}

That wraps up our analysis. I hope if you’ve managed to get this far, you can see the benefit of Macauley’s method for speeding up the deflection calculation process. In the days of computerised structural analysis, we rarely need to process a deflection like this – but it’s still nice to keep in touch with the fundamental principles. If you found this post helpful, you might like this project where we build a beam deflection calculator in Python – if you think Macauley’s Method is fast, building your own calculator makes calculating beam deflections lightning fast.

Speaking of Python and computerised structural analysis, if you want to continue exploring how to use programming within your structural analysis, take a look at the direct stiffness method course bundle below. 

The Direct Stiffness Method for Truss Analysis with Python | DegreeTutors.com
Play Video

The Direct Stiffness Method for Truss Analysis with Python

Build your own finite element truss analysis software using Python and tackle large scale structures.

After completing this course…

  • You’ll understand how to use the Direct Stiffness Method to build complete structural models that can be solved using Python.
  • You’ll have your own analysis programme to identify displacements, reactions and internal member forces for any truss.
  • You’ll understand how common models of elastic behaviour such as plane stress and plane strain apply to real-world structures.
Play Video

Beam & Frame Analysis using the Direct Stiffness Method in Python

Build a sophisticated structural analysis software tool that models beams and frames using Python.

After completing this course…

  • You’ll understand how to model beam elements that resist axial force, shear forces and bending moments within the Direct Stiffness Method.
  • You’ll have your own analysis software that can generate shear force diagrams, bending moment diagrams, deflected shapes and more.
  • You’ll understand how to model rotational pins within your structures, inter-nodal distributed loading and realistic flexible supports.
Play Video

3D Space Frame Analysis using Python and Blender

Develop tools to model and analyse complex 3D space frame structures using Python.

After completing this course…

  • You’ll understand how to apply the Direct Stiffness Method to solve 3D space frame structures.
  • You’ll have your own analysis programme to identify displacements, reactions and internal member forces for any 3D space frame.
  • You’ll be able to use Blender, a powerful open source 3D modelling software to build, visualise and export your structural models.

That’s all for now, perhaps for our next analysis, we’ll look at speeding up this calculation even more by evaluating the differential equation of the deflection curve numerically – using only the bending moment diagram…watch this space!

Tutorial by:

Dr Seán Carroll

Hi, I’m Seán, the founder and lead tutor at DegreeTutors.com. I hope your found whatever I wrote up there helpful! Feel free to get in touch or follow DegreeTutors on any of the social accounts.