1.10 Least squares

Back

In FE analysis, the movements of a structure are constrained, because the

structure is only allowed to assume shapes that can be expressed as piecewise

linear, piecewise quadratic, or similar shape functions.

We then find that the FE solution is that deformation of the structure

which among all remaining movements renders the potential energy a minimum.

This is equivalent to the statement that the associated load case ph is

3 From now on we simply say scalar product instead of L2-scalar product.

34 1 What are finite elements?

Fig. 1.26. Best fit with a

straight line

work-equivalent to the original load case p. Finally the distance (in terms of

energy) between the exact solution u and the FE solution uh is a minimum:

a(e, e) = a(u − uh, u − uh) ≤ a(u − vh, u − vh) vh ∈ Vh . (1.114)

No other function vh ∈ Vh comes closer to u in this sense than uh. But a

minimum in energy means (in the case of a beam for example)

a(e, e) =

_ l

0

(M −Mh)2

EI

dx → minimum, (1.115)

that the least-squares error of the internal actions is a minimum.

Least squares is a concept of numerical analysis. When a straight line

Mh(x) = ax + b is drawn through a number of points (Fig. 1.26)

ax1 + b = M(x1) ,

ax2 + b = M(x2) , (1.116)

. . . = . . . ,

axn + b = M(xn) ,

as to minimize the sum of squared errors

F =

_n

i=1

(Mh(xi) −M(xi))2 → minimum, (1.117)

∂ F

∂ a

= 0,

∂ F

∂ b

= 0, (1.118)

we obtain a least-squares solution. Solving (1.118) is equivalent to solving the

2 × 2 system of equations

AT

(2×n)A(n×2)

_

a

b

_

= AT

(2×n)m(n) , (1.119)

the normal equations, where A is the coefficient matrix of (1.116) and the

vector m = [M(x1),M(x2) . . .] is the right-hand side of (1.116).

In FE analysis, the bending moment is a function

Mh(x) =

_

i

wiMi(x), Mi(x) = −EI ϕ

__

i (x) , (1.120)

1.10 Least squares 35

and because the FE solution is minimized with respect to all points 0 < x < l

of a beam, the sum is replaced by an integral (see (1.115)).

The equation Ku = f is the associated normal equation [232]. This equation

is obtained if the infinitely many equations Mh(x) = M(x) (there are

infinitely many points x1, x2, . . . in the interval [0, l]) basically a matrix A∞×n

is multiplied by the transposed matrix AT

n×∞ from the left, and the diagonal

matrix C∞×∞ with weights Cii = EI is placed in between:

AT(

n×∞) C(∞×∞)A(∞×n) = K(n×n) . (1.121)

Weighted least squares

The term EI comes from the strain energy product

a(w, ˆ w) =

_ l

0

M ˆM

EI

dx , (1.122)

because in FE analysis the bending moments are scaled in such a way that in

the weighted least-squares sense the discrepancies between Mh and the exact

bending moment M are minimized

F :=

_ l

0

(M −Mh)2

EI

dx =

_ l

0

(M −Mh) (κκh) dx → min. (1.123)

which means that the error in the bending moments is multiplied by the error

in the curvature, and that this product is minimized.

Of course if the bending stiffness EI is constant, then there is no difference

between least squares and weighted least squares.

Global and local

In least squares the global picture—the error over the whole interval [0, l]—is

studied. But the algorithm also incorporates a local match, as we will see in

the following. To simplify the derivation it is assumed that EI = 1, and the

short-hand notation for integrals (see Chap. 7) is used.

For the integral

F = (M −Mh,M −Mh) = (M,M) − 2 (M,Mh) + (Mh,Mh)

= (M,M) − 2

_

i

(M,Mi)wi +

_

i,j

(Mi,Mj)wi wj (1.124)

to attain a minimum at w it is necessary that the gradient of F vanishes at

this point with respect to the nodal values wi:

∂F

∂wi

= 2

_

j

(Mi,Mj)wj − 2 (M,Mi) = 2(Mh,Mi) − 2(M,Mi) = 0.

(1.125)

36 1 What are finite elements?

Fig. 1.27. Cantilever beam and FE load case. At the Gauss points the error in the

bending moment is zero

Hence the solution of (1.123) is equivalent to making the error M − Mh orthogonal

to the functions Mi, with EI reinserted

_ l

0

(M −Mh)Mi

EI

dx = 0, i= 1, 2 . . . (local match) . (1.126)

Test range

The test range for the local match between M and Mh is not a single element,

but the support of the individual hat function Mi. (This is what the bending

moments of the unit deflections ϕi look like.) Outside of the support Ωi the

hat function is zero, and therefore the local test is an integral over the interval

Ωi = [xi−1, xi+1]:

_ l

0

(M −Mh)Mi

EI

dx =

_ xi+1

xi−1

(M −Mh)Mi

EI

dx = 0. (1.127)

In the case of a continuous beam the test range of an individual hat function

Mi consists of two consecutive elements. And the positive and negative errors

M−Mh must be distributed in such a way that in the weighted average sense

the error over any two neighboring elements disappears, see Fig. 1.27.

1.11 Distance inside = distance outside 37