```

PHOTON USE
p

gr ou x 1
MSG Velocity vectors
vec x 1 sh
msg
pause
cont p1 x 1 fil;.01
msg
msg Type e to End
ENDUSE

GROUP 1. Run title
TEXT(Flow over a symmetric aerofoil:     B528
TITLE
DISPLAY
This case represents one half of an isolated symmetrical aerofoil.
The grid is curved through a right angle about the leading edge
centre so that the front portion of the North boundary is a
circular quadrant. Coding in GXBFC is used to set the inflow
boundary conditions by calculating the resolutes of the
free-stream velocity on the grid directions.
ENDDIS
REAL(PI2,PWRY,PWRZ,YV)
PI2=1.5708
INTEGER(K0,KKA,NPR)
REAL(YYC,ZZC,ZTF,DZLF);INTEGER(NI,NJ,NK)
ARRAY(YP,REAL,25);ARRAY(ZP,REAL,6)
YP(1) =0.01646;YP(2) =0.01985;YP(3) =0.02517;YP(4) =0.03485
YP(5) =0.04871;YP(6) =0.05915;YP(7) =0.06769;YP(8) =0.08108
YP(9) =0.09095;YP(10)=0.09807;YP(11)=0.10269;YP(12)=0.10481
YP(13)=0.10431;YP(14)=0.10030;YP(15)=0.09404;YP(16)=0.08607
YP(17)=0.07678;YP(18)=0.06649;YP(19)=0.05549;YP(20)=0.04416
YP(21)=0.03287;YP(22)=0.02213;YP(23)=0.01245;YP(24)=0.00449
YP(25)=0.00000
ZP(1)=0.0050;ZP(2)=0.0075;ZP(3)=0.0125;ZP(4)=0.0250
ZP(5)=0.0500;ZP(6)=0.0750
GROUP 4. Y-direction grid specification
NY=10;YVLAST=1.0;YV=1.0;PWRY=2.0
GROUP 5. Z-direction grid specification
NZ=35;PWRZ=1.5;ZWLAST=3.0;NPR=25
GROUP 6. Body-fitted coordinates or grid distortion
BFC=T;NONORT=T
NI=NX+1;NJ=NY+1;NK=NZ+1
** Number of points on curved part of outer boundary
K0=15
SETPT(1,1,1,0.0,0.0,YV);SETPT(2,1,1,1.0,0.0,YV)
** Set points along horizontal low boundary
DO II=1,2
+ DO JJ=1,NJ
+   YC(II,JJ,1)=0.0
+ ENDDO
ENDDO
DO JJ=1,NJ
+ ZZC=YV*(1.0-((JJ-1)/NY)**PWRY)
+ DO II=1,2
+   ZC(II,JJ,1)=ZZC
+ ENDDO
ENDDO
** Set points along vertical high boundary
DO JJ=1,NJ
+ YYC=YV*((JJ-1)/NY)**PWRY
+ DO II=1,2
+   DO KK=NPR+1,NZ+1
+     YC(II,JJ,KK)=YYC
+   ENDDO
+ ENDDO
ENDDO
** Set points along curved boundary
DO KK=1,K0
+ YYC=PI2*(KK-1)/(K0-1);ZZC=YV*(1.0-COS(YYC));YYC=YV*SIN(YYC)
+ DO II=1,2
+   YC(II,NJ,KK)=YYC;ZC(II,NJ,KK)=ZZC
+ ENDDO
ENDDO
** Set points along horizontal north boundary
ZTF=ZC(1,NJ,K0);DZLF=(ZC(1,NJ,NPR)-ZTF)/(NPR-K0)
DO KK=K0,NPR
+ ZZC=ZTF+DZLF*(KK-K0)
+ DO II=1,2
+   ZC(II,NJ,KK)=ZZC
+ ENDDO
ENDDO
ZTF=ZC(1,NJ,NPR);DZLF=(ZC(1,NJ,NK)-ZTF)
DO KK=NPR,NK
+ ZZC=ZTF+DZLF*((KK-NPR)/(NK-NPR))**PWRZ
+ DO II=1,2
+   ZC(II,NJ,KK)=ZZC
+ ENDDO
ENDDO
** Set points on aerofoil surface
DO KK=2,7
+ KKA=KK-1
+ DO II=1,2
+   YC(II,1,KK)=YP(KKA);ZC(II,1,KK)=ZP(KKA)+YV
+ ENDDO
ENDDO
DO KK=8,26
+ KKA=KK-1
+ DO II=1,2
+   YC(II,1,KK)=YP(KKA);ZC(II,1,KK)=(KK-6)*0.05+YV
+ ENDDO
ENDDO
** Set power law grid on south boundary aft of aerofoil
KKA=KKA+1
ZTF=ZC(1,1,KKA);DZLF=(ZC(1,1,NK)-ZTF)
DO KK=KKA,NK
+ ZZC=ZTF+DZLF*((KK-KKA)/(NK-KKA))**PWRZ
+ DO II=1,2
+   ZC(II,1,KK)=ZZC
+ ENDDO
ENDDO
** Linear transfinite interpolation
DOMAIN(1,2,1,NY+1,1,NZ+1);MAGIC(T)
MSWP=5;MAGIC(L)
DOMAIN(1,2,1,NY+1,K0,NZ+1);MSWP=3;SLIDN=T;SLIDH=T;MAGIC(L)
GROUP 7. Variables stored, solved & named
** Solve for pressure (whole-field) and velocity.
SOLVE(P1,V1,W1);SOLUTN(P1,Y,Y,Y,N,N,N)
GROUP 9. Properties of the medium (or media)
RHO1=1.0;ENUL=1.0E-10
GROUP 11. Initialization of variable or porosity fields
PATCH(IBFC,INIVAL,1,NX,2,NY,1,NZ,1,1)
INIT(IBFC,V1,0.0,GRND1);INIT(IBFC,W1,0.0,GRND1)
INIT(IBFC,VCRT,0.0,0.0);INIT(IBFC,WCRT,0.0,10.0)
GROUP 13. Boundary conditions and special sources
** Inlet
PATCH(BFCINLET,NORTH,1,NX,NY,NY,1,NZ,1,1)
COVAL(BFCINLET,P1,FIXFLU,GRND1);COVAL(BFCINLET,W1,ONLYMS,GRND1)
COVAL(BFCINLET,V1,ONLYMS,GRND1)
COVAL(BFCINLET,VCRT,0.0,0.0);COVAL(BFCINLET,WCRT,0.0,10.0)
** Pass incoming density to GXBFC through BFCA
BFCA=RHO1
** Outlet
PATCH(OUTLET,HIGH,1,NX,1,NY,NZ,NZ,1,1);COVAL(OUTLET,P1,1000.,0.0)
COVAL(OUTLET,V1,ONLYMS,0.0);COVAL(OUTLET,W1,ONLYMS,0.0)
GROUP 15. Termination of sweeps
LSWEEP=100
GROUP 17. Under-relaxation devices
RELAX(P1,LINRLX,0.6)
RELAX(V1,FALSDT,0.01); RELAX(W1,FALSDT,0.01)
GROUP 21. Print-out of variables
OUTPUT(V1,N,N,N,N,N,N)
GROUP 22. Spot-value print-out
ITABL=3;NPLT=2;IYMON=3;IZMON=10;NZPRIN=2;TSTSWP=-1
SELREF=T; RESFAC=0.01
PATCH(SURFACE,PROFIL,1,1,1,1,1,NZ,1,1)
PLOT(SURFACE,P1,0.0,0.0);PLOT(SURFACE,W1,0.0,0.0)
```