߂
C
PROGRAM MAIN
C
C ******************************************************************** C
C C
C Graphics Software for HF3D-FEM-R C
C C
C ( Rectangular Parallelepiped Elements ) C
C C
C ( Version 3.8 ) C
C C
C 3D-View, Velocity distribution & Contour lines C
C C
C Copyright : Yasuhiro MATSUDA C
C C
C ******************************************************************** C
C
IMPLICIT REAL*8 (A-H,O-Z)
C
C ----- DATA ( 1) ----------------------------------------------------
C
*** PARAMETER ( NEX = 3, NEY = 3, NEZ = 3,
C
*** PARAMETER ( NEX = 4, NEY = 4, NEZ = 4,
C
PARAMETER ( NEX = 10, NEY = 10, NEZ = 10,
*** PARAMETER ( NEX = 20, NEY = 20, NEZ = 20,
*** PARAMETER ( NEX = 30, NEY = 30, NEZ = 30,
C
C ----- Attention : NEX etc. < = 50 ---
C
1 LX = NEX + 1, LY = NEY + 1, LZ = NEZ + 1, NNP = LX* LY* LZ,
2 NXY= NEX* NEY, NYZ= NEY* NEZ, NZX= NEZ* NEX )
C
C --------------------------------------------------------------------
C
DIMENSION X ( NNP), Y ( NNP), Z ( NNP), U ( NNP),
1 V ( NNP), W ( NNP), T ( NNP), XLX( LX),
2 YLY( LY), LXY( NXY,4), LYZ( NYZ,4), LZX( NZX,4)
C
DIMENSION VVR( 8000), VMX( 8000)
C
COMMON /PFAC/ ZFAC
COMMON /CNST/ IKY, MSH, LOP
C
C ----- DATA ( 2) ----------------------------------------------------
C
*** OPEN ( 11, FILE='HF3D_FEM_UM3_710_OUT.DAT' ) ! DT = 0.04
C
*** OPEN ( 11, FILE='HF3D_FEM_UM4_710_OUT.DAT' ) ! DT = 0.04
C
OPEN ( 11, FILE='HF3D_FEM_UM10_710_OUT.DAT' ) ! DT = 0.014
*** OPEN ( 11, FILE='HF3D_FEM_M20_7100_001_OUT' )
*** OPEN ( 11, FILE='HF3D_FEM_M30_10X6_00002_OUT' )
*** OPEN ( 11, FILE='HF3D_FEM_M30_10X6_000025_OUT' )
C
C --------------------------------------------------------------------
C
CALL TTLE
C
C -----------------------------------------------------------------------
C
CALL INPT ( U, V, W, T, VVR, VMX, X, Y, Z,
1 LXY, LYZ, LZX, XLX, YLY, NNP, LYX, LYY, NXY,
2 NYZ, NZX, RA )
C
C -----------------------------------------------------------------------
C
CALL VGA@
C
MSH = 1
C ----------------------------------------------------------------------
C
IF ( IKY.EQ.1) CALL MSHF ( XLX, YLY, LAX, LAY, NEX, MSH )
C
C ------------------------------------------------
IF ( IKY.EQ.2) THEN
C
CALL PLVR ( VVR, 1, RA )
C
CALL PLVR ( VMX, 2, RA )
C
CALL VW3D ( X, Y, Z, U, V, W, NNP )
C
END IF
C ------------------------------------------------
C
IF ( IKY.EQ.3)
C
1 CALL PVLTM ( U, V, W, T, X, Y, Z, LXY,
2 LYZ, LZX, NNP, NXY, NYZ, NZX )
C
C ---------------------------------------------------------------------
C
IF ( IKY.EQ.4) CALL CMPT ( U, V, X, Y, NNP, RA )
C
C ---------------------------------------------------------------------
CLOSE ( 11)
C
STOP
END
C **********************************************************************
C
SUBROUTINE CMPT ( UN, VN, X, Y, NNP, RA )
C
IMPLICIT REAL*8(A-H,O-Z)
C
INTEGER*2 ICR, IX1, IY1, IX2, IY2
C
DIMENSION UN( NNP), VN( NNP), X( NNP), Y( NNP)
C
NPL = 0
C
IX1 = 0
IY1 = 0
IX2 = 639
IY2 = 479
C
ICR = 7
CALL CLEAR_SCREEN_AREA@( IX1, IY1, IX2, IY2, ICR )
C
C -----------------------------------------------
C
CALL RFCP ( X, Y, UN, VN, NNP, RA )
C
C -----------------------------------------------
RETURN
END
C **********************************************************************
C
SUBROUTINE INPT ( U, V, W, T, VVR, VMX, X, Y,
1 Z, LXY, LYZ, LZX, XLX, YLY, NNP, LYX,
2 LYY, NXY, NYZ, NZX, RA )
C
IMPLICIT REAL*8 (A-H,O-Z)
C
DIMENSION U ( NNP), V( NNP), W( NNP), T( NNP), XLX( LYX),
1 YLY( LYY ), X( NNP), Y( NNP), Z( NNP)
C
DIMENSION VVR( 8000), VMX( 8000), LXY( NXY,4), LYZ( NYZ,4),
1 LZX( NZX,4)
C
COMMON /PFAC/ ZFAC
COMMON /CNST/ IKY, MSH, LOP
C
C -----------------------------------------------
C
READ(11,1000) RA, LYX, LYY, NEX, DT
1000 FORMAT(D10.2,3I3,D11.4)
WRITE(6,*) ' '
IF ( NEX.LE.3) THEN
WRITE(6,*) ' ** IKY = 1 : ( Meshes )'
WRITE(6,*) ' = 2 : ( Var. & Vel. & Perspective-View)'
END IF
C
C ----- Attention : When defined for X=0.5, Y=0.5 and Z=0.5 ---
C
IF ( NEX.GE.4) THEN
WRITE(6,*) ' '
WRITE(6,*) ' ** IKY = 1 : ( Meshes )'
WRITE(6,*) ' = 2 : ( Var. & Vel. & Perspective-View )'
WRITE(6,*) ' = 3 : ( X-Y face) ( Y-Z face) ( Z-X face)'
WRITE(6,*) ' = 4 : U - V Distr. on Z = 0.5'
END IF
C
READ(6,1100) IKY
1100 FORMAT(I3)
WRITE(6,2000) IKY
2000 FORMAT(' ** IKY =',I3)
C
READ(11,1200) ( XLX( I), I=1, LYX), ( YLY( I), I=1, LYY)
READ(11,1300) LOP, ZFAC
1300 FORMAT(I5,D11.4)
READ(11,1200) ( VVR( I), I=1, LOP)
READ(11,1200) ( VMX( I), I=1, LOP)
C
C ----- Attention : From " 2 ` LOP" DATA ---
READ(11,1500) NNP
1500 FORMAT( I5)
C
C -----------------------------------------------
WRITE(6,2100) RA, DT, LOP, ZFAC, NNP
2100 FORMAT(/,' ** Ra =',F9.2,' DT =',F6.3,' LOP =',I4,' Zfac =',
1 F5.2,' NNP =',I5 )
READ(11,1200) ( X( I), Y( I), Z( I), I=1, NNP)
READ(11,1200) ( U( I), V( I), W( I), I=1, NNP)
1200 FORMAT(10D11.4)
READ(11,1600) NXY, NYZ, NZX
1600 FORMAT(3I3)
READ(11,1200) ( T( I), I=1, NNP)
READ(11,1700) ( ( LXY( I,J), I=1, NXY ), J = 1, 4)
READ(11,1700) ( ( LYZ( I,J), I=1, NYZ ), J = 1, 4)
READ(11,1700) ( ( LZX( I,J), I=1, NZX ), J = 1, 4)
1700 FORMAT(16I5)
C
C ----- To set Node numbers for U-V Distribution on Z = 0.5 ---
C
CALL STND ( X, Y, LXY, NXY, NNP )
C
C -------------------------------------------------------------
C
RETURN
END
C **********************************************************************
C
SUBROUTINE MSHF ( XLX, YLY, LAX, LAY, NEX, MESH )
C
IMPLICIT REAL*8 (A-H,O-Z)
C
DIMENSION X ( 50), Y ( 50), GX ( 50), GY ( 50),
1 GXR( 50), GYR( 50), GZX( 50), GZY( 50),
2 ZXR( 50), ZYR( 50), XRA( 50), XRB( 50),
3 YRA( 50), YRB( 50), XLX( LAX), YLY( LAY)
C
COMMON /VWCT/ HR, VR, DIS, DS2, OX, OY, OZ, SCA, RR
C
INTEGER*2 ICR, IY0, IX0, IX1, IY1, IX2, IY2, IX3, IY3, IX4,
1 IY4, IY5, IY6, IWK1, IWK2, IWK3, IWK4
C
CALL SET_TEXT_ATTRIBUTE@( 3, 1.5, 0., 0. )
C
CALL DRAW_TEXT@('HF3D-FEM-R', 25, 430, 15 )
C
IF ( NEX.EQ. 2) CALL DRAW_TEXT@( ' ( 2 X 2 X 2 )', 140,430,15)
IF ( NEX.EQ. 3) CALL DRAW_TEXT@( ' ( 3 X 3 X 3 )', 140,430,15)
IF ( NEX.EQ. 4) CALL DRAW_TEXT@( ' ( 4 X 4 X 4 )', 140,430,15)
C
IF ( NEX.EQ.10) CALL DRAW_TEXT@( ' (10 X 10 X 10)',140,430,15)
IF ( NEX.EQ.20) CALL DRAW_TEXT@( ' (20 X 20 X 20)',140,430,15)
IF ( NEX.EQ.30) CALL DRAW_TEXT@( ' (30 X 30 X 30)',140,430,15)
C
IF ( NEX.LE.10) CALL DRAW_TEXT@( 'Ra = 710', 335,430,15)
IF ( NEX.EQ.20) CALL DRAW_TEXT@( 'Ra = 7100', 340,430,15)
IF ( NEX.EQ.30) CALL DRAW_TEXT@( 'Ra = 10** 6',340,430,15)
C
C ----------------------------
DO 100 I = 1, NEX
C
X( I) = XLX( I)
Y( I) = YLY( I)
C
100 CONTINUE
C ----------------------------
C
ICR = 15
IWD = 200
C
IY0 = 30
IX1 = 100
IY1 = IY0 + 100
IX2 = IX1 + 67
IY2 = 50 + IY0
C
IX3 = IX2 + IWD
IY3 = IY2
IX4 = IX1 + IWD
IY4 = IY1
IY5 = IY1 + IWD
IY6 = IY3 + IWD
C
CALL DRAW_LINE@( IX1, IY1, IX4, IY4, ICR )
CALL DRAW_LINE@( IX1, IY1, IX1, IY5, ICR )
C
CALL DRAW_LINE@( IX1, IY5, IX4, IY5, ICR )
CALL DRAW_LINE@( IX4, IY5, IX3, IY6, ICR )
C
CALL DRAW_LINE@( IX1, IY1, IX2, IY2, ICR )
CALL DRAW_LINE@( IX2, IY2, IX3, IY3, ICR )
C
CALL DRAW_LINE@( IX3, IY3, IX3, IY6, ICR )
CALL DRAW_LINE@( IX4, IY1, IX4, IY5, ICR )
C
CALL DRAW_LINE@( IX4, IY4, IX4, IY1, ICR )
CALL DRAW_LINE@( IX4, IY1, IX3, IY2, ICR )
C
CALL DRAW_LINE@( IX2, IY2, IX2, IY0, ICR )
C
C ----- Arrows ---------------------------------------
C
IX0 = IX2 - 5
IWK2 = 20 + IY0
CALL DRAW_LINE@( IX0, IWK2, IX2, IY0, ICR )
C
CALL SET_TEXT_ATTRIBUTE@( 3, 1.5, 0., 0. )
C
CALL DRAW_TEXT@( ' y', IX2-50, 30, 15 )
C
IWK4 = IX2 + 5
CALL DRAW_LINE@( IX2, IY0, IWK4, IWK2, ICR )
C
IWK3 = IX1 - 40
IWK4 = IY5 + 30
CALL DRAW_LINE@( IX1, IY5, IWK3, IWK4, ICR )
C
IWK1 = IY1 + IWD + 30
IWK2 = IX1 - 33
IWK4 = IX1 + IWD + 15 + IY0
CALL DRAW_LINE@( IWK3, IWK1, IWK2, IWK4, ICR )
C
CALL DRAW_TEXT@( ' z', IWK3-40, IWK1+20, 15 )
C
IWK2 = IX1 - 23
IWK4 = IY1 + IWD + 28
CALL DRAW_LINE@( IWK3, IWK1, IWK2, IWK4, ICR )
C
IWK1 = IY3 + IWD
IWK2 = IX1 + IWD + 140
CALL DRAW_LINE@( IX3, IWK1, IWK2, IWK1, ICR )
C
IWK1 = IX1 + IWD + 140
IWK2 = IY3 + IWD
IWK3 = IWK1 - 20
IWK4 = IWK2 + 5
CALL DRAW_LINE@( IWK1, IWK2, IWK3, IWK4, ICR )
C
IWK4 = IY3 + IWD - 5
CALL DRAW_LINE@( IWK1, IWK2, IWK3, IWK4, ICR )
C
CALL DRAW_TEXT@( ' x', IWK1-10, IWK2+15, 15 )
C
C --------------------------------------
DO 200 I = 1, NEX
C
L = I
GX( L) = 200.D0* X( L)
GY( L) = 200.D0* Y( L)
C
200 CONTINUE
C --------------------------------------
IWK1 = IY1 + IWD
DO 300 L = 2, NEX
C
GXR( L) = 100.D0 + GX( L)
IX1 = IDNINT( GXR( L))
CALL DRAW_LINE@( IX1, IWK1, IX1, IY1, ICR )
C
300 CONTINUE
C --------------------------------------
IWK1 = 100 + IY0
IWK2 = 50 + IY0
DO 400 I = 2, NEX
C
GXR( I)= 100.D0 + GX( I)
IX1 = IDNINT( GXR( I))
IX2 = 167 + GX( I)
CALL DRAW_LINE@( IX1, IWK1, IX2, IWK2, ICR )
C
400 CONTINUE
C --------------------------------------
IWK1 = 100
IWK2 = 100 + IWD
DO 500 I = 2, NEX
C
GYR( I) = 100.D0 + GY( I)
IY1 = IDNINT( 100.D0 + GY( I)) + IY0
CALL DRAW_LINE@( IWK1, IY1, IWK2, IY1, ICR )
C
500 CONTINUE
C --------------------------------------
C
IWK1 = 300
IWK2 = 367
C --------------------------------------
DO 600 I = 2, NEX
C
GYR( I) = 100.D0 + GY( I)
IY1 = IDNINT( 100.D0 + GY( I)) + IY0
IY2 = 47 + GY( I) + IY0
CALL DRAW_LINE@( IWK1, IY1, IWK2, IY2, ICR )
C
600 CONTINUE
C --------------------------------------
DO 700 I = 2, NEX
C
GZX( I) = 0.8D0* GX( I)
GZY( I) = 0.6D0* GX( I)
C
700 CONTINUE
C --------------------------------------
C
** WRITE(6,*) ' * IY0 =', IY0
** WRITE(6,*) ' * MESH =', MESH
C
C ----------------------------------------------------------
IF ( MESH.EQ.1) THEN
C
C ------------------------------------------------
DO 800 L = 2, NEX
C
ZXR( L) = 367.D0 - GZX( L)* 0.42D0
YRA( L) = 250.D0 + GZY( L)* 0.42D0
YRB( L) = 50.D0 + GZY( L)* 0.42D0
IX1 = ZXR( L)
IY1 = YRA( L) + IY0
IY2 = YRB( L) + IY0
CALL DRAW_LINE@( IX1, IY1, IX1, IY2, ICR )
C
800 CONTINUE
C ------------------------------------------------
DO 900 L = 2, NEX
C
XRA( L) = 167.D0 - GZX( L)* 0.42D0
XRB( L) = 367.D0 - GZX( L)* 0.42D0
ZYR( L) = 50.D0 + GZY( L)* 0.42D0
IX1 = IDNINT( XRA( L))
IY1 = IDNINT( ZYR( L)) + IY0
IX2 = IDNINT( XRB( L))
CALL DRAW_LINE@( IX1, IY1, IX2, IY1, ICR )
C
900 CONTINUE
C ------------------------------------------------
C
END IF
C ----------------------------------------------------------
IF ( MESH.EQ.2 ) THEN
C
C -----------------------------------------------
DO 920 L = 2, NEX
C
ZXR( L) = 367.D0 - Gzx( L)* 0.4D0
YRA( L) = 250.D0 + Gzy( L)* 0.4D0
YRB( L) = 50.D0 + Gzy( L)* 0.4D0
IX1 = ZXR( L)
IY1 = IDNINT( YRA( L)) + IY0
IY2 = IDNINT( YRB( L)) + IY0
CALL DRAW_LINE@( IX1, IY1, IX1, IY2, ICR )
C
920 CONTINUE
C ------------------------------------------------
DO 950 l = 2, NEX
C
XRA( L) = 167.D0 - GZX( L)* 0.4D0
XRB( L) = 367.D0 - GZX( L)* 0.4D0
ZYR( L) = 50.D0 + GZY( L)* 0.4D0
IX1 = IDNINT( XRA( L))
IY1 = IDNINT( ZYR( L)) + IY0
IX2 = IDNINT( XRB( L))
CALL DRAW_LINE@( IX1, IY1, IX2, IY1, ICR )
C
950 CONTINUE
C ------------------------------------------------
C
ENDIF
C ----------------------------------------------------------
RETURN
END
C **********************************************************************
C
SUBROUTINE PGEO ( X, Y, Z, NNP )
C
IMPLICIT REAL*8 (A-H,O-Z)
C
INTEGER*2 ICR, IX1, IY1, IX2, IY2, IX3, IY3, IX4, IY4, IX5, IY5,
1 IX6, IY6, IX7, IY7, IX8, IY8
C
COMMON /VWCT/ HR, VR, DIS, DS2, OX, OY, OZ, SCA, RR
C
DIMENSION X( NNP), Y( NNP), Z( NNP)
C
ICR = 15
C
XMIN = X( 1)
XMAX = X( 1)
C
YMIN = Y( 1)
YMAX = Y( 1)
C
ZMIN = Z( 1)
ZMAX = Z( 1)
C
C -------------------------------------
DO 100 I = 1, NNP
C
XMIN = DMIN1( XMIN, X( I))
XMAX = DMAX1( XMAX, X( I))
C
YMIN = DMIN1( YMIN, Y( I))
YMAX = DMAX1( YMAX, Y( I))
C
ZMIN = DMIN1( ZMIN, Z( I))
ZMAX = DMAX1( ZMAX, Z( I))
C
100 CONTINUE
C --------------------------------------
C
OX = ( XMIN + XMAX)/ 2.D0
OY = ( YMIN + YMAX)/ 2.D0
OZ = ( ZMIN + ZMAX)/ 2.D0
C
CALL VPNT ( XMIN, YMIN, ZMIN, X1, Y1 )
CALL VPNT ( XMAX, YMIN, ZMIN, X2, Y2 )
CALL VPNT ( XMAX, YMAX, ZMIN, X3, Y3 )
CALL VPNT ( XMIN, YMAX, ZMIN, X4, Y4 )
C
CALL VPNT ( XMIN, YMIN, ZMAX, X5, Y5 )
CALL VPNT ( XMAX, YMIN, ZMAX, X6, Y6 )
CALL VPNT ( XMAX, YMAX, ZMAX, X7, Y7 )
CALL VPNT ( XMIN, YMAX, ZMAX, X8, Y8 )
C
IY1 = 399 - Y1
IY2 = 399 - Y2
IY3 = 399 - Y3
IY4 = 399 - Y4
C
IY5 = 399 - Y5
IY6 = 399 - Y6
IY7 = 399 - Y7
IY8 = 399 - Y8
C
ICT = 150
C
IX1 = X1 + ICT
IX2 = X2 + ICT
IX3 = X3 + ICT
IX4 = X4 + ICT
C
IX5 = X5 + ICT
IX6 = X6 + ICT
IX7 = X7 + ICT
IX8 = X8 + ICT
C
CALL DRAW_LINE@( IX1, IY1, IX2, IY2, ICR )
CALL DRAW_LINE@( IX2, IY2, IX3, IY3, ICR )
CALL DRAW_LINE@( IX3, IY3, IX4, IY4, ICR )
CALL DRAW_LINE@( IX4, IY4, IX1, IY1, ICR )
C
CALL DRAW_LINE@( IX1, IY1, IX5, IY5, ICR )
CALL DRAW_LINE@( IX5, IY5, IX6, IY6, ICR )
CALL DRAW_LINE@( IX6, IY6, IX7, IY7, ICR )
CALL DRAW_LINE@( IX7, IY7, IX8, IY8, ICR )
CALL DRAW_LINE@( IX8, IY8, IX5, IY5, ICR )
C
CALL DRAW_LINE@( IX8, IY8, IX4, IY4, ICR )
CALL DRAW_LINE@( IX4, IY4, IX3, IY3, ICR )
CALL DRAW_LINE@( IX3, IY3, IX7, IY7, ICR )
CALL DRAW_LINE@( IX7, IY7, IX6, IY6, ICR )
CALL DRAW_LINE@( IX6, IY6, IX2, IY2, ICR )
C
C ----- PLOT X - Y - Z - > ---
C
XP = XMIN - 0.1D0
YP = YMIN - 0.1D0
ZP = ZMAX
C
XC = XP
YC = YP
ZC = ZP
C
CALL VPNT ( XP, YP, ZP, X1, Y1)
CALL VPNT ( XP + 0.2D0, YP, ZP, X2, Y2 )
CALL VPNT ( XP, YP + 0.2D0, ZP, X3, Y3 )
C
CALL VPNT ( XP, YP, ZP + 0.2D0, X4, Y4 )
CALL VPNT ( XP, YP, ZP-0.1D0, X5, Y5 )
C
IX1 = X1
IY1 = 399 - Y1
IX2 = X2
IY2 = 399 - Y2
C
IX3 = X3
IY3 = 399 - Y3
IX4 = X4
IY4 = 399 - Y4
C
IX1 = IX1 + ICT
IX2 = IX2 + ICT
IX3 = IX3 + ICT
IX4 = IX4 + ICT
C
CALL DRAW_LINE@( IX1, IY1, IX2, IY2, ICR )
CALL DRAW_LINE@( IX1, IY1, IX3, IY3, ICR )
CALL DRAW_LINE@( IX1, IY1, IX4, IY4, ICR )
C
C ----- PLOT 'X' ---
C
XP = XC + 0.25D0
YP = YC
ZP = ZC
C
CALL VPNT ( XP-0.02D0, YP+0.02D0, ZP, X1, Y1)
CALL VPNT ( XP+0.02D0, YP-0.02D0, ZP, X2, Y2 )
CALL VPNT ( XP-0.02D0, YP-0.02D0, ZP, X3, Y3 )
CALL VPNT ( XP+0.02D0, YP+0.02D0, ZP, X4, Y4 )
C
IX1 = X1
IY1 = 399-Y1
IX2 = X2
IY2 = 399-Y2
C
IX3 = X3
IY3 = 399-Y3
IX4 = X4
IY4 = 399-Y4
C
IX1 = IX1 + ICT
IX2 = IX2 + ICT
IX3 = IX3 + ICT
IX4 = IX4 + ICT
C
CALL DRAW_LINE@( IX1, IY1, IX2, IY2, ICR )
CALL DRAW_LINE@( IX3, IY3, IX4, IY4, ICR )
C
C ----- PLOT 'Y' ---
C
XP = XC
YP = YC + 0.25D0
ZP = ZC
C
CALL VPNT ( XP, YP, ZP, X1, Y1)
CALL VPNT ( XP-0.02D0, YP+0.03D0, ZP, X2, Y2 )
CALL VPNT ( XP+0.02D0, YP+0.03D0, ZP, X3, Y3 )
CALL VPNT ( XP, YP-0.02D0, ZP, X4, Y4 )
C
IY1 = 399 - Y1
IY2 = 399 - Y2
C
IY3 = 399 - Y3
IY4 = 399 - Y4
C
IX1 = X1 + ICT
IX2 = X2 + ICT
IX3 = X3 + ICT
IX4 = X4 + ICT
C
CALL DRAW_LINE@( IX2, IY2, IX1, IY1, ICR )
CALL DRAW_LINE@( IX1, IY1, IX3, IY3, ICR )
CALL DRAW_LINE@( IX3, IY3, IX1, IY1, ICR )
CALL DRAW_LINE@( IX1, IY1, IX4, IY4, ICR )
C
C ----- PLOT 'Z' ---
C
XP = XC
YP = YC
ZP = ZC + 0.25D0
C
CALL VPNT ( XP-0.02D0, YP+0.02D0, ZP, X1, Y1)
CALL VPNT ( XP+0.02D0, YP+0.02D0, ZP, X2, Y2 )
CALL VPNT ( XP-0.02D0, YP-0.02D0, ZP, X3, Y3 )
CALL VPNT ( XP+0.02D0, YP-0.02D0, ZP, X4, Y4 )
C
IY1 = 399 - Y1
IY2 = 399 - Y2
IY3 = 399 - Y3
IY4 = 399 - Y4
C
IX1 = X1 + ICT
IX2 = X2 + ICT
IX3 = X3 + ICT
IX4 = X4 + ICT
C
CALL DRAW_LINE@( IX1, IY1, IX2, IY2, ICR )
CALL DRAW_LINE@( IX2, IY2, IX3, IY3, ICR )
CALL DRAW_LINE@( IX3, IY3, IX4, IY4, ICR )
C
RETURN
END
C *********************************************************************
C
SUBROUTINE PLTN ( X, Y, UN, VN, NODE, NNP, NEL, IKY )
C
IMPLICIT REAL*8 (A-H,O-Z)
C
DIMENSION X( NNP), Y( NNP), UN( NNP), VN( NNP), NODE( NEL,4)
C
INTEGER*2 ICR, IX1, IY1, IX2, IY2, IX3, IY3
C
COMMON /PFAC/ ZFAC
C
ICR = 15
C
IWD = 160
IX1 = 20
IF ( IKY.EQ.4) IX1 = IX1 + 200
IF ( IKY.EQ.5) IX1 = IX1 + 400
IY1 = 20
IF ( IKY.EQ.4) IY1 = 20 + IDNINT( DFLOAT( IWD)* ( 1.D0 - ZFAC))
IX3 = IX1 + IWD
IY3 = IY1 + IWD
IF ( IKY.EQ.4) IY3 = IY1 + IDNINT( DFLOAT( IWD)* ZFAC)
IF ( IKY.EQ.5) IX3 = IX1 + IDNINT( DFLOAT( IWD)* ZFAC)
C
CALL RECTANGLE@( IX1, IY1, IX3, IY3, ICR )
C
ICR = 11
C
X0 = 0.1D0
Y0 = 0.1D0
C
SCX = 0.8D0
SCY = 0.8D0
SCL = 0.02D0
SCL = 0.04D0
C
IF ( NNP.GE. 9000) SCL = 0.005D0
IF ( NNP.GE.29000) SCL = 0.001D0 ! 31*31*31
C -----------------------------------------------------
C
EPS = 1.D-4
C
C ----------------------------------------------------------
DO 100 J = 1, NEL
C
I = NODE( J, 1)
GX = X0 + X( I)* SCX
GY = Y0 + Y( I)* SCY
IF ( DABS( GX - X0).LE.0.0001D0) GO TO 100
IF ( DABS( GY - Y0).LE.0.0001D0) GO TO 100
IF ( DABS( GX - X0 - SCX).LE.0.0001D0) GO TO 100
IF ( DABS( GY - Y0 - SCY).LE.0.0001D0) GO TO 100
C
AL = DSQRT( UN( I)**2 + VN( I)**2 )
C ------------------------------------------------
IF ( AL.GE.EPS) THEN
C
XS = GX
YS = GY
XE = XS + UN( I)* SCL
YE = YS + VN( I)* SCL
IX1 = 200* XS
IF ( IKY.EQ.4) IX1 = IX1 + 200
IF ( IKY.EQ.5) IX1 = IX1 + 400
IY1 = 200 - 200* YS
IX2 = 200* XE
IF ( IKY.EQ.4) IX2 = IX2 + 200
IF ( IKY.EQ.5) IX2 = IX2 + 400
IY2 = 200.D0 - 200.D0* YE
C
CALL DRAW_LINE@( IX1, IY1, IX2, IY2, ICR )
C
END IF
C ------------------------------------------------
C
100 CONTINUE
C ----------------------------------------------------------
C
RETURN
END
C **********************************************************************
C
SUBROUTINE PLTT ( X, Y, NODE, TEP, NNP, NEL, IKY )
C
C ----- To Plot Contour Lines for Rectangular Meshes ---
C
C --- ( Linear interpolation ) ---
C
IMPLICIT REAL*8 (A-H,O-Z)
C
DIMENSION X ( NNP), Y ( NNP), NODE( NEL,4), TEP( NNP),
1 PX( 4), PY( 4)
C
INTEGER*2 ICR, IX1, IY1, IX2, IY2, IX3, IY3
C
COMMON /PFAC/ ZFAC
COMMON /PLOT/ HVN( 20)
C
CHARACTER*11 TSB( 3)
DATA TSB / '( Z = 0.5 )', '( X = 0.5 )', '( Y = 0.5 )'/
C
ICR = 15
C
X0 = 0.1D0
Y0 = 0.1D0
SCX = 0.8D0
SCY = 0.8D0
C
NM = 9
TMIN = 0.D0
TMAX = 1.D0
C
C --------------------------------------------------------------
DO 100 I = 1, NM
C
HVN( I) = TMIN + ( TMAX-TMIN )* DFLOAT( I)/ DFLOAT( NM+1)
100 CONTINUE
C --------------------------------------------------------------
C
IWD = 160
IX1 = 20
IF ( IKY.EQ.4) IX1 = IX1 + 200
IF ( IKY.EQ.5) IX1 = IX1 + 400
IY1 = 220
IF ( IKY.EQ.4) IY1 = 220 + DFLOAT( IWD)* ( 1.D0-ZFAC)
IX3 = IX1 + IWD
IY3 = IY1 + IWD
IF ( IKY.EQ.4) IY3 = IY1 + DFLOAT( IWD)* ZFAC
IF ( IKY.EQ.5) IX3 = IX1 + DFLOAT( IWD)* ZFAC
C
CALL RECTANGLE@( IX1, IY1, IX3, IY3, ICR )
C
ICR = 12
IWD = 200
C
PXV = 0.D0
PYV = 0.D0
C ----------------------------------------------------------
DO 200 L = 1, NEL
C
C -----------------------------------------------
DO 220 IH = 1, NM
C
HV = HVN( IH)
ICL = 0
C
C --------------------------------------
DO 240 J = 1, 4
C
I1 = NODE( L,J)
IF ( J .NE. 4) THEN
I2 = NODE( L,J+1)
ELSE
I2 = NODE( L,1)
END IF
C
IF ( DABS( TEP( I1)-TEP( I2)).LE.0.0001D0) THEN
IF ( DABS( TEP( I1)-HV).LE.0.0001D0) THEN
C
ICL = ICL + 1
C
PX( ICL) = X( I1)* SCX
PY( ICL) = Y( I1)* SCY
C
END IF
C
ELSE IF ( ( ( TEP( I1).LE.HV).AND.( HV.LT.TEP( I2)) )
1 .OR. ( ( TEP( I1).GE.HV).AND.( HV.GT.TEP( I2)) ))
2 THEN
ICL = ICL + 1
RT = ( HV - TEP( I1))/( TEP( I2) - TEP( I1))
PX( ICL) = ( X( I1) + ( X( I2) - X( I1))* RT)* SCX
PY( ICL) = ( Y( I1) + ( Y( I2) - Y( I1))* RT)* SCY
END IF
C
240 CONTINUE
C --------------------------------------
C
IF ( ICL.GE.2 ) THEN
C
C --------------------------------------
DO 260 J = 1, ICL
C
JM = J
D2 = ( PXV - PX( J))** 2 + ( PYV - PY( J))** 2
C
C ----------------------------
C
DO 280 JJ = J+1, ICL
C
D2V = ( PXV - PX( JJ))** 2 + ( PYV - PY( JJ))** 2
IF ( D2V .LT. D2 ) THEN
JM = JJ
D2 = D2V
END IF
C
280 CONTINUE
C ----------------------------
PXV = PX( JM)
PYV = PY( JM)
PX(JM) = PX( J)
PY(JM) = PY( J)
PX( J) = PXV
PY( J) = PYV
C
260 CONTINUE
C --------------------------------------
PXV = PX( ICL)
PYV = PY( ICL)
END IF
C
C -----------------------------------------------------
IF ( ICL.GT.0) THEN
C ------------------------------------------------
DO 290 J = 1, ICL-1
C
IX1 = 200.D0* ( X0 + PX( J))
IF ( IKY.EQ.4) IX1 = IX1 + 200
IF ( IKY.EQ.5) IX1 = IX1 + 400
IY1 = 401 - 200.D0* ( Y0 + PY( J))
IX2 = 200.D0* ( X0 + PX( J+1))
IF ( IKY.EQ.4) IX2 = IX2 + 200
IF ( IKY.EQ.5) IX2 = IX2 + 400
IY2 = 401 - 200.D0*( Y0 + PY( J+1))
CALL DRAW_LINE@( IX1, IY1, IX2, IY2, ICR )
C
290 CONTINUE
C ------------------------------------------------
DO 295 J = 2, ICL
C
IX1 = 200.D0* ( X0 + PX( 1))
IF ( IKY.EQ.4) IX1 = IX1 + 200
IF ( IKY.EQ.5) IX1 = IX1 + 400
IY1 = 401 - 200.D0* ( Y0 + PY( 1))
IX2 = 200.D0* ( X0 + PX( J))
IF ( IKY.EQ.4) IX2 = IX2 + 200
IF ( IKY.EQ.5) IX2 = IX2 + 400
IY2 = 401 - 200.D0* ( Y0 + PY( J))
C
295 CONTINUE
C ------------------------------------------------
END IF
C -----------------------------------------------------
C
220 CONTINUE
C ------------------------------------------------
C
200 CONTINUE
C ----------------------------------------------------------
C
ICR = 15
C
IY1 = 240 + IWD
IX1 = 10
IX2 = IX1 + 40
IY2 = IY1 - 40
CALL DRAW_LINE@( IX1, IY1, IX2, IY1, ICR )
CALL DRAW_LINE@( IX1, IY1, IX1, IY2, ICR )
C
CALL SET_TEXT_ATTRIBUTE@( 3, 1.3, 0., 0. )
CALL DRAW_TEXT@( 'X', 60, 450, ICR )
CALL DRAW_TEXT@( 'Y', 20, 410, ICR )
IX1 = 10 + IWD
IY1 = IY1 - 40
IX2 = IX1 + 40
IY2 = IY1 + 40
CALL DRAW_LINE@( IX1, IY1, IX2, IY1, ICR )
CALL DRAW_LINE@( IX1, IY1, IX1, IY2, ICR )
CALL DRAW_TEXT@( 'Z', IX1+10, 450, ICR )
CALL DRAW_TEXT@( 'Y', IX1+50, 410, ICR )
C
IX1 = 10 + IWD + IWD + 140
IX2 = IX1 + 40
IY2 = IY1 + 40
CALL DRAW_LINE@( IX1, IY1, IX2, IY1, ICR )
CALL DRAW_LINE@( IX2, IY1, IX2, IY2, ICR )
CALL DRAW_TEXT@( 'Z', IX1-20, 410, ICR )
CALL DRAW_TEXT@( 'X', IX1+20, 450, ICR )
C
C --------------------------------------------------------
C
CALL SET_TEXT_ATTRIBUTE@( 3, 1.3, 0., 0. )
C
CALL DRAW_TEXT@( TSB( 1), 50, 210, ICR )
C
CALL DRAW_TEXT@( TSB( 2), 50 + IWD, 210, ICR )
C
CALL DRAW_TEXT@( TSB( 3), 50 + 2* IWD, 210, ICR )
C
C --------------------------------------------------------
C
RETURN
END
C **********************************************************************
C
SUBROUTINE PLVR ( XX, KEY, RA )
C
IMPLICIT REAL*8 (A-H,O-Z)
C
INTEGER*2 ICR, IX0, IY0, IX1, IY1
C
C ------ Attention : Size ---
C
DIMENSION XX( 8000)
C
INTEGER*2 IVAR, IXX
DIMENSION IVAR( 8000), IXX( 8000)
C
COMMON /CNST/ IKY, MSH, LOP
C
CHARACTER*80 STR0, STR4, STR5
CHARACTER*4 STR1( 3)
C
DATA STR0
1/'*** Result for HF3D-FEM-R ** ( Max_Var_Vel.(%)) & Max_Vel.***'/
C
DATA STR4
1/'*** Result for HF3D-FEM-R ** ( Max_Var_Vel.(%)) & Max_Vel.( x10
200)***'/
C
DATA STR5
1/'*** Result for HF3D-FEM-R ** ( Max_Var_Vel.(%)) & Max_Vel.( x10
2000)***'/
C
DATA STR1 /'0 ', 'LP/2', 'LOP'/
C
IX0 = 50
IY0 = 40
IF ( KEY.NE.1) GO TO 90
ICR = 15
C
CALL SET_TEXT_ATTRIBUTE@( 3, 1., 0., 0. )
C
IF ( DABS( RA-710.D0).LE.0.01D0)
1 CALL DRAW_TEXT@( STR0, 50, 10, ICR )
C
IF ( DABS( RA-7100.D0).LE.0.01D0)
1 CALL DRAW_TEXT@( STR4, 50, 10, ICR )
C
IF ( DABS( RA-1000000.D0).LE.0.01D0)
1 CALL DRAW_TEXT@( STR5, 50, 10, ICR )
C
CALL DRAW_TEXT@( STR1( 1), 50, 460, ICR )
CALL DRAW_TEXT@( STR1( 2), 320, 460, ICR )
CALL DRAW_TEXT@( STR1( 3), 580, 460, ICR )
C
CALL DRAW_TEXT@( '1000', 10, 30, ICR )
CALL DRAW_TEXT@( ' 100', 10, 100, ICR )
CALL DRAW_TEXT@( ' 10', 10, 170, ICR )
CALL DRAW_TEXT@( ' 1', 10, 240, ICR )
C
CALL DRAW_TEXT@( ' 0.1', 10, 310, ICR )
CALL DRAW_TEXT@( '0.01', 10, 380, ICR )
CALL DRAW_TEXT@( '0.001', 10, 450, ICR )
C
CALL SET_TEXT_ATTRIBUTE@( 3, 1.2, 0., 0. )
CALL DRAW_TEXT@( '(%)', 5, 255, ICR )
C
CALL SET_TEXT_ATTRIBUTE@( 3, 1., 0., 0. )
CALL DRAW_TEXT@( '5', 615, 31, ICR )
CALL DRAW_TEXT@( '4', 615, 115, ICR )
CALL DRAW_TEXT@( '3', 615, 199, ICR )
C
CALL DRAW_TEXT@( '2', 615, 283, ICR )
CALL DRAW_TEXT@( '1', 615, 367, ICR )
CALL DRAW_TEXT@( '0', 615, 451, ICR )
C
IWD = 560
IX1 = IX0 + IWD
IY1 = IY0 + 420
CALL RECTANGLE@( IX0, IY0, IX1, IY1, 15 )
C
90 CONTINUE
C
C ----- XX( I) : Max. Velocity ---
C
RIWD = DFLOAT( IWD)/ DFLOAT( LOP)
C -----------------------------------------------------------
C
IF ( KEY.EQ.2 ) THEN
IF ( DABS( RA- 710.D0) .LE.0.01D0) SCL = 84.D0
IF ( DABS( RA-7100.D0) .LE.0.01D0) SCL = 0.084D0
IF ( DABS( RA-1000000.D0).LE.0.01D0) SCL = 0.0084D0 ! Ra=10**6
END IF
C
C -----------------------------------------------------------
C
C -----------------------------
IF ( KEY.EQ.1) THEN
C
DO 100 I= 1, LOP
C
XX( I) = DLOG10( XX( I))
C
100 CONTINUE
C
SCL = 70.D0
C
END IF
C -----------------------------
C
C ---------------------------------------------------------------
DO 200 I = 1, LOP
C
IF ( KEY.EQ.1) IVAR( I) = 250 - IDNINT( XX( I)* SCL)
C
IF ( KEY.EQ.2 ) IVAR( I) = IY0 + 420 - IDNINT( XX( I)* SCL)
C
IXX( I) = IX0 + IDNINT( DFLOAT( I)* RIWD )
C ==========================================
C
200 CONTINUE
C --------------------------------------------------------------
IF ( KEY.EQ.1) ICR = 9
IF ( KEY.EQ.2) ICR = 11
C
CALL POLYLINE@( IXX, IVAR, LOP, ICR )
C
RETURN
END
C **********************************************************************
C
SUBROUTINE PVLTM ( U, V, W, T, X, Y, Z, LXY, LYZ,
1 LZX, NNP, NXY, NYZ, NZX )
C
IMPLICIT REAL*8 (A-H,O-Z)
C
DIMENSION U( NNP), V( NNP), W( NNP), T( NNP), X( NNP), Y( NNP),
1 Z( NNP)
C
DIMENSION LXY( 100,4), LYZ( 100,4), LZX( 100,4)
C
COMMON /CNST/ IKY, MSH, LOP
C
C ----------------------------------------------------------
C
CALL PLTT ( X, Y, LXY, T, NNP, NXY, 3 )
CALL PLTN ( X, Y, U, V, LXY, NNP, NXY, 3 )
C
CALL PLTT ( Y, Z, LYZ, T, NNP, NYZ, 4)
CALL PLTN ( Y, Z, V, W, LYZ, NNP, NYZ, 4)
C
CALL PLTT ( Z, X, LZX, T, NNP, NZX, 5 )
CALL PLTN ( Z, X, W, U, LZX, NNP, NZX, 5 )
C
C ----------------------------------------------------------
C
RETURN
END
C **********************************************************************
C
SUBROUTINE RFCP ( X, Y, U, V, NNP, RA )
C
IMPLICIT REAL*8(A-H,O-Z)
C
DIMENSION X( NNP), Y( NNP), U( NNP), V( NNP)
C
COMMON /CMPD/ KSY( 200), NY, KSX( 200), NX
C
INTEGER*2 X0, Y0, X1, Y1, IX1, IY1, IX2, IY2, ICR, IC2
C
CHARACTER*80 STR0
C
DATA STR0 /'** Results for HF3D-FEM ** ( Velocity U - V Distribut
1ion on Z = 0.5 )'/
C
NPL = 0
ICR = 11
IC2 = 16
C
CALL SET_TEXT_ATTRIBUTE@( 2.5, 1.5, 0., 0. )
C
CALL DRAW_TEXT@( STR0, 10, 5, IC2 )
C
X1 = 345
Y1 = 50
X0 = 80
Y0 = 50
C
IWD = 400
SCL = 10.D0
C
IX1 = X0
IY1 = Y0
IX2 = X0 + IWD
IY2 = Y0 + IWD
CALL RECTANGLE@( IX1, IY1, IX2, IY2, IC2 )
C
IX1 = X0 + IDNINT( DFLOAT( IWD)* 0.5D0)
IY1 = Y0 + IWD
IX2 = IX1
IY2 = Y0
CALL RECTANGLE@( IX1, IY1, IX2, IY2, IC2 )
C
IX1 = X0
IY1 = Y0 + IDNINT( DFLOAT( IWD)* 0.5D0)
IX2 = IX1 + IWD
IY2 = IY1
CALL DRAW_LINE@( IX1, IY1, IX2, IY2, IC2 )
C
C ----- For Ra < = 7100 -------------------------------
C
IF ( RA.LE.7200) THEN
C
CALL SET_TEXT_ATTRIBUTE@( 2.5, 1.5, 0., 0. )
C
CALL DRAW_TEXT@(' 20', 450, 475, IC2 )
CALL DRAW_TEXT@('-20', 50, 475, IC2 )
CALL DRAW_TEXT@(' 0', 250, 475, IC2 )
C
CALL DRAW_TEXT@(' 20', 40, 70, IC2 )
CALL DRAW_TEXT@(' 0', 40, 270, IC2 )
C
END IF
C --------------------------------------------------------------
IF ( DABS(RA-1000000.D0).LE.0.01D0) THEN ! Ra = 10**6
C
CALL SET_TEXT_ATTRIBUTE@( 2.5, 1., 0., 0. )
C
CALL DRAW_TEXT@(' 100', 450, 475, IC2 )
CALL DRAW_TEXT@('-100', 50, 475, IC2 )
CALL DRAW_TEXT@(' 0', 240, 475, IC2 )
C
CALL DRAW_TEXT@(' 300', 20, 60, IC2 )
CALL DRAW_TEXT@('-300', 20, 455, IC2 )
CALL DRAW_TEXT@(' 0', 20, 260, IC2 )
C
END IF
C ---------------------------------------------------------------
IX1 = 40
IY1 = 80 + NPL* 20
IX2 = 80
IY2 = IY1
*** WRITE(6,*) ' ** NX =', NX, ' NY =', NY
C
C ----- Computed U ---------------------------------------------------
C
IF ( DABS( RA-1000000.D0).LE.0.001D0) SCL= 2.D0
C ======================================================
C
C --------------------------------------------------------------------
DO 300 J = 1, NY
C
NN = KSY( J)
IX1 = X0 + IDNINT( DFLOAT( IWD)* 0.5D0) + IDNINT( U( NN)* SCL)
C -------------------------------------------------------------------
IY1 = Y0 + IWD - IDNINT( Y( NN)* DFLOAT( IWD))
IX2 = 2
IY2 = 2
IF ( J.EQ.1) GO TO 350
CALL ELLIPSE@( IX1, IY1, IX2, IY2, ICR )
C
350 CONTINUE
NN1 = KSY( J+1)
IX2 = X0 + 10 + IDNINT( U( NN1)* SCL)
IY2 = Y0 + IWD - IDNINT( Y( NN1)* DFLOAT( IWD))
IF ( J.EQ.NY) THEN
IX2 = X0 + IWD
IY2 = Y0
END IF
C
300 CONTINUE
C ---------------------------------------------------------------------
C
C ----- Computed V ----------------------------------------------------
C
IF ( DABS( RA-1000000.D0).LE.0.001D0) SCL = 0.7D0
C
DO 500 J = 1, NX
C
NN = KSX( J)
IX1 = X0 + IDNINT( X( NN)* DFLOAT( IWD))
IY1 = Y0 + IDNINT( DFLOAT( IWD)*0.5D0) - IDNINT( V( NN)* SCL)
C --------------------------------------------------------------------
IX2 = 2
IY2 = 2
IF ( J.EQ.1) GO TO 550
CALL FILL_ELLIPSE@( IX1, IY1, IX2, IY2, ICR )
C
550 CONTINUE
NN1 = KSX( J+1)
IX2 = X0 + IDNINT( X( NN1)* DFLOAT( IWD))
IY2 = Y0 + IDNINT( DFLOAT( IWD)* 0.5D0) - IDNINT( V( NN1)* SCL)
C ---------------------------------------------------------------------
IF ( J.EQ.NX) THEN
IX2 = X0 + IWD
IY2 = Y0 + IDNINT( DFLOAT( IWD)* 0.5D0)
END IF
C
500 CONTINUE
C --------------------------------------------------------------------
C
RETURN
END
C **********************************************************************
C
SUBROUTINE STND ( X, Y, LXY, NXY, NNP )
C
IMPLICIT REAL*8 (A-H,O-Z)
C
DIMENSION X( NNP), Y( NNP), LXY( NXY, 4)
C
COMMON /CMPD/ KSY( 200), NY, KSX( 200), NX
C
C ----- To set Node numbers for U-V Distribution on Z = 0.5 ---
C =========
C ==> LXY( NXY,4) <==> NODE( NEL,4) ---
C
C -----------------------------------------
C
C ----- Attention : All point information kept ---
C
NY = 0
C ------------------------------------------------
DO 100 I = 1, NXY
C
NI = LXY( I,1) ! or 2
C ----- Attention ---------------
C
IF ( DABS( X( NI) - 0.5D0).LE.0.01D0) THEN
NY = NY + 1 ! NY < = 200
KSY( NY) = NI
END IF
C
100 CONTINUE
C ------------------------------------------------
DO 200 I = 1, NXY
C
NI = LXY( I,2) ! or 2
C ----- Attention ---------------
C
IF ( DABS( X( NI) - 0.5D0).LE.0.01D0) THEN
NY = NY + 1 ! NY < = 200
KSY( NY) = NI
END IF
C
200 CONTINUE
C ------------------------------------------------
DO 300 I = 1, NXY
C
NI = LXY( I,3) ! or 2
C ----- Attention ------------
C
IF ( DABS( X( NI) - 0.5D0).LE.0.01D0) THEN
NY = NY + 1 ! NY < = 200
KSY( NY) = NI
END IF
C
300 CONTINUE
C ------------------------------------------------
DO 400 I = 1, NXY
C
NI = LXY( I,4) ! or 2
C ----- Attention -------------
C
IF ( DABS( X( NI) - 0.5D0).LE.0.01D0) THEN
NY = NY + 1 ! NY < = 200
KSY( NY) = NI
END IF
C
400 CONTINUE
C ------------------------------------------------
C
NX = 0
C
DO 500 I = 1, NXY
C
NI = LXY( I,1)
C ----- Attention ---
C
IF ( DABS( Y( NI) - 0.5D0).LE.0.01D0) THEN
NX = NX + 1 ! NX < = 200
KSX( NX) = NI
END IF
C
500 CONTINUE
C ------------------------------------------------
DO 600 I = 1, NXY
C
NI = LXY( I,2)
C ----- Attention ---
C
IF ( DABS( Y( NI) - 0.5D0).LE.0.01D0) THEN
NX = NX + 1 ! NX < = 200
KSX( NX) = NI
END IF
C
600 CONTINUE
C ------------------------------------------------
DO 700 I = 1, NXY
C
NI = LXY( I,3)
C ----- Attention ---
C
IF ( DABS( Y( NI) - 0.5D0).LE.0.01D0) THEN
NX = NX + 1 ! NX < = 200
KSX( NX) = NI
END IF
C
700 CONTINUE
C ------------------------------------------------
DO 800 I = 1, NXY
C
NI = LXY( I,4)
C ----- Attention ---
C
IF ( DABS( Y( NI) - 0.5D0).LE.0.01D0) THEN
NX = NX + 1 ! NX < = 200
KSX( NX) = NI
END IF
C
800 CONTINUE
C ------------------------------------------------
* WRITE(6,*) ' * NX =', NX,' NY =', NY
* WRITE(6,*) ' * KSY=', ( KSY(I), I=1, NY )
* WRITE(6,*) ' * KSX=', ( KSX(I), I=1, NX)
C
RETURN
END
C **********************************************************************
C
SUBROUTINE TTLE
C
IMPLICIT REAL*8 (A-H,O-Z)
C
WRITE(6,*) ' '
WRITE(6,*) 'C ************************************************* C'
WRITE(6,*) 'C C'
WRITE(6,*) 'C HF3D - FEM - R - FIG C'
WRITE(6,*) 'C C'
WRITE(6,*) 'C ( Version 3.8 ) C'
WRITE(6,*) 'C C'
WRITE(6,*) 'C ( Graphics output ) C'
WRITE(6,*) 'C C'
WRITE(6,*) 'C ( Rectangular Parallelepiped Elements ) C'
WRITE(6,*) 'C C'
WRITE(6,*) 'C Copyright 2011 : Yasuhiro MATSUDA C'
WRITE(6,*) 'C C'
WRITE(6,*) 'C ************************************************* C'
C
RETURN
END
C **********************************************************************
C
SUBROUTINE VIEW ( U, V, W, X, Y, Z, NNP )
C
IMPLICIT REAL*8 (A-H,O-Z)
C
C ----- Attention : XXX** ( Integer) ---
C
DIMENSION X( NNP), Y( NNP), Z( NNP), U( NNP), V( NNP), W( NNP)
C
COMMON /VWCT/ HR, VR, DIS, DS2, OX, OY, OZ, SCA, RR
C
INTEGER*2 ICR, IX3, IY3, IX4, IY4, IXS, IYS, IXE, IYE
C
ICR = 15
C
SCL = 0.05D0
IF ( NNP.GE.9000 ) SCL = 0.005D0
IF ( NNP.GE.29000) SCL = 0.0001D0 ! 31 x 31 x 31 = 29791
C
C -------- Attention ----------------------------------------------
C
ARSC = 0.01D0
C ----------------------------------------------------------
C
C ----- To Plot 3-D Geometry with x-y-z coordinates ---
C
CALL PGEO ( X, Y, Z, NNP )
C
C ----------------------------------------------------------
C
ICT = 150
S = DSIN( 0.2D0)
C = DCOS( 0.2D0)
EPS = 1.D-4
C
C ----- To Plot 3D-Velocity Perspectives -------------------
C
DO 100 I = 1, NNP
C
AL = DSQRT( U( I)** 2 + V( I)** 2 + W( I)** 2 )
C =================================================
C
C ------------------------------------------------
IF ( AL.GT.EPS ) THEN
C
X1 = X( I)
Y1 = Y( I)
Z1 = Z( I)
CALL VPNT ( X1, Y1, Z1, XS, YS )
C
IXS = XS
IYS = 399 - YS
X1 = X1 + U( I)* SCL
Y1 = Y1 + V( I)* SCL
Z1 = Z1 + W( I)* SCL
CALL VPNT ( X1, Y1, Z1, XE, YE )
C
IXE = XE
IYE = 399 - YE
IXS = IXS + ICT
IXE = IXE + ICT
CALL DRAW_LINE@( IXS, IYS, IXE, IYE, ICR )
C
AL = DSQRT( ( XE-XS )** 2 + ( YE-YS )** 2 )
C ============================================
C
C --------------------------------------
IF ( AL.GT.EPS ) THEN
C
T = ARSC / AL
IF ( T.GT.1.D0) T = 1.D0
X2 = XE + T* ( XS - XE )
Y2 = YE + T* ( YS - YE )
C
DX = X2 - XE
DY = Y2 - YE
X3 = XE + DX* C + DY* S
Y3 = YE - DX* S + DY* C
X4 = XE + DX* C - DY* S
Y4 = YE + DX* S + DY* C
IY3 = 399 - Y3
IY4 = 399 - Y4
C
IXE = XE + ICT
IYE = 399 - YE
IX4 = X4 + ICT
IX3 = X3 + ICT
C
CALL DRAW_LINE@( IX4, IY4, IXE, IYE, ICR )
CALL DRAW_LINE@( IXE, IYE, IX3, IY3, ICR )
C
END IF
C --------------------------------------
C
END IF
C ------------------------------------------------
C
100 CONTINUE
C ----------------------------------------------------------
C
RETURN
END
C **********************************************************************
C
SUBROUTINE VPNT ( X, Y, Z, VX, VY )
C
IMPLICIT REAL*8 (A-H,O-Z)
C
COMMON /VWCT/ HR, VR, DIS, DS2, OX, OY, OZ, SCA, RR
C
XX = X - OX
YY = Y - OY
ZZ = Z - OZ
C
PI = 3.141592653589793D0/ 180.D0
RH = PI* HR
RV = PI* VR
RL = PI* RR
C
X1 = XX* DCOS( RV) + YY* DSIN( RV)
Y2 = - XX* DSIN( RV) + YY* DCOS( RV)
X2 = X1* DCOS( RH) + ZZ* DSIN( RH)
C
Z1 = - X1* DSIN( RH) + ZZ* DCOS( RH)
X3 = X2* DCOS( RL) + Y2* DSIN( RL)
Y3 = - X2* DSIN( RL) + Y2* DCOS( RL)
D = DIS - Z1
C
VX = 450.D0* ( X3* DS2/ ( D + DS2)* SCA + 0.5D0) + 100.D0
VY = 450.D0* ( Y3* DS2/ ( D + DS2)* SCA + 0.5D0)
C
RETURN
END
C ************************************************************************
C
SUBROUTINE VW3D ( X, Y, Z, U, V, W, NNP )
C
IMPLICIT REAL*8 (A-H,O-Z)
C
COMMON /VWCT/ HR, VR, DIS, DS2, OX, OY, OZ, SCA, RR
C
DIMENSION X( NNP), Y( NNP), Z( NNP), U( NNP), V( NNP), W( NNP)
C
SCA = 0.5D0
DIS = 2.D0
DS2 = 6.D0
C
C ----- View Point ---
C
VR = 50.D0
HR = - 40.D0
RR = - 40.D0
C
C -------------------------------------------------
C
CALL VIEW ( U, V, W, X, Y, Z, NNP )
C
C -------------------------------------------------
C
RETURN
END
C **********************************************************************