TALK=T;RUN(1,1) 749.htm

  DISPLAY
  
  This file illustrates how to create simple shapes
  by using the SPHERE In-Form object.

  The shape and grid can be selected by setting CASENO 
  from 1 to 11 according to the following table
  
  grid              <----------- shape ---------------->
  ----               sphere   cone   spiral  curved duct !poros ?
  cartesian             1      4       7       10          yes
  cylindrical-polar     2      5       8                   no
  body-fitted           3      6       9       11          no
  -----------------------------------------------------

  The temperature equation is solved with, only conduction active
  The Q1 contains PHOTON USE commands
  ENDDIS  
  PHOTON USE
  p;;;
 
  up z
  vi 3 2 1
  gr ou x m col 6; gr ou y m col 6; gr ou z m col 6
  gr ou x 1 col 6; gr ou y 1 col 6; gr ou z 1 col 6
  set prop off
  surf mark x .99 col 2
  surf mark y .99 col 4
  surf mark z .99 col 6
  dump 249
  *end
  ENDUSE
 
 ************************************************************
  Group 1. Run Title
 LIBREF=749
 TITLE
boolean(alterbfc,altercrt)
alterbfc=f   ! set one of these =t so as the make
altercrt=f   ! the grid non-uniform

INTEGER(CASENO)
mesg(CASENO=1 .. Simple SPHERE for cartesian grid
mesg(CASENO=2 .. Simple SPHERE for BFC
mesg(CASENO=3 .. Simple SPHERE for polar grid
mesg(CASENO=4 .. Simple cone for cartesian grid
mesg(CASENO=5 .. Simple cone for BFC
mesg(CASENO=6 .. Simple cone for polar grid
mesg(CASENO=7 .. Simple spiral for cartesian grid
mesg(CASENO=8 .. Simple spiral for BFC
mesg(CASENO=9 .. Simple spiral for polar grid
mesg(CASENO=10 .. Simple curve duct for cartesian grid
mesg(CASENO=11 .. Simple curve duct for BFC
mesga(Please select case number
readvdu(caseno,int,1)

REAL(PI,PI2); PI=3.14159; PI2=2.*PI  ! Pythagoras's constant
  ************************************************************
  Group 2. Transience
 STEADY=T
 ************************************************************
  Groups 3, 4, 5  Grid Information
 NX=50;NY=50;NZ=50
 XULAST=1.0; YVLAST=1.0; ZWLAST=1.0
IF(CASENO.EQ.3.OR.CASENO.EQ.6.OR.CASENO.EQ.9) THEN
 CARTES=F
 XULAST=PI2
 altercrt=f
 alterbfc=f
ENDIF
TEXT(Shapes made bu sphere; caseno=:caseno:
 
grdpwr(x,nx   ,xulast,1.0)
grdpwr(y,ny   ,yvlast,1.0)
grdpwr(z,nz   ,zwlast,1.0)
if(altercrt) then
grdpwr(x,nx   ,xulast,0.5)
grdpwr(y,ny   ,yvlast,-0.5)
grdpwr(z,-nz   ,zwlast,0.5)
endif
IF(CASENO.EQ.2.OR.CASENO.EQ.5.OR.CASENO.EQ.8.OR.CASENO.EQ.11) THEN
 BFC=T     ! 
 if(alterbfc) then
 nx=20;ny=20;nz=20
 grdpwr(x,nx   ,xulast,1.0)
 grdpwr(y,ny   ,yvlast,1.0)
 grdpwr(z,nz   ,zwlast,1.0)
 real(xchanged,dummy1,dummy2)
 BFC=T
 dummy2=nx
 dummy2=1./dummy2
 dummy2
 do ixx=2,nx+1
+  dummy1=ixx-1
+  dummy1=(dummy1*dummy2)**.5  ! distort x-direction grid
+  do iyy=1,ny+1
+    do izz=1,nz+1
+      xc(:ixx:,:iyy:,:izz:)=dummy1
+    enddo
+  enddo
 enddo
 endif
ENDIF
 
 
if(cartes) then
 mesg(uniform cartesian grid
else
 mesg(uniform polar grid
endif
  Group 7. Variables: STOREd,SOLVEd,NAMEd
    * Solved variables list
 SOLVE(TEM1)
 Liter(tem1)=1000
    * Stored variables list
 STORE(PRPS,MARK)
IF(CASENO.EQ.1.OR.CASENO.EQ.4.OR.CASENO.EQ.7.OR.CASENO.EQ.10) THEN
 STORE(VPOR); FIINIT(VPOR)=1.
ENDIF
 ************************************************************
  Group 11.Initialise Var/Porosity Fields
 patch(.patch1,inival,0,1000,0,1000,0,1000,0,1000)

FIINIT(PRPS) =  0.000000E+00 ! 
  inform11begin
(stored of mark at .patch1 is 1.0 with infob_1)

(initial of prps is 100 with infob_1)
char(x0,y0,z0,rad)
IF(CASENO.EQ.1) THEN ! Cartesian; fixed sphere
 x0=xulast/2
 y0=yvlast/2
 z0=zwlast/2
 rad=xulast/4

 (infob at .patch1 is sphere(x0,y0,z0,rad) with infob_1!poros)
endif
IF(CASENO.EQ.2) THEN ! bfc; fixed sphere
xulast
 x0=xulast/2
 y0=yvlast/2
 z0=zwlast/2
 rad=xulast/4
 (infob at .patch1 is sphere(x0,y0,z0,rad) with infob_1)
ENDIF
IF(CASENO.EQ.3) THEN ! polar; fixed sphere
 x0=yvlast
 y0=yvlast
 z0=zwlast/2
 rad=yvlast/2
 (infob at .patch1 is sphere(x0,y0,z0,rad) with infob_1)
ENDIF
IF(CASENO.EQ.4) THEN ! cartesian; cone
 x0=xulast/2
 y0=yvlast/2
 z0=zg
 rad=0.25*(:zwlast:-zg)  ! radius diminishes linearly with z
 (infob at .patch1 is sphere(x0,y0,z0,rad) with infob_1!poros)
ENDIF
IF(CASENO.EQ.5) THEN ! bfc; cone
 x0=xulast/2
 y0=yvlast/2
 z0=zg
 rad=0.25*(:zwlast:-zg)  ! radius diminishes linearly with z
 (infob at .patch1 is sphere(x0,y0,z0,rad) with infob_1)
ENDIF
IF(CASENO.EQ.6) THEN ! polar; cone
 x0=yvlast
 y0=yvlast
 z0=zg
 rad=0.5*(:zwlast:-zg)   ! radius diminishes linearly with z
 (infob at .patch1 is sphere(x0,y0,z0,rad) with infob_1)
ENDIF
IF(CASENO.EQ.7) THEN ! Cartesian; spiral
 x0=:xulast/2:*(1+0.8*cos(:PI2:*zg)) ! centre cordinates vary
 y0=:yvlast/2:*(1+0.8*sin(:PI2:*zg)) ! with z
 z0=zg                    !
 rad=zwlast/10                  ! radius is constant
 (infob at .patch1 is sphere(x0,y0,z0,rad) with infob_1!poros)
ENDIF
IF(CASENO.EQ.8) THEN ! bfc; spiral
 x0=:xulast/2:*(1+0.8*cos(:PI2:*zg)) ! centre cordinates vary
 y0=:yvlast/2:*(1+0.8*sin(:PI2:*zg)) ! with z                
 z0=zg                    !                       
 rad=zwlast/10                  ! radius is constant    
 (infob at .patch1 is sphere(x0,y0,z0,rad) with infob_1)
ENDIF
IF(CASENO.EQ.9) THEN ! polar; spiral
 rad=yvlast/4                       ! radius is constant    
 x0=:yvlast:*(1+2*:rad:*cos(:PI2:*zg)) ! centre cordinates vary
 y0=:yvlast:*(1+2*:rad:*sin(:PI2:*zg)) ! with z                
 z0=zg                        !                       
 (infob at .patch1 is sphere(x0,y0,z0,rad) with infob_1)
ENDIF
IF(CASENO.EQ.10) THEN ! Cartesian; curved duct
 x0=xg                        ! centre cordinates vary
 y0=:yvlast/2:                           ! with z                
 z0=:zwlast/2:+0.25*sin(:pi2:*xg/xulast) !                       
 rad=yvlast/4                       ! radius is constant    
 (infob at .patch1 is sphere(x0,y0,z0,rad) with infob_1!poros)
ENDIF
IF(CASENO.EQ.11) THEN ! BFC; curved duct
 x0=xg                        ! centre cordinates vary
 y0=:yvlast/2:                           ! with z                
 z0=:zwlast/2:+0.25*sin(:pi2:*xg/xulast) !                       
 rad=yvlast/4                       ! radius is constant    
 (infob at .patch1 is sphere(x0,y0,z0,rad) with infob_1)
ENDIF
  inform11end
 ************************************************************
  Group 13. Boundary & Special Sources
patch(low,lwall,1,nx,1,ny,1,1,1,1)
coval(low,tem1,1,-1)
patch(high,hwall,1,nx,1,ny,nz,nz,1,1)
coval(high,tem1,1,1)
 ************************************************************
  Group 15. Terminate Sweeps
 LSWEEP  =      3
 
yzpr=t
nxprin=1;ixprf=nx/2;ixprl=nx/2+1
spedat(set,material,100,l,t)
tstswp=-1
isg50=1       ! endpause
isg52=2       ! maxabs
store(xcen,ycen,zcen)
stop
 ************************************************************