pois.f90

FDS为火灾动力学模拟软件源代码,该软件为开源项目,代码语言主要为FORTRAN,可在WINDOWS和LINUX下编译运行,详细说明可参考http://fire.nist.gov/fds/官方网址

  DO  I=1,L
    XT(I,K,MS2+1) = 2._EB*X(I,MS2+1,K)
  END DO
END DO
404 CONTINUE
RETURN
END SUBROUTINE VCOST1


SUBROUTINE VCOSTA(M,LDIMX,MM1,PL,X,XT)

! +--------------------------------------------------------------------+
! |                                                                    |
! |                       COPYRIGHT (C) 1989 BY                        |
! |                          ROLAND A. SWEET                           |
! |                        ALL RIGHTS RESERVED                         |
! |                                                                    |
! +--------------------------------------------------------------------+

!     PACKAGE VFFTPAK, VERSION 1, JUNE 1989



INTEGER                       :: M
INTEGER                       :: LDIMX
INTEGER                   :: MM1
REAL(EB)   X(LDIMX,M), PL(LDIMX),XT(LDIMX,M)

DO  I=1,LDIMX
  X(I,1) = XT(I,1)
  X(I,2) = PL(I)
END DO
DO  J=4,M,2
  DO  I=1,LDIMX
    X(I,J) = X(I,J-2)-XT(I,J-1)
    X(I,J-1) = XT(I,J-2)
  END DO
END DO
IF (MOD(M,2) == 0) GO TO 409
DO  I=1,LDIMX
  X(I,M) = XT(I,MM1)
END DO
409 CONTINUE
RETURN
END SUBROUTINE VCOSTA


SUBROUTINE VCOSTI(N,C,WSAVE)

! +--------------------------------------------------------------------+
! |                                                                    |
! |                       COPYRIGHT (C) 1989 BY                        |
! |                          ROLAND A. SWEET                           |
! |                        ALL RIGHTS RESERVED                         |
! |                                                                    |
! +--------------------------------------------------------------------+

!     PACKAGE VFFTPAK, VERSION 1, JUNE 1989


!     *******   NOTE   ******

!     SPECIAL DOCUMENTATION FOR USE OF INTERNAL PARAMETER TPOSE.  IF
!     TPOSE = .TRUE., THEN ALL PREPROCESSING ROUTINES TRANSPOSE SECOND
!                     AND THIRD INDICES IN PREPARATION FOR TRANSFORM OF
!                     ORIGINAL SECOND INDEX.  IF
!     TPOSE = .FALSE., THEN ALL PREPROCESSING ROUTINES PREPARE ONLY FOR
!                      TRANSFORM OF THIRD INDEX.  TO USE IN THIS CASE,
!                      CALL ROUTINE WITH M AND N INTERCHANGED.



INTEGER                       :: N
REAL(EB)    C(N), WSAVE(N+15)

!                               INITIALIZE NOCOPY AND TPOSE TO DEFAULT
!                               VALUES

NOCOPY = .FALSE.
TPOSE = .TRUE.

PI = FSH19S()
IF (N <= 3) RETURN
NP1 = N+1
NS2 = N/2
DT = PI/(N-1)

DO  K=2,NS2
  KC = NP1-K
  C(K) = 2._EB*SIN((K-1)*DT)
  C(KC) = 2._EB*COS((K-1)*DT)
END DO
CALL VRFFTI (N-1,WSAVE)
RETURN
END SUBROUTINE VCOSTI


SUBROUTINE VRFFTF (M,N,R,MDIMR,RT,WSAVE)

! +--------------------------------------------------------------------+
! |                                                                    |
! |                       COPYRIGHT (C) 1989 BY                        |
! |                          ROLAND A. SWEET                           |
! |                        ALL RIGHTS RESERVED                         |
! |                                                                    |
! +--------------------------------------------------------------------+

!     PACKAGE VFFTPAK, VERSION 1, JUNE 1989



INTEGER                   :: M
INTEGER                   :: N
INTEGER                   :: MDIMR
REAL(EB)   RT(M,N)
REAL(EB)   WSAVE(N+15)
REAL(EB)         R(MDIMR,N)
IF (N == 1) RETURN
CALL VRFTF1 (M,N,R,MDIMR,RT,WSAVE(1),WSAVE(N+1))
RETURN
END SUBROUTINE VRFFTF


SUBROUTINE VRFFTI (N,WSAVE)

! +--------------------------------------------------------------------+
! |                                                                    |
! |                       COPYRIGHT (C) 1989 BY                        |
! |                          ROLAND A. SWEET                           |
! |                        ALL RIGHTS RESERVED                         |
! |                                                                    |
! +--------------------------------------------------------------------+

!     PACKAGE VFFTPAK, VERSION 1, JUNE 1989



INTEGER                   :: N
REAL(EB)         WSAVE(N+15)

!                               INITIALIZE NOCOPY AND TPOSE TO DEFAULT
!                               VALUES

NOCOPY = .FALSE.
TPOSE = .TRUE.

IF (N <= 1) RETURN
CALL VRFTI1 (N,WSAVE(1),WSAVE(N+1))
RETURN
END SUBROUTINE VRFFTI


SUBROUTINE VRFTF1 (M,N,C,MDIMC,CH,WA,FAC)

! +--------------------------------------------------------------------+
! |                                                                    |
! |                       COPYRIGHT (C) 1989 BY                        |
! |                          ROLAND A. SWEET                           |
! |                        ALL RIGHTS RESERVED                         |
! |                                                                    |
! +--------------------------------------------------------------------+

!     PACKAGE VFFTPAK, VERSION 1, JUNE 1989



INTEGER                       :: M
INTEGER                       :: N
INTEGER                   :: MDIMC
REAL(EB)    CH(M,N)
REAL(EB)   WA(N)
REAL(EB)  FAC(15)
REAL(EB)         C(MDIMC,N)

NF = FAC(2)
NA = 1
L2 = N
IW = N
DO  K1=1,NF
  KH = NF-K1
  IP = FAC(KH+3)
  L1 = L2/IP
  IDO = N/L2
  IDL1 = IDO*L1
  IW = IW-(IP-1)*IDO
  NA = 1-NA
  IF (IP /= 4) GO TO 102
  IX2 = IW+IDO
  IX3 = IX2+IDO
  IF (NA /= 0) GO TO 101
  CALL VRADF4 (M,IDO,L1,C,MDIMC,CH,M,WA(IW),WA(IX2),WA(IX3))
  GO TO 110
  101    CALL VRADF4 (M,IDO,L1,CH,M,C,MDIMC,WA(IW),WA(IX2),WA(IX3))
  GO TO 110
  102    IF (IP /= 2) GO TO 104
  IF (NA /= 0) GO TO 103
  CALL VRADF2 (M,IDO,L1,C,MDIMC,CH,M,WA(IW))
  GO TO 110
  103    CALL VRADF2 (M,IDO,L1,CH,M,C,MDIMC,WA(IW))
  GO TO 110
  104    IF (IP /= 3) GO TO 106
  IX2 = IW+IDO
  IF (NA /= 0) GO TO 105
  CALL VRADF3 (M,IDO,L1,C,MDIMC,CH,M,WA(IW),WA(IX2))
  GO TO 110
  105    CALL VRADF3 (M,IDO,L1,CH,M,C,MDIMC,WA(IW),WA(IX2))
  GO TO 110
  106    IF (IP /= 5) GO TO 108
  IX2 = IW+IDO
  IX3 = IX2+IDO
  IX4 = IX3+IDO
  IF (NA /= 0) GO TO 107
  CALL VRADF5(M,IDO,L1,C,MDIMC,CH,M,WA(IW),WA(IX2),WA(IX3),WA(IX4))
  GO TO 110
  107 CALL VRADF5 (M,IDO,L1,CH,M,C,MDIMC,WA(IW),WA(IX2),WA(IX3),WA(IX4))
  GO TO 110
  108    IF (IDO == 1) NA = 1-NA
  IF (NA /= 0) GO TO 109
  CALL VRADFG (M,IDO,IP,L1,IDL1,C,C,C,MDIMC,CH,CH,M,WA(IW))
  NA = 1
  GO TO 110
  109    CALL VRADFG (M,IDO,IP,L1,IDL1,CH,CH,CH,M,C,C,MDIMC,WA(IW))
  NA = 0
  110    L2 = L1
END DO
OUTARY=.TRUE.
IF (NOCOPY) THEN
  SCALE=SCALE*SQRT(1.0_EB/REAL(N,EB))
  IF (NA==0) THEN
    OUTARY=.FALSE.
  END IF
ELSE
  SCALE=SQRT(1.0_EB/REAL(N,EB))
  IF (NA == 1) GO TO 113
  DO  J=1,N
    DO  I=1,M
      C(I,J) = SCALE*CH(I,J)
    END DO
  END DO
  RETURN
  113    DO  J=1,N
    DO  I=1,M
      C(I,J)=SCALE*C(I,J)
    END DO
  END DO
END IF
RETURN
END SUBROUTINE VRFTF1


SUBROUTINE VRFTI1 (N,WA,FAC)

! +--------------------------------------------------------------------+
! |                                                                    |
! |                       COPYRIGHT (C) 1989 BY                        |
! |                          ROLAND A. SWEET                           |
! |                        ALL RIGHTS RESERVED                         |
! |                                                                    |
! +--------------------------------------------------------------------+

!     PACKAGE VFFTPAK, VERSION 1, JUNE 1989


INTEGER                       :: N
REAL(EB)   FAC(15)
REAL(EB)         WA(N)      , NTRYH(4)
DATA NTRYH(1),NTRYH(2),NTRYH(3),NTRYH(4)/4,2,3,5/

NL = N
NF = 0
J = 0
101 J = J+1
IF (J-4<0) GO TO 102
IF (J-4==0) GO TO 102
IF (J-4>0) GO TO 103
102 NTRY = NTRYH(J)
GO TO 104
103 NTRY = NTRY+2
104 NQ = NL/NTRY
NR = NL-NTRY*NQ
IF (NR<0) GO TO 101
IF (NR==0) GO TO 105
IF (NR>0) GO TO 101
105 NF = NF+1
FAC(NF+2) = NTRY
NL = NQ
IF (NTRY /= 2) GO TO 107
IF (NF == 1) GO TO 107
DO  I=2,NF
  IB = NF-I+2
  FAC(IB+2) = FAC(IB+1)
END DO
FAC(3) = 2
107 IF (NL /= 1) GO TO 104
FAC(1) = N
FAC(2) = NF
TPI = 2._EB*FSH19S()
ARGH = TPI/REAL(N,EB)
IS = 0
NFM1 = NF-1
L1 = 1
IF (NFM1 == 0) RETURN
DO  K1=1,NFM1
  IP = FAC(K1+2)
  LD = 0
  L2 = L1*IP
  IDO = N/L2
  IPM = IP-1
  DO  J=1,IPM
    LD = LD+L1
    I = IS
    ARGLD = REAL(LD,EB)*ARGH
    FI = 0._EB 
    DO  II=3,IDO,2
      I = I+2
      FI = FI+1._EB 
      ARG = FI*ARGLD
      WA(I-1) = COS(ARG)
      WA(I) = SIN(ARG)
    END DO
    IS = IS+IDO
  END DO
  L1 = L2
END DO
RETURN
END SUBROUTINE VRFTI1


SUBROUTINE VSCOSB(F,L,M,N,LDIMF,FT,C1,C2,WORK)

! +--------------------------------------------------------------------+
! |                                                                    |
! |                       COPYRIGHT (C) 1989 BY                        |
! |                          ROLAND A. SWEET                           |
! |                        ALL RIGHTS RESERVED                         |
! |                                                                    |
! +--------------------------------------------------------------------+

!     PACKAGE VFFTPAK, VERSION 1, JUNE 1989

INTEGER                   :: L
INTEGER                       :: M
INTEGER                       :: N
INTEGER                       :: LDIMF
REAL(EB)   FT(LDIMF*N,M)
REAL(EB)  C1(M)
REAL(EB)    C2(M)
REAL(EB)   WORK(M+15)
REAL(EB)   F(LDIMF*N,M)

!     TPOSE = .TRUE. IF TRANSFORMING SECOND INDEX (SO TRANSPOSE)
!           = .FALSE. IF TRANSFORMING THIRD INDEX


!     PREPROCESSING

IF (TPOSE) THEN
  CALL VSCSB1(L,M,N,LDIMF,F,FT,C1,C2)
ELSE
  DO  I=1,LDIMF*N
    FT(I,1)=.5_EB*F(I,1)
  END DO
DO  J=2,M
  DO  I=1,LDIMF*N
    FT(I,J) =(C1(J)*F(I,J)+C2(J)*F(I,M-J+2))
  END DO
END DO
END IF

!     REAL(EB),PERIODIC ANALYSIS

CALL VRFFTF(LDIMF*N,M,FT,LDIMF*N,F,WORK)

!     POSTPROCESSING

SCALE=SQRT(2.0_EB)*SCALE
IF (OUTARY) THEN
  CALL VSCSBA(M,N*LDIMF,FT,F)
ELSE
  CALL VSCSBA(M,N*LDIMF,F,FT)
END IF
RETURN
END SUBROUTINE VSCOSB


SUBROUTINE VSCOSF(F,L,M,N,LDIMF,FT,C1,C2,WORK)

! +--------------------------------------------------------------------+
! |                                                                    |
! |                       COPYRIGHT (C) 1989 BY                        |
! |                          ROLAND A. SWEET                           |
! |                        ALL RIGHTS RESERVED                         |
! |                                                                    |
! +--------------------------------------------------------------------+

!     PACKAGE VFFTPAK, VERSION 1, JUNE 1989


INTEGER                   :: L
INTEGER                       :: M
INTEGER                       :: N
INTEGER                       :: LDIMF
REAL(EB)    FT(LDIMF*N,M)
REAL(EB)  C1(M)
REAL(EB)   C2(M)
REAL(EB)   WORK(M+15)
REAL(EB)   F(LDIMF*N,M)

!     TPOSE = .TRUE. IF TRANSFORMING SECOND INDEX (SO TRANSPOSE)
!           = .FALSE. IF TRANSFORMING THIRD INDEX


!     PREPROCESSING

SCALE=SQRT(2.0_EB)*SCALE
IF (TPOSE) THEN
  CALL VSCSF1(L,M,N,LDIMF,F,FT)
ELSE
  DO  I=1,LDIMF*N
    FT(I,1)=F(I,1)
  END DO
CONTINUE
IF (MOD(M,2)==0) THEN
  DO  I=1,LDIMF*N
    FT(I,M)=F(I,M)
  END DO
END IF
DO  J=2,M-1,2
  DO  I=1,LDIMF*N
    FT(I,J)   =  0.5_EB*(F(I,J)+F(I,J+1))
    FT(I,J+1) = -0.5_EB*(F(I,J)-F(I,J+1))
  END DO
END DO
END IF

!     REAL(EB),PERIODIC SYNTHESIS

CALL VRFFTB(LDIMF*N,M,FT,LDIMF*N,F,WORK)

!     POSTPROCESSING

IF (OUTARY) THEN
  CALL VSCSFA(M,N*LDIMF,F,FT,C1,C2)
ELSE
  CALL VSCSFA(M,N*LDIMF,FT,F,C1,C2)
END IF
RETURN
END SUBROUTINE VSCOSF


SUBROUTINE VSCOSI(N,C1,C2,WSAVE)

! +--------------------------------------------------------------------+
! |                                                                    |
! |                       COPYRIGHT (C) 1989 BY                        |
! |                          ROLAND A. SWEET                           |
! |                        ALL RIGHTS RESERVED                         |
! |                                                                    |
! +--------------------------------------------------------------------+

!     PACKAGE VFFTPAK, VERSION 1, JUNE 1989

!      ENTRY VSSINI(N,C1,C2,WSAVE)


INTEGER                       :: N
REAL(EB)    WSAVE(N+15)
REAL(EB)   C1(N),C2(N)

PI=FSH19S()
DX=PI/(2*N)

!     GENERATE A(I)+-B(I)

DO  I=1,N
  C=COS((I-1)*DX)
  S=SIN((I-1)*DX)
  C1(I)=.5_EB*(S+C)
  C2(I)=.5_EB*(S-C)
END DO

!     INITIALIZE VRFFTPK ROUTINE

CALL VRFFTI(N,WSAVE)
RETURN
END SUBROUTINE VSCOSI


SUBROUTINE VSCOSQ(F,L,M,N,LDIMF,FT,C1,C2,C3,C4,WORK)

! +--------------------------------------------------------------------+
! |                                                                    |
! |                       COPYRIGHT (C) 1989 BY                        |
! |                          ROLAND A. SWEET                           |
! |                        ALL RIGHTS RESERVED                         |
! |                                                                    |
! +--------------------------------------------------------------------+

!     PACKAGE VFFTPAK, VERSION 1, JUNE 1989


INTEGER                   :: L
INTEGER                       :: M
INTEGER                       :: N
INTEGER                       :: LDIMF
REAL(EB)  FT(LDIMF*N,M)
REAL(EB)  C1(M)
REAL(EB)  C2(M)
REAL(EB)  C3(M)
REAL(EB)  C4(M)
REAL(EB)  WORK(M+15)
REAL(EB)   F(LDIMF*N,M)

!     TPOSE = .TRUE. IF TRANSFORMING SECOND INDEX (SO TRANSPOSE)
!           = .FALSE. IF TRANSFORMING THIRD INDEX


!     PREPROCESSING

IF (TPOSE) THEN
   CALL VSCSQ1(L,M,N,LDIMF,F,FT,C1,C2)
ELSE
   DO  I=1,LDIMF*N
      FT(I,1)=F(I,1)
   END DO
   IF (MOD(M,2)==0) THEN
      DO  I=1,LDIMF*N
         FT(I,M)=-F(I,M)
      END DO