1 
jgs 
110 
% $Id$ 
2 
jgs 
121 
\section{Bending Beam under Uniform Load} 
3 
jgs 
110 
\label{BEAM CHAP} 
4 


\sectionauthor{Jannine Eisenmann}{} 
5 



6 


Given is a twodimension bending beam which is fixed in all directions 
7 


at the left end and free at the other, see \fig{BEAM FIG 1}. Furthermore the beam is loaded 
8 


with a uniform load $p$. 
9 



10 


\begin{figure} 
11 


% \centerline{\includegraphics[width=\figwidth]{DiffusionRes1}} 
12 


\caption{Loaded Beam.} 
13 


\label{BEAM FIG 1} 
14 


\end{figure} 
15 



16 



17 



18 


For emphasizing the displacement this picture is drawn (uebertrieben), 
19 


cause since we use the linear theory otherwise no displacements would 
20 


be visible. 
21 


There are two ways of solving this problem: on the one hand one can 
22 


use the differential equation for the deformed neutral fibre of the 
23 


beam. This classical differential equation is based on several simplifying 
24 


assumptions concerning the statics and kinematics of the problem. 
25 


However the results are known to be highly accurate at least for slender 
26 


beams with length to hight ratios $> 10$. Alternatively, in connection 
27 


with finite element based differential equation toolkits one may simply 
28 


consider the beam as a special case of a 2D or 3D elastic continuum 
29 


problem and solve the stress equilibrium equations combined with Hooke's 
30 


law for the specific boundary conditions of a beam. Both cases assume 
31 


isotropic and linear elastic material properties. 
32 



33 


The beam equation is easily solved analytically. The analytical solutions 
34 


can be used for benchmarkung finite element solutions. In section 
35 


1.2 we formluate a finite element code for general isotropic elasticity 
36 


problems including thin and deep beams under arbitrary loading conditiond 
37 


in 2D or 3D. 
38 



39 



40 


\section{2dimensional} 
41 


As the stress equilibrium equation is a partial differential equation 
42 


(PDE), we choose to use Finley to solve it, since Finley is a finite 
43 


element kernel library, that has been incorporated as a differential 
44 


equation solver into the python based numerical toolkit called escript. 
45 


We divided the beam into ten square elements of the size 1x1. Each 
46 


element consists of 8 nodes, which leads to a quadratic interpolation 
47 


of the model point displacements \\ 
48 



49 


The key ingredients is \textbf{Hooks Law}. We use Hooks Law because 
50 


we are dealing with \textbf{linear elastic material} \textbf{behaviour}. 
51 


We have \\ 
52 



53 



54 


$\sigma_{ik}=2G\left(\varepsilon_{ik}+\frac{\nu}{12\nu}\cdot e\cdot\delta_{ik}\right)$\hfill{}(1)\\ 
55 


where the engineering strain$\varepsilon_{ik}$is defined as: 
56 



57 


$\varepsilon_{ik}=\frac{1}{2}\cdot\left(u_{k,i}+u_{i,k}\right)$\hfill{}(2)\\ 
58 



59 



60 


with 
61 



62 


\begin{enumerate} 
63 


\item e= Volume strain = $\varepsilon_{kk}$ 
64 


\item $\delta_{ik}$= Kronecker symbol 
65 


\end{enumerate} 
66 


Inserting equation (2) in (1) (and with further mathematical conversions) 
67 


leads to the following partial differential equation :\\ 
68 



69 



70 


$\sigma_{ij}=\left[\mu\left(\delta_{ik}\delta_{jl}+\delta_{il}\delta_{jk}\right)+\lambda\left(\delta_{ij}\delta_{kl}\right)\right]u_{k,l}$\\ 
71 



72 



73 


For tension equilibrium we require:\\ 
74 



75 



76 


$\sigma_{ij,j}=0$\\ 
77 



78 



79 


which leads to the following expression:\\ 
80 



81 



82 


\[ 
83 


\left(\left[\mu\left(\delta_{ik}\delta_{jl}+\delta_{il}\delta_{jk}\right)+\lambda\left(\delta_{ij}\delta_{kl}\right)\right]u_{k,l}\right)_{,j}=0\] 
84 



85 



86 


with the natural boundary condition: 
87 



88 


\[ 
89 


n_{j}\sigma_{ij}=p_{i}\] 
90 


\\ 
91 


$p_{i}$ representing a uniform load on top of the beam. 
92 



93 


A Dirichlet Boundary condition is assumed on the left end of the beam 
94 


for which no displacements can occure.\\ 
95 


\\ 
96 


\includegraphics[% 
97 


width=0.60\linewidth,bb = 0 0 200 100, draft, type=eps]{/home/jeannine/sandbox/report/draws/dir_cond_beam.eps}\\ 
98 


This is described in the code with the setting a mask for the left 
99 


end and setting values to that mask: 
100 



101 


\begin{python} 
102 


q = xNodes{[}0{]}.whereZero(){*}{[}1.0,1.0{]} 
103 



104 


r = Vector({[}0.0, 0.0{]}, where = nodes) 
105 


\end{python} 
106 


The Finley template PDE reads:\\ 
107 



108 



109 


\[ 
110 


(A_{ijkl}u_{k,l})_{,j}(B_{ijk}u_{k})_{,j}+C_{ikl}u_{k,l}+D_{ik}u_{k}=X_{ij,j}+Y_{i}\] 
111 


\\ 
112 


with the natural boundary condition: 
113 



114 


\[ 
115 


n_{j}(A_{ijkl}u_{k,l}+B_{ijkl}u_{k})+d_{ik}u_{k}=n_{j}X_{ij}+y_{i}on\Gamma_{i}^{D}\] 
116 



117 



118 


Yields by comparing the coefficients : 
119 



120 


\begin{enumerate} 
121 


\item $A_{ijkl}$= $\left[\mu\left(\delta_{ik}\delta_{ij}+\delta_{jl}\delta_{il}\right)+\lambda\left(\delta_{ij}\delta_{kl}\right)\right]$ 
122 


\item $y_{i}$= $p_{i}$ 
123 


\item $u_{k}$= displacement $u$ 
124 


\end{enumerate} 
125 


$B_{ijk,}=C_{ikl}=D_{ik}=X_{ij}=Y_{i}=d_{ik}=0$\\ 
126 



127 



128 


Where 0 in the last line is taken as a scalar, vector or tensor, depending 
129 


on what the belonging coefficient is. 
130 



131 


These equations are the base for the \textbf{Finley Script}: 
132 



133 


\begin{python} 
134 


from ESyS import {*} 
135 


import Finley 
136 



137 



138 



139 


\#mu lamda 
140 



141 


def mu (E, nu): \#= shear modul G 
142 



143 


return E/(2{*}(1+nu)) 
144 



145 


def lamda (E, nu): 
146 



147 


return (nu{*}E)/((12{*}nu){*}(1+nu)) 
148 



149 



150 



151 


def main() 
152 



153 


\#material, beam PROPERTIES 
154 



155 


L = 10.0 \#length of beam {[}m{]} 
156 



157 


h = 1 \#height of beam {[}m{]} 
158 



159 


p = 1 \#outer uniform load {[}kN/m{]} 
160 



161 


E0 = 210000 \#Young's modulus {[}kN/m\textasciicircum{}2{]} 
162 



163 


nu = 0.4 \#Poisson ratio {[}{]} 
164 



165 



166 



167 


print \char`\"{}L=\char`\"{}, L 
168 



169 


print \char`\"{}h=\char`\"{}, h 
170 



171 


print \char`\"{}p=\char`\"{}, p 
172 



173 


print \char`\"{}E=\char`\"{}, E0 
174 



175 


print \char`\"{}nu=Poisson =\char`\"{}, nu 
176 



177 


print \char`\"{}mu=\char`\"{}, mu (E0,nu) 
178 



179 


print \char`\"{}lamda=\char`\"{}, lamda (E0,nu) 
180 



181 



182 



183 


\#SET MESH for FE 
184 



185 


mesh= Finley.Rectangle(n0=10 ,n1=1 ,order=2, l0=L, l1=h) 
186 



187 


nodes = mesh.Nodes() 
188 



189 


xNodes = nodes.getX() 
190 



191 


elements = mesh.Elements() 
192 



193 


faceElements = mesh.FaceElements() 
194 



195 


xFaceElements = faceElements.getX() 
196 



197 



198 



199 


\#DISPLACEMENT boundary 
200 



201 


q = xNodes{[}0{]}.whereZero(){*}{[}1.,1.0{]} \#setting the mask for r 
202 



203 


r = Vector({[}0.0, 0.0{]}, where = nodes) \#fixed end 
204 



205 



206 



207 


\#STRESS boundary 
208 



209 


ymask = xFaceElements{[}1{]}.whereEqualTo(h) 
210 



211 


y = Vector({[}0, p{]}, where = faceElements) 
212 



213 


y = y{*}ymask 
214 



215 



216 



217 


\#Finley coeff. 
218 



219 


A = Tensor4(0, where = elements) 
220 



221 


for i in range (2) : 
222 



223 


for j in range (2) : 
224 



225 


A{[}i,i,j,j{]}+= lamda (E0,nu) 
226 



227 


A{[}j,i,j,i{]}+= mu (E0,nu) 
228 



229 


A{[}j,i,i,j{]}+= mu (E0,nu) 
230 



231 



232 



233 


M,b = mesh.assemble(A=A, B=B, q=q, 
234 



235 


y=y,r=r,type=CSR,num\_equations=2) 
236 



237 


u= M.iterative(b, tolerance=1e8,iter\_max=50000, 
238 



239 


iterative\_method=PCG) 
240 



241 


print \char`\"{}u{[}0{]}=\char`\"{},u{[}0{]} 
242 



243 


print \char`\"{}u{[}1{]}=\char`\"{},u{[}1{]} 
244 



245 


main() 
246 


\end{python} 
247 


The finer the mesh the more exact is the solution. E.g. a 20x2 mesh 
248 


is more exact than a 10x1 mesh. 
249 


