SUBROUTINE TFINLR(XB,YB,XL,YL,XT,YT,XR,YR,NB,NT, * NUMI,NLR,X,Y) C***BEGIN PROLOGUE TFINLR C***DATE WRITTEN 000614 (YYMMDD) C***REVISION DATE 000614 (YYMMDD) C***CATEGORY NO. C***KEYWORDS TRANSFINITE BILINEAR INTERPOLATION, BLENDING FUNCTION C INTERPOLATION, ALGEBRAIC GRID GENERATION, BOUNDARY-FITTED C GRID GENERATION, NUMERICAL GRID GENERATION, CURVILINEAR C COORDINATE SYSTEM C C***AUTHOR SAUNDERS,B.V., (NIST) C***PURPOSE Generates a two dimensional boundary-fitted grid by C using transfinite bilinear interpolation. Grid points C on left and right boundaries are fixed to be the points C input by the user. C C***DESCRIPTION C C TFINLR generates a two dimensional boundary-fitted grid by C transfinite bilinear interpolation. In effect, TFINLR creates C a curvilinear coordinate system that maps the unit square onto C a region whose boundary is input by the user. The boundary must be C divided into four sides: bottom, left, top, right. Grid points C on the left and right boundaries are fixed to be the points C input by the user so that the left and right boundaries are C interpolated accurately. The number of points on the left and right C boundaries must be the same. C (Also see TFINTB, TFINBLIN,TFINXACT) C C Advantage: C Left and right boundary grid points interpolate user inputted C left and right boundary points exactly. C C Disadvantage: C 1. Number of grid points along left and right boundaries restricted. C C INPUT C C XB,YB REAL(NB), REAL(NB) C Real vectors containing abcissas, ordinates C for bottom boundary points. C C XL,YL REAL(NLR),REAL(NLR) C Real vectors containing abcissas, ordinates C for left boundary points. C C XT,YT REAL(NT),REAL(NT) C Real vectors containing abcissas, ordinates C for top boundary points. C C XR,YR REAL(NLR),REAL(NLR) C Real vectors containing abcissas, ordinates C for right boundary points. C C NB INTEGER C Number of bottom boundary points. NB must be C less than or equal to NMAXBP. C C NLR INTEGER C Number of left(right) boundary points. C NLR must be less than or equal C to NMAXBP. C C NT INTEGER C Number of top boundary points. NT must be C less than or equal to NMAXBP. C C NUMI INTEGER C NUMI x NLR equal the dimensions of the outputted C grid. C C C C OUTPUT C C X,Y REAL(NUMI,NLR),REAL(NUMI,NLR) C Matrices containing abcissas, ordinates of C grid points. First row of grid points would C be (X(1,1),Y(1,1)),....,(X(NUMI,1),Y(NUMI,1)). C Last row of grid points would be C (X(1,NLR),Y(1,NLR)),....,(X(NUMI,NLR), C Y(NUMI,NLR)) C C C PARAMETER C CONSTANT C C NMAXBP INTEGER C Maximum number of boundary points allowed C as input on any side of the domain. If size C insufficient for the needs of the user, user C can globally change parameter statement and C recompile. C C C EXAMPLE PROGRAM C C PROGRAM TESTLR CC CC Sample driver for TFINLR. Computes NUMI x NLR grid for nonconvex figure CC whose boundary data is in file input.data . CC C PARAMETER(NUMI=5,NLR=4) C REAL X(NUMI,NLR),Y(NUMI,NLR),XB(100),YB(100),XL(100), C * YL(100),XT(100),YT(100),XR(100),YR(100) C OPEN(10, FILE='input.data') C OPEN(17, FILE='output.data') CC CC Number of points on left and right curves must be CC the same value, NLR. Consequently, only NB and CC NT need to be in the input file, but placing CC values for all 4 sides allows the same input CC format to be used for other transfinite routines CC TFINBLIN AND TFINTB. Routine TFINXACT can also CC use the same format, but the first read ignores CC NT and NR. CC CC CC INPUT BOUNDARY POINTS CC C READ(10,*) NB,NL,NT,NR C READ(10,*) (XB(J),J=1,NB) C READ(10,*) (YB(J),J=1,NB) C READ(10,*) (XL(J),J=1,NL) C READ(10,*) (YL(J),J=1,NL) C READ(10,*) (XT(J),J=1,NT) C READ(10,*) (YT(J),J=1,NT) C READ(10,*) (XR(J),J=1,NR) C READ(10,*) (YR(J),J=1,NR) C C CHECK TO MAKE SURE NL = NR = NLR. IF NOT PRINT C ERROR MESSAGE AND STOP. C C IF(NL.NE.NLR.OR.NR.NE.NLR) THEN C PRINT*,'NUMBER OF POINTS ON LEFT, RIGHT BOUNDARIES MUST BE NLR' C STOP C ENDIF C CALL TFINLR(XB,YB,XL,YL,XT,YT,XR,YR,NB,NT,NUMI, C * NLR,X,Y) CC CC OUTPUT GRID POINTS CC C DO 500 J=1,NLR C DO 400 I=1,NUMI C WRITE(17,*) X(I,J),Y(I,J) C 400 CONTINUE C WRITE(17,*) C 500 CONTINUE C STOP C END C C C input.data C C 3 4 3 4 C -2 0 4 C -3 0 2 C -2 -3 -3.5 -4 C -3 0 2 6 C -4 0 6 C 6 6 6 C 4 4.5 5 6 C 2 4 5 6 C C C output.data C C -2.000000 -3.000000 C -0.8798263 -1.319740 C 0.3875484 0.1937742 C 2.193774 1.096887 C 4.000000 2.000000 C C -3.000000 0.0000000E+00 C -1.378217 1.286840 C 0.3416991 2.462516 C 2.420849 3.231258 C 4.500000 4.000000 C C -3.500000 2.000000 C -1.501609 2.893420 C 0.5458498 3.731258 C 2.772924 4.365629 C 5.000000 4.999999 C C -4.000000 6.000000 C -1.500000 6.000000 C 1.000000 6.000000 C 3.500000 6.000000 C 6.000000 6.000000 C***REFERENCES Knupp,P., Steinberg, S., Fundamentals of Grid C Generation, CRC Press, 1993. C Thompson, J.F., Warsi, Z.U.A., Mastin, C.W., C Numerical Grid Generation: Foundations and C Applications, Elsevier, 1985. C Gordon, W.J. and Hall, C.A., "Constructions of C Curvilinear Coordinate Systems and Applications C to Mesh Generation," Int. Journal for Numer. C Methods in Eng., Vol. 7, 461-477, 1963. C C C***END PROLOGUE TFINLR C PARAMETER(NMAXBP=500) COMMON/TOPBOT/S(2,NMAXBP),XTB(2,NMAXBP),YTB(2,NMAXBP) REAL XB(NB),YB(NB),XL(NLR),YL(NLR),XT(NT),YT(NT),XR(NLR), * YR(NLR),X(NUMI,NLR),Y(NUMI,NLR),A(NMAXBP),B(NMAXBP) C C FLEXIBLE BOTTOM, TOP BOUNDARIES C DO 5 J=1,NB XTB(1,J) = XB(J) YTB(1,J) = YB(J) 5 CONTINUE DO 10 J=1,NT XTB(2,J) = XT(J) YTB(2,J) = YT(J) 10 CONTINUE C C CONSTRUCT ARCLENGTH PARAMETRIZATIONS OF BOTTOM, TOP BOUNDARIES C DO 20 I=1,2 S(I,1)=0.0 20 CONTINUE DO 30 J=2,NB S(1,J)= SQRT((XB(J)-XB(J-1))**2 + * (YB(J)-YB(J-1))**2) +S(1,J-1) 30 CONTINUE DO 40 J=2,NT S(2,J)= SQRT((XT(J)-XT(J-1))**2 + * (YT(J)-YT(J-1))**2) +S(2,J-1) 40 CONTINUE DO 50 J=2,NB-1 S(1,J) = S(1,J)/S(1,NB) 50 CONTINUE DO 60 J=2,NT-1 S(2,J) = S(2,J)/S(2,NT) 60 CONTINUE S(1,NB)=1.0 S(2,NT)=1.0 C C COMPUTE VECTORS A AND B C A(1)=0.0 A(NUMI)=1.0 DO 160 I=2,NUMI-1 A(I)=FLOAT(I-1)/FLOAT(NUMI-1) 160 CONTINUE B(1)=0.0 B(NLR)=1.0 DO 170 J=2,NLR-1 B(J)=FLOAT(J-1)/FLOAT(NLR-1) 170 CONTINUE C C C******************TRANSFINITE MAPPING ******************** C DO 500 J=1,NLR ETA = B(J) DO 400 I=1,NUMI XI = A(I) X(I,J) = (1.0-XI)*XL(J)+XI*XR(J) * + (1.0-ETA)*XBOT(XI)+ETA*XTOP(XI) * -(1.0-XI)*(1.0-ETA)*XBOT(0.0)-(1.0-XI)*ETA*XTOP(0.0) * -XI*(1.0-ETA)*XBOT(1.0)-XI*ETA*XTOP(1.0) Y(I,J) = (1.0-XI)*YL(J)+XI*YR(J) * +(1.0-ETA)*YBOT(XI)+ETA*YTOP(XI) * -(1.0-XI)*(1.0-ETA)*YBOT(0.0)-(1.0-XI)*ETA*YTOP(0.0) * -XI*(1.0-ETA)*YBOT(1.0)-XI*ETA*YTOP(1.0) 400 CONTINUE 500 CONTINUE RETURN END REAL FUNCTION XBOT(U) PARAMETER(NMAXBP=500) COMMON/TOPBOT/S(2,NMAXBP),XTB(2,NMAXBP),YTB(2,NMAXBP) J=1 10 CONTINUE IF(U.GE.S(1,J).AND.U.LE.S(1,J+1)) THEN XBOT=(XTB(1,J)*(S(1,J+1)-U)+XTB(1,J+1)*(U-S(1,J)))/ * (S(1,J+1)-S(1,J)) RETURN ELSE J=J+1 ENDIF GOTO 10 END C REAL FUNCTION YBOT(U) PARAMETER(NMAXBP=500) COMMON/TOPBOT/S(2,NMAXBP),XTB(2,NMAXBP),YTB(2,NMAXBP) J=1 10 CONTINUE IF(U.GE.S(1,J).AND.U.LE.S(1,J+1)) THEN YBOT=(YTB(1,J)*(S(1,J+1)-U)+YTB(1,J+1)*(U-S(1,J)))/ * (S(1,J+1)-S(1,J)) RETURN ELSE J=J+1 ENDIF GOTO 10 END C C REAL FUNCTION XTOP(U) PARAMETER(NMAXBP=500) COMMON/TOPBOT/S(2,NMAXBP),XTB(2,NMAXBP),YTB(2,NMAXBP) J=1 10 CONTINUE IF(U.GE.S(2,J).AND.U.LE.S(2,J+1)) THEN XTOP=(XTB(2,J)*(S(2,J+1)-U)+XTB(2,J+1)*(U-S(2,J)))/ * (S(2,J+1)-S(2,J)) RETURN ELSE J=J+1 ENDIF GOTO 10 END C C REAL FUNCTION YTOP(U) PARAMETER(NMAXBP=500) COMMON/TOPBOT/S(2,NMAXBP),XTB(2,NMAXBP),YTB(2,NMAXBP) J=1 10 CONTINUE IF(U.GE.S(2,J).AND.U.LE.S(2,J+1)) THEN YTOP=(YTB(2,J)*(S(2,J+1)-U)+YTB(2,J+1)*(U-S(2,J)))/ * (S(2,J+1)-S(2,J)) RETURN ELSE J=J+1 ENDIF GOTO 10 END