The Analysis of Rectangular solid slabs by Finite Differences

I have a BASIC program which provides an elastic analysis of a rectangular slab/beam system with variable
section properties and any combination of supported nodes, the program also treats non-linear cases where
lift-off occurs at supports which allows foundation rafts and slabs with unrestrained corners to be analysed.
This sheet is a trial to see how the analysis might transfer to Mathcad format and, in the first instance, is limited
to the analysis of rectangular slabs using a square mesh (the program allows a rectangular mesh), Poisson's
ratio is zero and beam supports are not yet included.
The loads are specified as point loads at any node, by apportioning distributed loads appropriately any
combination of point, whole- or part-UD or hydrostatic loads can be accepted. The stiffness matrix is assembled
from the finite-difference operator patterns and then inverted, the deflexions are obtained by multiplying the
inverted matrix by the load vector. The stress resultants are calculated from the deflexions using the appropriate
finite difference formulae.

Finite differences have been pretty well been superseded by finite elements for slab design but, when the
plan-form is rectangular, they still have a lot to offer by way of simplicity and accuracy.

1 2 3 N-2 N-1 N

h
N+1 N+2 N+3 2N-2 2N-1 2N

2N+1 2N+2 2N+3 3N-2 3N-1 3N

No. of longit.
grid lines, N

No. of transverse
grid lines, M

(M-3)*N+1 (M-3)*N+2 (M-3)*N+3 (M-2)*N-2 (M-2)*N-1 (M-2)*N

(M-2)*N+1 (M-2)*N+2 (M-2)*N+3 (M-1)*N-2 (M-1)*N-1 (M-1)*N

(M-1)*N+1 (M-1)*N+2 (M-1)*N+3 M*N-2 M*N-1 M*N

TYPICAL RECTANGULAR SLAB

Enter the data:-

width
Width of slab (mm)...................................................... width 7000 Number of longit grid lines, N 1
Grid line spacing (mm)...................................................... h 500 h
Number of transverse grid lines.......................................... M 15 Length of slab, (mm), length h ( M 1)
Depth of slab (mm)........................................................... d 200
3
Young's modulus (kN/mm^2)............................................. E 28 d
distributed inertia, D E
12

Note:- N and M not to be less than 5 length N

1
 The stiffness matrix is assembled from the appropriate finite difference patterns for a rectangular slab with all edges
free/free in the collapsed section below. For information on the method see Structural Analysis, Ghali & Neville,
Chapman & Hall, or some similar text

Infill below leading diagonal (at the same time correcting for the half-value leading diagonal coefficient
input above).

T
A A A

Supply some supports, initially, say a line support on the left hand side and two point supports near the other
two corners:- (If the slab edges are not tied down, lift-off from the line support can be engineered by trial and
error, after the first run any support with a tensile reaction is freed. After several runs, repeating this procedure
this results in only the seven central wall nodes being supported, as follows),

Dimension the array: Support 0 Set the supported nodes to a high number:
MN
12
i 4N 1 5N 1 ( M 5) N 1 Support 10
i
12
i N 4 Support 10
i
12
i ( M N) 4 Support 10
i

Add to leading diagonal of A matrix: i 1 M N A A Support

i i i i i

Then part of the stiffness matrix looks like:-

2
Supply a load vector, say a UDL at the rate of W kN/sq.m. W 10 W length width
Panel_load
6
( M 1) ( N 1) 1 10

initially set the whole load matrix so that each node is loaded by Panel_load i 1 2 M N Load Panel_load
i
then modify the edge and corner nodes:-
Load Load
i i
i 1 2 N Load i N 2 N M N Load
i 2 i 2

Load
i
i 1 N 1 ( M 1) N 1 Load Load
i 2 i
i ( M 1) N 1 ( M 1) N 2 M N Load
i 2

In the BASIC version of the program the heavily banded stiffness matrix is solved using Choleski decomposition.
With Mathcad it seems much simpler to invert the matrix

2
1 h
Solve the problem: Deflexions A Load
D

Deflns p0
for i 1 N M
Rearrange deflexions into a rectangular array, (for plotting).
pp1

1
q ceil p
N
r p ( q 1) N
D Deflexions
q r p

B Deflns

3
 A table containing the deflexions at the nodes, in mm units, is as follows:-

Deflns

and a graph of deflexions appears thus:-

DEFLEXIONS

Deflns

(Notice how the unrestrained rear corners pick themselves up off the end support wall)

4
The support reactions are found by multiplying the support stiffnesses by the deflexions at the nodes and the
appropriate constant, thus:

D
i 1 M N Rn Deflexions Support
i i i 2
h

Check that sum of reactions is equal to sum of applied loads:

Sum of reactions= Sum_reactions

Rn Sum_reactions
(If these two are equal we have
a balance of vertical forces,
Sum of Loads= Sum_loads
Load Sum_loads if not, we have an error!.)

The two column reactions are:- Rn and Rn

N 4 ( MN) 4

Rearrange reactions into a rectangular

array, for plotting: Rect( Rn) p0
for i 1 N M
pp1

1
q ceil p
N
r p ( q 1) N
A plot of the reactions is as follows:-
D Rn
r q p

D
REACTIONS at SUPPORTED NODES

Rect( Rn)

5
 Calculate some bending moments:
Mx p0
for i 1 M
pp1
q1
for j 1 N 2
qq1
M1
p q pq1 2Bpq Bpq1
B

M1 0
M N

D
M1 M1
2
h
Table of Mx values in kN:- M1

Mx

Bending Moment Mx

Mx

6
 My p1
for i 1 M 2
pp1
q0
for j 1 N
qq1
M1
p q p1q 2Bpq Bp1q
B

M1 0
M N

D
M1 M1
2
h
M1
My in kN. :-

My

Bending Moment My

My

7
 Calculate the twisting moments
Mxy p1
for i 2 M 1
pp1
q1
for j 2 N 1
qq1

M1
Bp1q1 Bp1q1 Bp1q1 Bp1q1
p q 4
M1
M N MN BM N1 BM1N1 BM1N
B

B
M 1 2
M1 B B B
M 1 M 2 M 1 M 1 1
M1 B B B B
1 1 2 2 2 1 1 1 1 2
M1 B B B B
1 N 2 N 2 N 1 1 N 1 1 N
p1
for i 2 N 1
pp1
M1
1 p 2p1 B2p1 B1p1 B1p1.5
B

B
M 1 p 1
M1 B B B .5
M p M p 1 M p 1 M 1 p 1

p1
for i 2 3 M 1
pp1
M1
p 1 p12 Bp11 Bp11 Bp12.5
B

B
p 1 N
M1 B B B .5
p N p 1 N p 1 N 1 p 1 N 1

D
Table of Mxy, in kN:- M1 M1
2
h
M1

Mxy

8
and a plot of Mxy:-

Mxy