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

revision 578 by gross, Mon Mar 6 06:12:04 2006 UTC revision 579 by gross, Tue Mar 7 01:28:23 2006 UTC
# Line 4  Line 4
4  In this section we want to discuss the deformation of a linear elastic body caused by expansion through a heat distribution. We want  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  to displacement field $u\hackscore{i}$ which solves the momentum equation
6  \index{momentum equation}:  \index{momentum equation}:
7  \begin{eqnarray}\label{BEAM general problem}  \begin{eqnarray}\label{HEATEDBLOCK general problem}
8   - \sigma\hackscore{ij,j}=0   - \sigma\hackscore{ij,j}=0
9  \end{eqnarray}  \end{eqnarray}
10  where the stress $\sigma$ is given by  where the stress $\sigma$ is given by
11  \begin{eqnarray}\label{BEAM linear elastic}  \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})   \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} \;.   - (\lambda+\frac{2}{3} \mu)  \; \alpha  \;  (T-T\hackscore{ref})\delta\hackscore{ij} \;.
14  \end{eqnarray}  \end{eqnarray}
15  In this formula $\lambda$ and $\mu$ are the Lame coefficients, $\alpha$ is the  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  temperature expansion coefficient, $T$ is the temperature distribution and $T_{ref}$ a reference temperature. Note that
17  \eqn{BEAM general problem} is similar to eqn{WAVE general problem} introduced in section~\Sec{WAVE CHAP} but the  \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  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}  comparison to the \eqn{WAVE stress}
20  definition of stress $\sigma$ in \eqn{BEAM linear elastic} an extra term is introduced  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.    to bring in stress due to volume changes trough temperature dependent expansion.
22
23  Our domain is the unit cube  Our domain is the unit cube
24  \begin{eqnarray} \label{BEAM natural}  \begin{eqnarray} \label{HEATEDBLOCK natural}
25  \Omega=\{(x\hackscore{i} | 0 \le x\hackscore{i} \le 1 \}  \Omega=\{(x\hackscore{i} | 0 \le x\hackscore{i} \le 1 \}
26  \end{eqnarray}  \end{eqnarray}
27  On the boundary the normal stress component is set to zero  On the boundary the normal stress component is set to zero
28  \begin{eqnarray} \label{BEAM natural}  \begin{eqnarray} \label{HEATEDBLOCK natural}
29  \sigma\hackscore{ij}n\hackscore{j}=0  \sigma\hackscore{ij}n\hackscore{j}=0
30  \end{eqnarray}  \end{eqnarray}
31  and on the face with $x\hackscore{i}=0$ we set the $i$-th component of the displacement to $0$  and on the face with $x\hackscore{i}=0$ we set the $i$-th component of the displacement to $0$
32  \begin{eqnarray} \label{BEAM constraint}  \begin{eqnarray} \label{HEATEDBLOCK constraint}
33  u\hackscore{i}(x)=0 & \mbox{ where } & x\hackscore{i}=0 \;  u\hackscore{i}(x)=0 & \mbox{ where } & x\hackscore{i}=0 \;
34  \end{eqnarray}  \end{eqnarray}
35  For the temperature distribution we use  For the temperature distribution we use
36  \begin{eqnarray} \label{BEAM temperature}  \begin{eqnarray} \label{HEATEDBLOCK temperature}
37  T(x)= \frac{\beta}{\|x-x^{c}\|};  T(x)= T\hackscore{0} e^{-\beta \|x-x^{c}\|};
38  \end{eqnarray}  \end{eqnarray}
39  with a given positive constant $\beta$ and location $x^{c}$ in the domain\footnote{This selection of $T$ corresponds to  with a given positive constant $\beta$ and location $x^{c}$ in the domain. Later in \Sec{MODELFRAME} we will use
40  a temperature distribution in an indefinite domain created by a nodal heat source at $x^{c}$. Later in \Sec{X} we will calculate  $T$ from an time-dependent temperature diffusion problem as discussed in \Sec{DIFFUSION CHAP}.
the $T$ from an time-dependent temperature diffusion problem as discussed in \Sec{DIFFUSION CHAP}.}
41
42  When we insert~\eqn{BEAM linear elastic} we get a second oder system of linear PDEs for the displacements $u$ which is called  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  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  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  takes PDEs of the form
# Line 66  A\hackscore{ijkl}=\lambda \delta\hacksco Line 65  A\hackscore{ijkl}=\lambda \delta\hacksco
65  X\hackscore{ij}=(\lambda+\frac{2}{3} \mu) \;  \alpha \; (T-T\hackscore{ref})\delta\hackscore{ij} \\  X\hackscore{ij}=(\lambda+\frac{2}{3} \mu) \;  \alpha \; (T-T\hackscore{ref})\delta\hackscore{ij} \\
66  \end{eqnarray}  \end{eqnarray}
67  The characteristic function $q$ defining the locations and components where constraints are set is given by:  The characteristic function $q$ defining the locations and components where constraints are set is given by:
69  q\hackscore{i}(x)=\left\{  q\hackscore{i}(x)=\left\{
70  \begin{array}{cl}  \begin{array}{cl}
71  1  & x\hackscore{i}=0  \\  1  & x\hackscore{i}=0  \\
# Line 77  q\hackscore{i}(x)=\left\{ Line 76  q\hackscore{i}(x)=\left\{
76  Under the assumption that $\lambda$, $\mu$, $\beta$ and $T\hackscore{ref}$  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,  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  this choice would lead to a different natural boundary condition which does not set the normal stress component as defined
79  in~\eqn{BEAM linear elastic} to zero.  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  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}  \index{symmetric PDE}
# Line 88  D\hackscore{ik}=D\hackscore{ki} \\ Line 87  D\hackscore{ik}=D\hackscore{ki} \\
87  d\hackscore{ik}=d\hackscore{ki} \  d\hackscore{ik}=d\hackscore{ki} \
88  \end{eqnarray}  \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  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 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
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
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
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

Legend:
 Removed from v.578 changed lines Added in v.579