Encyclopaedia Index
``` TALK=F;RUN( 1, 1)

TEXT( Library case Y622: Thermal stresses in composite solid block

>>>>>>>>>>>>>>>>>>>>>> Comment begins >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
PLANT information :
* Data input groups used: 11, 13, 19
* Ground groups planted : 11, 13, 19-6
* Headings used  : INIT??, SORC??, SCO6??
* Functions used : None
* Commands used  : None
<<<<<<<<<<<<<<<<<<<<<<<<<  Comment ends <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<

REAL(LENGTH,WIDTH,ALFA,TISO,STRNY,STRNX,STIFFN,STRSX,STRSY)
REAL(NORTEM,SOUTEM,HEATSOR,POISSN)
BOOLEAN(ISOTHER,TEMLIN,TEMSOR)

TISO=100.0
LENGTH=1.0;WIDTH =1.0

STIFFN=2.e11
STRSY =-1.e8;STRSX=0.0
ALFA  =1.e-05          ; RG(1)=ALFA
STRNY =STRSY           ; RG(2)=STRNY
STRNX =STRSX           ; RG(3)=STRNX
POISSN=0.3333          ; RG(4)=POISSN

NX  =10 ; NY  =10

ISOTHER=t;TEMLIN=f;TEMSOR=f

GROUP 3. X-direction grid specification
GRDPWR(X,NX,LENGTH,1.0)
GROUP 4. Y-direction grid specification
GRDPWR(Y,NY,WIDTH,1.0)
GROUP 7. Variables stored, solved & named
SOLVE(P1,U1,V1,TEM)
Provide the storage for thermal expansion coefficient
,THEX, stiffnes, STIF, and Poisson ratio to be initialised
for the presence of different materials

STORE(THEX,STIF,POIS)

STORE(CON1,CON2,CON3,CON4)
GROUP 8. Terms (in differential equations) & devices
TERMS(U1,Y,N,Y,Y,Y,Y)
TERMS(V1,Y,N,Y,Y,Y,Y)
terms(tem,n,n,y,p,p,p)
GROUP 9. Properties of the medium (or media)
ENUL=1.
GROUP 11. Initialization of variable or porosity fields
IF(ISOTHER) THEN
PATCH(STELP,INIVAL,1,NX,1,NY,1,1,1,1)
INIT (STELP,POIS,0.0,POISSN)

PATCH(STELT,INIVAL,1,NX,1,NY,1,1,1,1)
INIT (STELT,THEX,0.0,ALFA)

PATCH(STELS,INIVAL,1,NX,1,NY,1,1,1,1)
INIT (STELS,STIF,0.0,STIFFN)

PATCH(STELC,INIVAL,1,NX,1,NY,1,1,1,1)
VAL=(1.-RG(4))/((1.+RG(4))*(1-2.*RG(4)))
INIT (STELC,CON1,0.0,GRND)
VAL=RG(4)/((1.+RG(4))*(1-2.*RG(4)))
INIT (STELC,CON2,0.0,GRND)
VAL=1./(1.-2.*RG(4))
INIT (STELC,CON3,0.0,GRND)
VAL=1./(2.*(1+RG(4)))
INIT (STELC,CON4,0.0,GRND)

fiinit(tem)=TISO
#
PATCH(ALMT1,INIVAL,1,3,NY/2+1,NY,1,1,1,1)
INIT (ALMT1,THEX,0.0,10.*ALFA)
INIT (ALMT1,STIF,0.0,0.1*STIFFN)
PATCH(ALMT2,INIVAL,8,NX,NY/2+1,NY,1,1,1,1)
INIT (ALMT2,THEX,0.0,10.*ALFA)
INIT (ALMT2,STIF,0.0,0.1*STIFFN)
ENDIF
GROUP 13. Boundary conditions and special sources
** Linear temperature field
IF(TEMLIN) THEN
PATCH(TEMPN,NWALL,1,NX,NY,NY,1,1,1,1)
COVAL(TEMPN,TEM,1.0,0.0)
PATCH(TEMPS,SWALL,1,NX,1,1,1,1,1,1)
COVAL(TEMPS,TEM,1.0,1.0)
ENDIF
** Source generated temperature field
IF(TEMSOR) THEN
PATCH(TEMPN,NWALL,1,NX,NY,NY,1,1,1,1)
COVAL(TEMPN,TEM,1.0,0.0)
PATCH(TEMSOR,VOLUME,1,NX,1,NY,1,1,1,1)
COVAL(TEMSOR,TEM,FIXFLU,2.0)
ENDIF
Stress - in -Solid BC
---------------- WEST ----------------------------
** Zero displacement at the west boundary
PATCH(WESZD,WWALL,1,1,1,NY,1,1,1,1)
COVAL(WESZD,U1,1.,0.0)
---------------- EAST ----------------------------
** Normal stress at the east solid boundary
PATCH(EASZN,EAST,NX-1,NX-1,1,NY,1,1,1,1)
VAL=U1[-1,,]+0.333*(P1+4.*(THEX*TEM+RG(3)/STIF))*DXU2D
COVAL(EASZN,U1,FIXVAL,GRND)
** Zero shear stress at the east boundary
PATCH(EASZS,EAST,NX,NX,1,NY-1,1,1,1,1)
VAL=-(U1[-1,+1,]-U1[-1,,])/DYG2D
COVAL(EASZS,V1,FIXFLU,GRND)
** Dilatation at the east solid boundary
PATCH(FXDEAS,EAST,NX,NX,2,NY-1,1,1,1,1)
VAL=2.*(THEX*TEM-RG(3)/STIF)-1.5*(V1-V1[,-1,])/DYV2D
COVAL(FXDEAS,P1,FIXVAL, grnd)

PATCH(FXDEAC,EAST,NX,NX,1,1,1,1,1,1)
VAL=2.*(THEX*TEM-RG(3)/STIF)-\$
1.5*(V1[,+1,]-V1)/DYV2D[,+1,]
COVAL(FXDEAC,P1,FIXVAL, grnd)

---------------- NORTH ----------------------------
** Normal stress at the free north boundary
PATCH(NORZN,CELL,1,NX,NY-1,NY-1,1,1,1,1)
VAL=V1[,-1,]+0.333*(P1+4.*(THEX*TEM+RG(2)/STIF))*DYV2D
COVAL(NORZN,V1,FIXVAL,grnd)
** Zero shear stress at the free north boundary
PATCH(NORZS,NORTH,1,NX-1,NY,NY,1,1,1,1)
VAL=-(V1[+1,-1,]-V1[,-1,])/DXG2D
COVAL(NORZS,U1,FIXFLU,GRND)
** Dilatation at the free north boundary
PATCH(FDNW,NORTH,1,1,NY,NY,1,1,1,1)
VAL=2.*(THEX*TEM-RG(2)/STIF)-\$
1.5*(U1[+1,,]-U1)/DXU2D[+1,,]
COVAL(FDNW,P1,FIXVAL, grnd)
PATCH(FDNOR,NORTH,2,NX-1,NY,NY,1,1,1,1)
VAL=2.*(THEX*TEM-RG(2)/STIF)-1.5*(U1-U1[-1,,])/DXU2D
COVAL(FDNOR,P1,FIXVAL, grnd)
PATCH(FDNEC,NORTH,NX,NX,NY,NY,1,1,1,1)
VAL=2.*(THEX*TEM-RG(2)/STIF)-\$
1.5*(U1[-1,,]-U1[-2,,])/DXU2D[-1,,]
COVAL(FDNEC,P1,FIXVAL, grnd)
---------------- SOUTH ----------------------------
** Zero displacement at the south symmetry
PATCH(SOUZD,SWALL,1,NX,1,1,1,1,1,1)
COVAL(SOUZD,V1,1.,0.0)
==================================================
** Volumetric mass source for dilatation
PATCH(DILAT,VOLUME,1,nx,1,ny,1,1,1,1)
VAL=-P1/CON3+THEX*TEM/CON4
COVAL(DILAT,P1,FIXFLU,GRND)
===================================================
OUTPUT INFORMATION
** Dilatation
STORE(DIL)
DIL=-P1/CON3+THEX*TEM/CON4
** Stress and strain calculations
STORE(STRX,STRY,TAXY)
STORE(EPSX,EPSY,EPST)
* Thermal strains
EPST=THEX*TEM
REGION(1,NX,1,NY,1,1)
* Normal strains
EPSX=(U1-WEST(U1))/DXU2D
REGION(2,NX-1,1,NY,1,1)
EPSX=(U1[+1,,]-U1)/DXU2D[+1,,]
REGION(1,1,1,NY,1,1)
EPSX=(U1[-1,,]-U1[-2,,])/DXU2D[-1,,]
REGION(NX,NX,1,NY,1,1)
EPSY=(V1-SOUTH(V1))/DYV2D
REGION(1,NX,2,NY-1,1,1)
EPSY=(V1[,+1,]-V1)/DYV2D[,+1,]
REGION(1,NX,1,1,1,1)
EPSY=(V1[,-1,]-V1[,-2,])/DYV2D[,-1,]
REGION(1,NX,NY,NY,1,1)
* Normal stresses
STRX=STIF*(CON1*EPSX+CON2*EPSY-CON3*EPST)
REGION(1,NX,1,NY,1,1)
STRY=STIF*(CON1*EPSY+CON2*EPSX-CON3*EPST)
REGION(1,NX,1,NY,1,1)
* Shear stresses
STORE(DDYU,DDXV)
== ddy.u ==
DDYU=0.25*((U1-U1[,-1,])/DYG2D[,-1,]+ \$
(U1[,+1,]-U1)/DYG2D+(U1[-1,,]-\$
U1[-1,-1,])/DYG2D[,-1,]+\$
(U1[-1,+1,]-U1[-1,,])/DYG2D)
REGION(2,NX-1,2,NY-1,1,1)

DDYU=0.5*((U1-U1[,-1,])/DYG2D[,-1,]+\$
(U1[,+1,]-U1)/DYG2D)
REGION(1,1,2,NY-1,1,1)

DDYU=(U1[,+1,]-U1)/DYG2D
REGION(1,1,1,1,1,1)

DDYU=(U1-U1[,-1,])/DYG2D[,-1,]
REGION(1,1,NY,NY,1,1)

DDYU=0.5*((U1[-1,,]-U1[-1,-1,])/DYG2D[,-1,]+\$
(U1[-1,+1,]-U1[-1,,])/DYG2D)
REGION(NX,NX,2,NY-1,1,1)

DDYU=(U1[-1,+1,]-U1[-1,,])/DYG2D
REGION(NX,NX,1,1,1,1)

DDYU=(U1[-1,,]-U1[-1,-1,])/DYG2D[,-1,]
REGION(NX,NX,NY,NY,1,1)

DDYU=0.5*((U1-U1[,-1,])/DYG2D[,-1,]+\$
(U1[-1,,]-U1[-1,-1,])/DYG2D[,-1,])
REGION(2,NX-1,NY,NY,1,1)

DDYU=0.5*((U1[,+1,]-U1)/DYG2D+\$
(U1[-1,+1,]-U1[-1,,])/DYG2D)
REGION(2,NX-1,1,1,1,1)

== ddx.v ==
DDXV=0.25*((V1[+1,-1,]-V1[,-1,]  )/DXG2D+\$
(V1[+1,,]-V1)/DXG2D+(V1[,-1,]-\$
V1[-1,-1,])/DXG2D[-1,,]+\$
(V1-V1[-1,,])/DXG2D[-1,,])
REGION(2,NX-1,2,NY-1,1,1)

DDXV=(V1[+1,,] -V1)/DXG2D
REGION(1,1,1,1,1,1)

DDXV=0.5*((V1[+1,-1,]-V1[,-1,])/DXG2D+\$
(V1[+1,,]-V1)/DXG2D)
REGION(1,1,2,NY-1,1,1)

DDXV=(V1[+1,-1,]-V1[,-1,]  )/DXG2D
REGION(1,1,NY,NY,1,1)

DDXV=(V1-V1[-1,,])/DXG2D[-1,,]
REGION(NX,NX,1,1,1,1)

DDXV=0.5*((V1[,-1,]-V1[-1,-1,])/DXG2D[-1,,]+\$
(V1-V1[-1,,])/DXG2D[-1,,])
REGION(NX,NX,2,NY-1,1,1)

DDXV=(V1[,-1,] -V1[-1,-1,])/DXG2D[-1,,]
REGION(NX,NX,NY,NY,1,1)

DDXV=0.5*((V1[+1,,]-V1)/DXG2D+\$
(V1-V1[-1,,])/DXG2D[-1,,])
REGION(2,NX-1,1,1,1,1)

DDXV=0.5*((V1[+1,-1,]-V1[,-1,])/DXG2D+\$
(V1[,-1,]-V1[-1,-1,])/DXG2D[-1,,])
REGION(2,NX-1,NY,NY,1,1)

TAXY=STIF*CON4*(DDYU+DDXV)
REGION(1,NX,1,NY,1,1)

RELAX(P1,LINRLX,0.005)
RELAX(U1,FALSDT,0.0005)
RELAX(V1,FALSDT,0.0005)

NAMSAT=MOSG

GROUP 15. Termination of sweeps
LSWEEP=3000
selref=f
RESREF(P1)=1.E-8
RESREF(U1)=1.E-8;RESREF(V1)=1.E-8
tstswp=-1
GROUP 22. Spot-value print-out
IYMON=NY-1
GROUP 23. Field print-out and plot control
NXPRIN=1;NYPRIN=1
LIBREF=458
PHOTON USE
p

gr ou z 1
con thex z 1 fil;.001
gr ou z 1 x 1 3 y 6 m
gr ou z 1 x 8 m y 6 m
msg Linear thermal expansion coefficient
pause
con cl;red
con stif z 1 fil;.001
gr ou z 1 x 1 3 y 6 m
gr ou z 1 x 8 m y 6 m
msg Stiffnes ( Young's module ) contours
pause
con cl;red
vec z 1 sh
msg Velocity and displacement vectors
pause
vec cl;red
con u1 z 1 fil;.001
gr ou z 1 x 1 3 y 6 m
gr ou z 1 x 8 m y 6 m
msg X-component displacement contours
pause
con cl;red
con v1 z 1 fil;.001
gr ou z 1 x 1 3 y 6 m
gr ou z 1 x 8 m y 6 m
msg Y-component displacement contours
pause
con cl;red
con dil z 1 fil;.001
gr ou z 1 x 1 3 y 6 m
gr ou z 1 x 8 m y 6 m
msg Dilatation contours
pause
con cl;red
con epsx z 1 fil;.001
gr ou z 1 x 1 3 y 6 m
gr ou z 1 x 8 m y 6 m
msg X-strain contours
pause
con cl;red
con epsy z 1 fil;.001
gr ou z 1 x 1 3 y 6 m
gr ou z 1 x 8 m y 6 m
msg Y-strain contours
pause
con cl;red
con epst z 1 fil;.001
gr ou z 1 x 1 3 y 6 m
gr ou z 1 x 8 m y 6 m
msg Thermal strain contours
pause
con cl;red
con strx z 1 fil;.001
gr ou z 1 x 1 3 y 6 m
gr ou z 1 x 8 m y 6 m
msg X-normal stress contours
pause
con cl;red
con stry z 1 fil;.001
gr ou z 1 x 1 3 y 6 m
gr ou z 1 x 8 m y 6 m
msg Y-normal stress contours
pause
con cl;red
con taxy z 1 fil;.001
gr ou z 1 x 1 3 y 6 m
gr ou z 1 x 8 m y 6 m
msg Shear stress contours
pause
msg Press E to end
enduse
STOP
```