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

Revision 586 - (show annotations)
Fri Mar 10 01:24:27 2006 UTC (17 years ago) by lkettle
File MIME type: application/x-tex
File size: 7936 byte(s)
corrected error in natural b.c.


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

## Properties

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