bridge.html

Deformation of a Beam

This worksheet shows the calculation of the deformation of a simply supported beam

The author expects that this worksheet will only be used for teaching and educational purposes and not for commerical profit without contacting the author for a licenced agreement.

> restart;

[Maple OLE 2.0 Object]

Introduction

Consider an elastic both side simply supported beam. This beam can be for example a bridge. It is loaded by a moving load. This might be a car driving over the bridge. The situation is considered to be quasi static. So there are no dynamic effects respected.

NOTE: Influence of shearing is not considered.

The following variables are used:

F: Force

L: Length of the beam

EI: Product of the elastic modulus and the moment of inertia

xF: Position of the force ralated to the left end of the beam

x: Position where the deflection of the beam is calculated

Cross section of the beam:

[Maple OLE 2.0 Object]

Further variables are defined at there first appearence.

All units are m and N

Initialisation

> with(plots):with(plottools):

Warning, the name changecoords has been redefined

Length of the beam

> L:=20:

Force

> F:=1*10**4:

Number of steps where the load is positioned

> nx:=25;

Length of the steps between the positions of the load

> dL:=L/nx;

Dimension of the cross section

> a:=1/2:b:=1/2:d:=1:e:=1:f:=1/2:

Number of parts of the cross section

> n:=3;

The Cross section is for later calculations divided in three parts, the plate and the two webs.

The dimensions of the parts of the cross section are collected in sets.

> width:=[a+b+d+b+a,b,b];height:=[f,e,e];

boundaries of the cross section in z-direction

> zmin:=0;zmax:=e+f;

Collect the coordinats centers of gravity of every part in sets (arbitrary with respect to the left upper corner, other points are also possible)

> cogy:=[a+b+d/2,a+b/2,a+b+d+b/2];cogz:=[f/2,f+e/2,f+e/2];

The material in this example should be concret, so the lastic modulus is

> Emod:=3*10**10;

Procedures

All the main calculations and all graphics are programmed in form of procedures.

Procedures for the Calculation

Area of the Cross Section

of a beam which is compounded by rectangle pieces

> crosssec:=proc(n,B,H)

n is the number of the parts of the cross section. In this example all part have to be rectangles.

B is the width and H the height of every rectangle.

> local i,i1,A;

A[i] is the area of the part i of the cross section.

> for i from 1 by 1 to n do

> A[i]:=B[i]*H[i]:

> od:

Cross section

> A:=sum(A[i1],i1=1..n):

End of the procedure

> end:

Center of Gravity

of a cross section which is compounded by rectangle pieces

> cog:=proc(n,B,H,ya,za)

n is the number of the parts of the cross section. In this example all part have to be rectangles.

B is the width and H the height of every rectangle.

ya and za are the coordinates of the center of gravity of every rectangle with respect to a arbitrary but fixed point P.

> local i,i1,yc,zc,A,Ay,Az,AA,AAy,AAz,yz;

yz=[yc,zc]

yc and zc are the coordinates of the center of gravity of the total cross section with respect to the point P.

A[i] is the area of the part i of the cross section.

yc=sum(ya[i]*A[i])/sum(A[i]), zc=sum(za[i]*A[i])/sum(A[i]).

> for i from 1 by 1 to n do

> A[i]:=B[i]*H[i]:

> Ay[i]:=A[i]*ya[i]:

> Az[i]:=A[i]*za[i]:

> od:

> AA:=sum(A[i1],i1=1..n):

> AAy:=sum(Ay[i1],i1=1..n):

> AAz:=sum(Az[i1],i1=1..n):

> yc:=AAy/AA:

> zc:=AAz/AA:

> yz:=[yc,zc];

End of the procedure

> end:

Moment of Inertia of a Beam

For calculating the moment of inertia of the beam a procedure is used.

> inertia:=proc(n,B,H,ya,za,yzc)

n is the number of the parts of the cross section. In this example all part have to be rectangles.

B is the width and H the height of every rectangle.

ya and za are the coordinates of the center of gravity of every rectangle with respect to a arbitrary but fixed point P.

yzc=[yc,zc] is the set of the coordinates of the center of gravity of the whole cross section with respect to point P.

> local i,i1,yc,zc,yz,A,Ay,Az,IIy,IIz,Iy,Iz,IAy,IAz,AA,AAy,AAz,II;

center of gravity

> yc:=yzc[1];

> zc:=yzc[2];

A[i] is the area of the part i of the cross section.

> for i from 1 by 1 to n do

> A[i]:=B[i]*H[i];

> od:

Iy[i]=(B[i]*H[i]^3)/12 is the moment of inertia of every rectangle with respect to the y-axis, Iz[i]=(H[i]*B[i]^3)/12 is the moment of inertia of every rectangle with respect to the z-axis.

IAy[i]=A[i]*dz[i]^2 and IAz[i]=A[i]*dy[i]^2 are the additional part of the moment of inertia which results by the influence of the difference of the center of gravity of every part of the cross section respectively the total center of gravity, dy[i]=ya[i]-yc, dz[i]=za[i]-zc.

> for i from 1 by 1 to n do

> Iy[i]:= B[i]*H[i]**3/12:

> Iz[i]:= H[i]*B[i]**3/12:

> IAy[i]:=A[i]*(za[i]-zc)**2:

> IAz[i]:=A[i]*(ya[i]-yc)**2:

> od:

The moment of inertia of the total cross section is the sum of all moments of inertia of all parts of the cross section plus the additional part which result of the difference of the centers of gravity.

> IIy:=sum(Iy[i1]+IAy[i1],i1=1..n):

> IIz:=sum(Iz[i1]+IAz[i1],i1=1..n):

> II:=[IIy,IIz]:

End of the procedure

> end:

Bending Moment of a Beam

which is both sided simply supported and loaded by single load

> moment:=proc(x,xF,L,F)

xF is the position of the load with respect to the left end of the beam

L ist the length of the beam

F is the load

x is the independent variable

M(x) means the bending moment at point x

> local M;

The bending moment left from the position of the load is given by

F/L*(L-xF)*x,

on the right side the bending moment is given by

F/L*(L-xF)*x-F*(x-xF)

> M(x):=piecewise(x<=xF,F/L*(L-xF)*x,x>xF,F/L*(L-xF)*x-F*(x-xF));

End of the procedure

> end:

Deformation of a Beam

which is both sided simply supported and loaded by single load

(shear deformations are neglectable in this example)

> deform:=proc(x,xF,L,EI,M)

xF is the position of the load with respect to the left end of the beam

L ist the length of the beam

EI is the product of the elastic modulus and the moment of inertia

x is the independent variable

M means the bending moment, depending on x

R(x) means the gradient at point x

> local R,w,eq1,eq2,RB1,RB2,lsg,dfrm;

The problem can be written as a system of differential equations.

The derivation of the gradient mutiplied by the elastic modulus and the moment of inertia yields the bending moment.

The derivation of the deflection yields the gradient.

> eq1:=diff(R(x),x)*EI=M;

> eq2:=diff(w(x),x)=R(x);

There is no deflection on both ends of the beam at the points of simply supporting.

> RB1:=w(0)=0;

> RB2:=w(L)=0;

The solution of the system of differential equation eq1 and eq2 with respect to the Conditions RB1 and RB2 yields the deflection of the beam

> lsg:=dsolve({eq1,eq2,RB1,RB2},{R(x),w(x)});

> assign(lsg);

> dfrm:=[(w(x)),(R(x))];

End of the procedure

> end:

Normal Stress of a Beam

which is both sided simply supported and loaded by single load

(shear deformations are neglectable in this example)

> stress:=proc(x,z,xF,L,Iy,M)

xF is the position of the load with respect to the left end of the beam

L ist the length of the beam

x is the independent variable for the position of the regarded cross section

z is the independent variable for the position within the cross section

M means the bending moment, depending on x

Iy is the moment of inertia

> local sigma,tau;

The normal stress is given by sigma=M/Iy*z.

> sigma:=M/Iy*z;

End of the procedure

> end:

Procedures for the Plots

Rotation and Following Translation of the Coordinate of a Point

This procedure is used for the later graphical representation of the solution

> transrot:=proc(XY,dxy,phi)

XY is the set [X,Y] of the originally Coordinates of the point

dxy is the set [dx,dy] of the translation of the point in x- and y-direction

phi is the angular of the rotation

xy is the set [x,y] of the new Coordinates of the point

> local X,Y,dx,dy,x,y,xy;

> X:=XY[1];Y:=XY[2];

> dx:=dxy[1];dy:=dxy[2];

> x:= X*cos(phi)+Y*sin(phi)+dx;

> y:=-X*sin(phi)+Y*cos(phi)+dy;

> xy:=[x,y];

> end:

Drawing a Car in a x-y-Coordinate System

This procedure is used to draw a car which symbolizes the moving load in the later graphical representation of the solution.

> carplot:=proc(dx,dy,phi)

dx and dy are the coordinates where the car must be drawn. The car is rotated in the Figure by the angular phi.

> local XY,x,y,car,wheel1,wheel2,win1,win2,win3,win4,P,i;

> x[1]:=-1; y[1]:= 0.25;

> x[2]:= 1; y[2]:= 0.25;

> x[3]:= -1.5; y[3]:= 0.5;

> x[4]:= 1.5; y[4]:= 0.5;

> x[5]:= 1.3; y[5]:= 1;

> x[6]:= 0.5; y[6]:= 1;

> x[7]:= 0.3; y[7]:= 1.5;

> x[8]:= -1.3; y[8]:= 1.5;

> x[9]:= -1.5; y[9]:= 0.5;

> x[10]:=-1.4; y[10]:=1;

> x[11]:=0.3; y[11]:=1;

> x[12]:=0.5; y[12]:=0.5;

> x[13]:=-0.5; y[13]:=1.5;

> x[14]:=-0.5; y[14]:=0.5;

Here the procedure transrot is used to displace the points to the actual position of the load and of the corresponding deflection of the beam.

> for i from 1 by 1 to 14 do

> XY:=transrot([x[i],y[i]],[dx,dy],phi);

> x[i]:=XY[1]:

> y[i]:=XY[2]:

> od;

Drawing the car by using the points (x[i],y[i])

> win1:=curve([[x[6],y[6]],[x[10],y[10]]],color=blue,thickness=2);

> win2:=curve([[x[7],y[7]],[x[11],y[11]]],color=blue,thickness=2);

> win3:=curve([[x[6],y[6]],[x[12],y[12]]],color=blue,thickness=2);

> win4:=curve([[x[13],y[13]],[x[14],y[14]]],color=blue,thickness=2);

> wheel1:=circle([x[1],y[1]],0.25,color=black,thickness=2);

> wheel2:=circle([x[2],y[2]],0.25,color=black,thickness=2);

> car:=curve([seq([x[i],y[i]],i=3..9)],color=blue,thickness=2);

> P:=display({car,wheel1,wheel2,win1,win2,win3,win4});

> end:

Draw the Supports of the Beam

No further comment necessary

> support:=proc(C,h)

C is the point which has to be supported, h is the dimension of the drawing of the symbol

> local P0,P1,P2,P3,PL1,PL2,P;

> P0:=[C[1]+h*sqrt(2)/4,C[2]-h*sqrt(2)/4];

> P1:=[C[1]+h,C[2]-h];

> P2:=[C[1]-h,C[2]-h];

> P3:=[C[1]-h*sqrt(2)/4,C[2]-h*sqrt(2)/4];

> PL1:=circle(C,h/2,color=magenta,thickness=2);

> PL2:=curve([P0,P1,P2,P3],color=magenta,thickness=2);

> P:=display([PL1,PL2]);

> end:

Draw the Bridge Head on Both Sides

by using the procedure support

> bridgeplot:=proc(L)

L is the length of the bridge

> local x,y,Q,i,P1,P2,P3,P4,P;

> P1:=support([0,0],L/50);

> P2:=support([L,0],L/50);

> P3:=curve([[-L/10,0],[0,0]],color=black,thickness=2);

> P4:=curve([[L,0],[L+L/10,0]],color=black,thickness=2);

> P:=display([P1,P2,P3,P4]);

> end:

Plot the Deformation

In this procedure the deflection of the beam will be plotted. To show the current position of the load there are a little car an an arrow drawn.

> plotsol:=proc(x,xF,L,w,R)

xF is the position of the load with respect to the left end of the beam

L ist the length of the beam

x is the independent variable

w is the deflection of the beam between the points of support

R is the gradient of the beam

W(x) is the deflection of the beam but 50 000 times inflated and including a small piece at both ends (only needed for better plotting)

V(x) is the corresponding gradient

Pfeil is the plot of an arrow

Car is the plot of the car

Bridge plots the supports

Beam plots the deformed beam

P is the plot of the total figure

> local W,V,Pfeil,Car,Bridge,Beam,P;

Inflate 50 000 times and append a smal piece at every end of the beam

> W(x):=piecewise(x<0,0,x>=0 and x<=L,w*5*10**4,x>L,0);

Gradient of the beam, also inflated 50 000 times. This is needed that the car in the drawing stands with all wheels on the deflected beam

> V(x):=piecewise(x<0,0,x>=0 and x<=L,-arctan(R*5*10**4),x>L,0);

Plot the car at the corrent point of the load by using the procedure carplot

> Car:=carplot(xF,subs(x=xF,W(x)),subs(x=xF,V(x)));

Plot the supports

> Bridge:=bridgeplot(L);

Plot an arrow at the corrent point of the load

> Pfeil:=arrow([xF,subs(x=xF,W(x)+2)], [xF,subs(x=xF,W(x))], L/200, L/50, L/100, color=red):

> Beam:=plot(W(x),x=-L/10..L+L/10,color=black,thickness=2);

Display the whole plot

> P:=display(Bridge,Car,Pfeil,Beam);

End of the procedure

> end:

Plot the Normal Stress in the Beam

In this procedure the normal forces in the beam will be plotted. To show the current position of the load there is an arrow drawn.

> stressplot:=proc(x,z,xF,L,sigma,zmin,zmax)

xF is the position of the load with respect to the left end of the beam

L ist the length of the beam

x is the independent variable for the position of the regarded cross section

z is the independent variable for the position within the cross section whithin the boundaries zmin and zmax

Pfeil is the plot of an arrow

Stress plots the stress

P is the plot of the total figure

> local S,Pfeil,Stress,P;

Plot an arrow at the corrent point of the load

> Pfeil:=arrow([xF,-zmax-(zmax-zmin)/2], [xF,-zmax], L/200, L/50, L/100, color=red):

Plot the normal stress

> S(x,z):=sigma;

> Stress:=contourplot(S(x,z),x=0..L,z=-zmax..-zmin,grid=[50,10],filled=true,contours=5,coloring=[red,green]);

Display the whole plot

> P:=display(Pfeil,Stress);

End of the procedure

> end:

Execution

Calculation

Calculate the center of gravity by using the procedure cog.

> yzc:=cog(n,width,height,cogy,cogz);

Boundaries of the cross section in z-direction in relation to the center of gravity

> zmins:=zmin-yzc[2];zmaxs:=zmax-yzc[2];

Calculate the moment of inertia by using the procedure inertia.

In this example only the moment of inertia with respect to the horizontal axis is necessary.

> Iy:=inertia(n,width,height,cogy,cogz,yzc)[1];

Area of the Cross section of the beam

> A:=crosssec(n,width,height);

For the further calculation the product of elastic modulus and moment of inertia is interresting.

> EI:=Emod*Iy;

Calculate the bending moment depending on x and xF.

> M:=convert(moment(x,xF,L,F),Heaviside):

Calculate the deformation of the beam depending on x and xF

> dfrm:=deform(x,xF,L,EI,M):

Separate the deflection

> w:=eval(dfrm[1]):

Separate the gradient

> R:=eval(dfrm[2]):

Calculate the normal stess in the beam depending on x and xF.

> sigma:=stress(x,z,xF,L,Iy,M):

Create the Animation

Calculate the line of deflection for every position of the load and return the solution in form of a plot by using the procedure plotsol . Thereby the load is moving over the bridge like a driving car. The deflection of the bridge is 50 000 times inflated. Note: there are no real dynamic effects regarded. When the car reaches the right end of the bridge, it drives back to the left end. The result is shown by an little animation.

> for i from 0 by 1 to 2*nx do

> if (i<=nx) then X:=eval(L*i/nx) else X:=eval(L-(i-nx)/nx*L) fi:

> WP:=subs(xF=X,w):

> P[i]:=display([bridgeplot(L),subs(xF=X,plotsol(x,X,L,WP,subs(xF=X,R)))]):

> od:

Calculate the normal stress for every position of the load and return the solution in form of a plot by using the procedure stressplot .

> for i from 0 by 1 to nx do

> X:=eval(L*i/nx):

> SP:=subs(xF=X,stress(x,z,xF,L,Iy,M)):

> S[i]:=display([subs(xF=X,stressplot(x,z,X,L,SP,zmaxs,zmins))]):

> od:

View the Result

Deformation

> display([seq(P[i],i=0..2*nx)],insequence=true,scaling=constrained,axes=none);

[Maple Plot]

The maximum deflection of the beam is reached when the load is positioned in the middle of the beam. The maximum deflected point is then the center of the beam. The maximum deflection is thus

> fmax:=evalf(subs(x=L/2,deform(x,L/2,L,EI,moment(x,L/2,L,F))[1]));

Normal Stress

> display([seq(S[i],i=0..nx)],insequence=true,scaling=unconstrained);

[Maple Plot]

The stress for the situation when the load is positioned in the middle of the beam at the center of the beam is

> Smax:=evalf(subs(x=L/2,stress(x,zmaxs,L/2,L,Iy,moment(x,L/2,L,F))));

>

Disclaimer: While every effort has been made to validate the solutions in this worksheet, Waterloo Maple Inc. and the contributors are not responsible for any errors contained and are not liable for any damages resulting from the use of this material.