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

Revision 737 - (hide annotations)
Wed Jun 14 04:00:41 2006 UTC (13 years, 8 months ago) by ksteube
File MIME type: application/x-tex
File size: 7906 byte(s)
Added a caution to be sure to set OMP_NUM_THREADS,
if unset you will get too many threads.

Also fixed a few spelling mistakes.


 1 jgs 110 % $Id$ 2 gross 625 % 3 % Copyright © 2006 by ACcESS MNRF 4 % \url{http://www.access.edu.au 5 % Primary Business: Queensland, Australia. 6 % Licensed under the Open Software License version 3.0 7 8 % 9 gross 578 \section{Elastic Deformation} 10 \label{ELASTIC CHAP} 11 In this section we want to discuss the deformation of a linear elastic body caused by expansion through a heat distribution. We want 12 to displacement field $u\hackscore{i}$ which solves the momentum equation 13 \index{momentum equation}: 14 gross 579 \begin{eqnarray}\label{HEATEDBLOCK general problem} 15 gross 578 - \sigma\hackscore{ij,j}=0 16 \end{eqnarray} 17 where the stress $\sigma$ is given by 18 gross 579 \begin{eqnarray}\label{HEATEDBLOCK linear elastic} 19 gross 578 \sigma\hackscore{ij}= \lambda u\hackscore{k,k} \delta\hackscore{ij} + \mu ( u\hackscore{i,j} + u\hackscore{j,i}) 20 - (\lambda+\frac{2}{3} \mu) \; \alpha \; (T-T\hackscore{ref})\delta\hackscore{ij} \;. 21 \end{eqnarray} 22 In this formula $\lambda$ and $\mu$ are the Lame coefficients, $\alpha$ is the 23 temperature expansion coefficient, $T$ is the temperature distribution and $T_{ref}$ a reference temperature. Note that 24 gross 579 \eqn{HEATEDBLOCK general problem} is similar to eqn{WAVE general problem} introduced in section~\Sec{WAVE CHAP} but the 25 gross 578 inertia term $\rho u\hackscore{i,tt}$ has been dropped as we assume a static scenario here. Moreover, in 26 comparison to the \eqn{WAVE stress} 27 gross 579 definition of stress $\sigma$ in \eqn{HEATEDBLOCK linear elastic} an extra term is introduced 28 ksteube 737 to bring in stress due to volume changes through temperature dependent expansion. 29 gross 568 30 gross 578 Our domain is the unit cube 31 gross 593 \begin{eqnarray} \label{HEATEDBLOCK natural location} 32 gross 578 \Omega=\{(x\hackscore{i} | 0 \le x\hackscore{i} \le 1 \} 33 \end{eqnarray} 34 On the boundary the normal stress component is set to zero 35 gross 579 \begin{eqnarray} \label{HEATEDBLOCK natural} 36 gross 578 \sigma\hackscore{ij}n\hackscore{j}=0 37 \end{eqnarray} 38 and on the face with $x\hackscore{i}=0$ we set the $i$-th component of the displacement to $0$ 39 gross 579 \begin{eqnarray} \label{HEATEDBLOCK constraint} 40 gross 578 u\hackscore{i}(x)=0 & \mbox{ where } & x\hackscore{i}=0 \; 41 \end{eqnarray} 42 For the temperature distribution we use 43 gross 579 \begin{eqnarray} \label{HEATEDBLOCK temperature} 44 T(x)= T\hackscore{0} e^{-\beta \|x-x^{c}\|}; 45 gross 578 \end{eqnarray} 46 gross 579 with a given positive constant $\beta$ and location $x^{c}$ in the domain. Later in \Sec{MODELFRAME} we will use 47 $T$ from an time-dependent temperature diffusion problem as discussed in \Sec{DIFFUSION CHAP}. 48 gross 578 49 gross 579 When we insert~\eqn{HEATEDBLOCK linear elastic} we get a second oder system of linear PDEs for the displacements $u$ which is called 50 gross 578 the Lame equation\index{Lame equation}. We want to solve 51 this using the \LinearPDE class to this. For a system of PDEs and a solution with several components the \LinearPDE class 52 takes PDEs of the form 53 gross 625 \begin{equation}\label{LINEARPDE.SYSTEM.1 TUTORIAL} 54 -A\hackscore{ijkl} u\hackscore{k,l}){,j}=-X\hackscore{ij,j} \; . 55 gross 568 \end{equation} 56 gross 625 $A$ is of \RankFour and $X$ is of \RankTwo. We show here the coefficients relevant 57 for the we traying to solve. The full form is given in~\eqn{LINEARPDE.SYSTEM.1}. 58 gross 568 The natural boundary conditions \index{boundary condition!natural} take the form: 59 gross 625 \begin{equation}\label{LINEARPDE.SYSTEM.2 TUTORIAL} 60 n\hackscore{j} A\hackscore{ijkl} u\hackscore{k,l}=n\hackscore{j}X\hackscore{ij} \;. 61 gross 568 \end{equation} 62 gross 625 Constraints \index{constraint} take the form 63 \begin{equation}\label{LINEARPDE.SYSTEM.3 TUTORIAL} 64 gross 568 u\hackscore{i}=r\hackscore{i} \mbox{ where } q\hackscore{i}>0 65 \end{equation} 66 gross 578 $r$ and $q$ are each \RankOne. 67 gross 625 We can easily identify the coefficients in~\eqn{LINEARPDE.SYSTEM.1 TUTORIAL}: 68 gross 578 \begin{eqnarray}\label{LINEARPDE ELASTIC COEFFICIENTS} 69 A\hackscore{ijkl}=\lambda \delta\hackscore{ij} \delta\hackscore{kl} + \mu ( 70 +\delta\hackscore{ik} \delta\hackscore{jl} 71 \delta\hackscore{il} \delta\hackscore{jk}) \\ 72 X\hackscore{ij}=(\lambda+\frac{2}{3} \mu) \; \alpha \; (T-T\hackscore{ref})\delta\hackscore{ij} \\ 73 \end{eqnarray} 74 The characteristic function $q$ defining the locations and components where constraints are set is given by: 75 gross 579 \begin{equation}\label{HEATEDBLOCK MASK} 76 gross 578 q\hackscore{i}(x)=\left\{ 77 \begin{array}{cl} 78 1 & x\hackscore{i}=0 \\ 79 0 & \mbox{otherwise} \\ 80 \end{array} 81 \right. 82 \end{equation} 83 Under the assumption that $\lambda$, $\mu$, $\beta$ and $T\hackscore{ref}$ 84 are constant setting $Y\hackscore{i}=\lambda+\frac{2}{3} \mu) \; \alpha \; T\hackscore{i}$ seems to be also possible. However, 85 this choice would lead to a different natural boundary condition which does not set the normal stress component as defined 86 gross 579 in~\eqn{HEATEDBLOCK linear elastic} to zero. 87 gross 578 88 gross 625 Analogously to concept of symmetry for a single PDE, we call the PDE defined by~\eqn{LINEARPDE.SYSTEM.1 TUTORIAL} symmetric if 89 gross 578 \index{symmetric PDE} 90 gross 625 \begin{eqnarray}\label{LINEARPDE.SYSTEM.SYMMETRY TUTORIAL} 91 gross 578 A\hackscore{ijkl} =A\hackscore{klij} \\ 92 gross 568 \end{eqnarray} 93 gross 578 Note that different from the scalar case now the coefficients $D$ and $d$ have to be inspected. It is easy to see that 94 gross 579 the Lame equation in fact is symmetric. The \LinearPDE class is notified by this fact by calling its \method{setSymmetryOn} method. 95 gross 568 96 gross 579 After we have solved the Lame equation we want to analyse the actual stress distribution. Typically the von--Mises stress\index{von--Mises stress} defined by 97 \begin{equation} 98 \sigma\hackscore{mises} = \sqrt{ 99 \frac{1}{6} ((\sigma\hackscore{00}-\sigma\hackscore{11})^2+ 100 (\sigma\hackscore{11}-\sigma\hackscore{22})^2 101 (\sigma\hackscore{22}-\sigma\hackscore{00})^2) 102 + \sigma\hackscore{01}^2+\sigma\hackscore{12}^2+\sigma\hackscore{20}^2} 103 \end{equation} 104 is used to detect material damages. Here we want to calculate the von--Mises and write the stress to a file for visualization. 105 gross 568 106 gross 579 \index{scripts!\file{diffusion.py}}: 107 The following script which is available in \file{heatedbox.py} in the \ExampleDirectory solves the Lame equation 108 and writes the displacements and the von--Mises stress\index{von--Mises stress} into a file \file{deform.xml} in the \VTK file format: 109 \begin{python} 110 from esys.escript import * 111 from esys.escript.linearPDEs import LinearPDE 112 from esys.finley import Brick 113 #... set some parameters ... 114 lam=1. 115 mu=0.1 116 alpha=1.e-6 117 xc=[0.3,0.3,1.] 118 beta=8. 119 T_ref=0. 120 T_0=1. 121 #... generate domain ... 122 mydomain = Brick(l0=1.,l1=1., l2=1.,n0=10, n1=10, n2=10) 123 x=mydomain.getX() 124 #... set temperature ... 125 T=T_0*exp(-beta*length(x-xc)) 126 #... open symmetric PDE ... 127 mypde=LinearPDE(mydomain) 128 mypde.setSymmetryOn() 129 #... set coefficients ... 130 C=Tensor4(0.,Function(mydomain)) 131 for i in range(mydomain.getDim()): 132 for j in range(mydomain.getDim()): 133 C[i,i,j,j]+=lam 134 C[j,i,j,i]+=mu 135 C[j,i,i,j]+=mu 136 msk=whereZero(x[0])*[1.,0.,0.] \ 137 +whereZero(x[1])*[0.,1.,0.] \ 138 +whereZero(x[2])*[0.,0.,1.] 139 sigma0=(lam+2./3.*mu)*alpha*(T-T_ref)*kronecker(mydomain) 140 mypde.setValue(A=C,X=sigma0,q=msk) 141 #... solve pde ... 142 u=mypde.getSolution() 143 #... calculate von-Misses stress 144 g=grad(u) 145 sigma=mu*(g+transpose(g))+lam*trace(g)*kronecker(mydomain)-sigma0 146 sigma_mises=sqrt(((sigma[0,0]-sigma[1,1])**2+(sigma[1,1]-sigma[2,2])**2+ \ 147 (sigma[2,2]-sigma[0,0])**2)/6. \ 148 +sigma[0,1]**2 + sigma[1,2]**2 + sigma[2,0]**2) 149 #... output ... 150 saveVTK("deform.xml",disp=u,stress=sigma_mises) 151 \end{python} 152 153 \begin{figure} 154 gross 599 \centerline{\includegraphics[width=\figwidth]{figures/HeatedBlock.eps}} 155 gross 579 \caption{von--Mises Stress and Displacement Vectors.} 156 \label{HEATEDBLOCK FIG 2} 157 \end{figure} 158 159 Finally the the results can be visualize by calling 160 \begin{python} 161 mayavi -d deform.xml -f CellToPointData -m VelocityVector -m SurfaceMap & 162 \end{python} 163 Note that the filter \text{CellToPointData} is applied to create smooth representation of the 164 von--Mises stress. \fig{HEATEDBLOCK FIG 2} shows the results where the vertical planes showing the 165 von--Mises stress and the horizontal plane shows the vector representing displacements. 166

## Properties

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