/[escript]/trunk/doc/user/heatedblock.tex
ViewVC logotype

Annotation of /trunk/doc/user/heatedblock.tex

Parent Directory Parent Directory | Revision Log Revision Log


Revision 121 - (hide annotations)
Fri May 6 04:26:16 2005 UTC (14 years, 9 months ago) by jgs
Original Path: trunk/esys2/doc/user/beam.tex
File MIME type: application/x-tex
File size: 6606 byte(s)
Merge of development branch back to main trunk on 2005-05-06

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 two-dimension 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{2-dimensional}
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}{1-2\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)/((1-2{*}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=1e-8,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    

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26