grid.f90

This program incorporates the FV method for solving the Navier-Stokes equations using 2D, Cartesian

module runparam
    integer::NI,NJ
    real*8,dimension(:),allocatable::XALL,YALL
end module


!C###################################################################      
PROGRAM EXPA

implicit none
!C###################################################################
!C#                                                                 #
!C#    PROGRAM FOR THE 2D LID DRIVEN CAVETY PROBLEM                 #
!C#                                                                 #
!C#                                   Jesper Winther St鎟dahl, 2008 #
!C###################################################################
	CHARACTER FILIN*12,FILOUT*12
	LOGICAL LPLOT
!C
!C----------------------------------------
!C.....READ INPUT FILE NAME AND OPEN FILES
!C----------------------------------------
!C
	1 FORMAT(A12)
	FILIN = 'grid.inp'
	FILOUT = 'grid.out'
	Write(*,'("INPUT FILE NAME: ",A10)') FILIN
	Write(*,'("OUTPUT FILE NAME: ",A10)') FILOUT
	
	OPEN (UNIT=5,FILE=FILIN)
	OPEN (UNIT=3,FILE=FILOUT)
	
	REWIND 3 ! places the curser at the beginning of unit 3 (filout)
	

	READ(5,*) LPLOT
    write(*,'("Plot (T=true,F=false)",L10)') LPLOT
	
	CALL GRID
    if (LPLOT) then
		CALL GRIDPL
	end if
	
end program
!C
!C####################################################################
subroutine GRID
use runparam
implicit none
!C####################################################################
    integer i,j
    read(5,*) NI,NJ
    allocate(XALL(NI),YALL(NJ))
    
    WRITE(3,*) NI   ,'  Number of grid in x-direction '
	WRITE(3,*) NJ   ,'  Number of grid in y-direction'
	
    
    do i = 1,NI
        XALL(i) = (dble(i)-1d0)/(dble(NI)-1d0)
    end do
    write(3,*) XALL
    do j = 1,NJ
        YALL(j) = (dble(j)-1d0)/(dble(NJ)-1d0)
    end do
    write(3,*) YALL

	

end subroutine 


!C
!C#############################################################
subroutine GRIDPL
!C#############################################################
use runparam
integer J,I
real*8 XMIN,XMAX,YMIN,YMAX,XL,YL,SCFX,SCFY,SCF
integer IXALL(NI),IYALL(NJ)
	CHARACTER FILOUT*12
!C
!C-----------------------------------------------------------
!C         CREATE POSTSCRIPT FILE WITH PLOT OF GRID LINES. 
!C         ROUTINE 'PSHEAD' WRITES HEADER OF THE POSCRIPT FILE. 
!C         PLOT AREA IS 7 x 7 INCHES, ORIENTATION IS PORTRAIT;
!C         COORDINATES ARE SCALED TO ALLOW ACCURATE PRINTING
!C         ON PRINTERS WITH 1200 x 1200 DPI 
!C-----------------------------------------------------------

	XMIN=XALL(1)
	XMAX=XALL(NI)
	YMIN=YALL(1)
	YMAX=YALL(NJ)

	XL=XMAX-XMIN
	YL=YMAX-YMIN
	SCFX=8400.0/XL
	SCFY=8400.0/YL
	SCF=MIN(SCFX,SCFY)
	XMIN=50.
	YMIN=50.
	XMAX=XMIN+XL*SCF*0.06
	YMAX=YMIN+YL*SCF*0.06
!C
	WRITE(FILOUT,'(4Hgrid,I1,3H.ps)') 1
	OPEN(UNIT=7,FILE=FILOUT)
	REWIND 7
	CALL PSHEAD(XMIN,XMAX,YMIN,YMAX)
!!!C
	DO I=1,NI
		IXALL(I)=INT(XALL(I)*SCF)
	END DO
	DO J=1,NJ
		IYALL(J)=INT(YALL(J)*SCF)
	END DO
!!!C---------------------------------------------------------
!!!C         PLOT BOUNDARY OF SOLUTION DOMAIN WITH THICK LINE
!!!C---------------------------------------------------------
	WRITE(7,*) '50 w'
	WRITE(7,71) IXALL(1),IYALL(1),' m ',IXALL(NI),IYALL(1),' l '
	WRITE(7,71) IXALL(NI),IYALL(NJ),' l ',IXALL(1),IYALL(NJ),' l '
	WRITE(7,72) IXALL(1),IYALL(1),' l  s '
	71 FORMAT(2I5,A3,2I5,A3)
	72 FORMAT(2I5,A6)
	73 FORMAT(2I5,A3,2I5,A6)
!!!C---------------------------------------------------------
!!!C         PLOT INNER GRID LINES WITH THIN LINES
!!!C---------------------------------------------------------
	WRITE(7,*) '1 w'
	DO J=2,NJ-1
		WRITE(7,73) IXALL(1),IYALL(J),' m ',IXALL(NI),IYALL(J),' l  s '
	END DO
	DO I=2,NI-1
		WRITE(7,73) IXALL(I),IYALL(1),' m ',IXALL(I),IYALL(NJ),' l  s '
	END DO
!   73 FORMAT(2I5,A3,2I5,A6)
	WRITE(7,*) 'p '
!!C
	CLOSE(7)
!      RETURN
END subroutine
!C
!C#######################################################
SUBROUTINE PSHEAD(XMIN,XMAX,YMIN,YMAX)
implicit none
!C#######################################################
!C
real*8 XMIN,XMAX,YMIN,YMAX

	WRITE(7,*) '%!PS-Adobe-2.0'
	WRITE(7,*) '%%Creator: GRIDPL'
	WRITE(7,*) '%%BoundingBox: ',XMIN,YMIN,XMAX,YMAX
	WRITE(7,*) '%%EndComments'
	WRITE(7,*) '/c {currentpoint} def /f {fill} def '
	WRITE(7,*) '/gr {grestore} def /gs {gsave} def /l {lineto} def '
	WRITE(7,*) '/m {moveto} def /n {newpath} def /p {showpage} def '
	WRITE(7,*) '/s {stroke} def /sg {setgray} def '
	WRITE(7,*) '/w {setlinewidth} def '
	WRITE(7,*) '50 50 translate 0.06 0.06 scale '
	WRITE(7,*) '1 setlinecap 1 setlinejoin '
	
END subroutine