Singularity functions are a mathematical tool used to describe systems where discontinuity occurs. These functions are written using a special notation using singularity brackets: . The general functions is written as
where
is the point of discontinuity and
is the function degree. A detailed explanation of the properties of these functions can be found on Wikipedia.
The singularity functions method is sometimes refereed to as “Macaulay’s method” and the brackets as “Macaulay brackets”.
The main advantage of using singularity functions for solving the continuous beam problem is the ability to describe and solve the system analytically with any kind of load. This can be done for complex loading systems such as having the load described as a high order function. The symbolic solution can be reached without any approximation which gives us a pure analytical solution to any continuous beam.
To demonstrate this, the following continuous beam is used as an example:
![](https://www.rajjoub.net/wp-content/uploads/2021/03/0001.png)
The load is described as a function of the length starting at
to
from left to right. Equation (1) describes the load function; where:
and
.
(1)
First we should determine the reaction forces Ra, Rb and Rc.
![](https://www.rajjoub.net/wp-content/uploads/2021/03/0003-1024x269.png)
This is done by simplifying the load as follows:
![](https://www.rajjoub.net/wp-content/uploads/2021/03/0002.png)
(2)
Using simple static we can determine the system equations 3 and 4 .
(3)
(4)
We start the solution by writing the equation of normal forces as shown in 5. To get the equation, a couple of rules should be known: for single loads we use a singularity function with a degree of , constant line loads get a degree of
and so on. Rules in applying singularity functions can be seen in this link.
The challenge occurs when having loads with higher order (such as the one in our example). In this case the load should be described as a Taylor expansion. The reason for this is explained in this discussion.
(5)
Now we actually finished the most important part on the solution and the rest is simple mathematics!
To get the shear force equation, we integrate the force equation 5 as seen in equation 6
(6)
Similarly we get the bending moment equation 7.
(7)
Rotation is calculated in 8. Here we should take an integration constant into account. This is due to the fact, that the rotation at the
doesn’t equal 0.
(8)
Finally we can derive the deflection equation as in 9.
(9)
Determining : when
then
therefore by substituting in the deflection equation we determine the following:
When then
therefore by substituting in the deflection equation we determine 10.
(10)
Now by solving equations 3,4 and 10 we can determine the reaction forces:
(11)
![](https://www.rajjoub.net/wp-content/uploads/2021/03/0004-1024x216.png)
Since we already determined the equations of the internal forces and calculated the reactions, we can simply plot the equations by varying the value of .
![](https://www.rajjoub.net/wp-content/uploads/2021/02/image-8.png)
![](https://www.rajjoub.net/wp-content/uploads/2021/02/image-9.png)
![](https://www.rajjoub.net/wp-content/uploads/2021/02/image-10.png)
![](https://www.rajjoub.net/wp-content/uploads/2021/02/image-11.png)
The following is the VBA code used to plot the solution. It can be used as a reference on the proper way to implement singularity functions in VBA/Excel.
Public Const Ra = 22.524
Public Const Rb = -135.048
Public Const Rc = 472.524
Public Const q0 = 60
Public Const E = 2.1 * 10 ^ 8
Public Const I = 3945001.5 * 10 ^ -8
Public Const EI = E * I
Function Shear(x)
y1 = Abs(x >= 0) * Ra * (x - 0) ^ 0
y2 = Abs(x >= 9) * Rb * (x - 9) ^ 0
y3 = Abs(x >= 18) * Rc * (x - 18) ^ 0
y4 = Abs(x >= 15.75) * (-4 * q0 / 18) * (x - 15.75) ^ 2
y5 = Abs(x >= 15.75) * (4 * q0 / 243) * (x - 15.75) ^ 3
Shear = Abs(x <= 24.75) * (y1 + y2 + y3 + y4 + y5)
End Function
Function Moment(x)
y1 = Abs(x >= 0) * Ra * (x - 0) ^ 1
y2 = Abs(x >= 9) * Rb * (x - 9) ^ 1
y3 = Abs(x >= 18) * Rc * (x - 18) ^ 1
y4 = Abs(x >= 15.75) * (-4 * q0 / 54) * (x - 15.75) ^ 3
y5 = Abs(x >= 15.75) * (4 * q0 / 972) * (x - 15.75) ^ 4
Moment = -Abs(x <= 24.75) * (y1 + y2 + y3 + y4 + y5)
End Function
Function Rotation(x)
If x > 24.75 Then
Rotation = Rotation(24.75)
Exit Function
End If
y1 = Abs(x >= 0) * (1 / 2) * Ra * (x - 0) ^ 2
y2 = Abs(x >= 9) * (1 / 2) * Rb * (x - 9) ^ 2
y3 = Abs(x >= 18) * (1 / 2) * Rc * (x - 18) ^ 2
y4 = Abs(x >= 15.75) * (-4 * q0 / 216) * (x - 15.75) ^ 4
y5 = Abs(x >= 15.75) * (4 * q0 / 4860) * (x - 15.75) ^ 5
y6 = -13.5 * Ra
EItheta = Abs(x <= 24.75) * (y1 + y2 + y3 + y4 + y5 + y6)
Rotation = EItheta / EI
End Function
Function Deflection(x)
If x > 24.75 Then
Deflection = Deflection(24.75) + Rotation(24.75) * (x - 24.75)
Exit Function
End If
y1 = Abs(x >= 0) * (1 / 6) * Ra * (x - 0) ^ 3
y2 = Abs(x >= 9) * (1 / 6) * Rb * (x - 9) ^ 3
y3 = Abs(x >= 18) * (1 / 6) * Rc * (x - 18) ^ 3
y4 = Abs(x >= 15.75) * (-4 * q0 / 1080) * (x - 15.75) ^ 5
y5 = Abs(x >= 15.75) * (4 * q0 / 29160) * (x - 15.75) ^ 6
y6 = -13.5 * Ra * x
EIy = Abs(x <= 24.75) * (y1 + y2 + y3 + y4 + y5 + y6)
Deflection = EIy / EI
End Function