调试沉积物动力学的旧 Fortran 代码

Question

我正在查看旧扫描论文中的一些 Fortran 代码。扫描质量不太好，所以我可能复印错了。我尝试使用在线 Fortran 编译器运行它，但它失败了。由于不熟悉 Fortran，我想知道是否有人可以指出语法没有意义的地方？该代码来自一篇关于沉积物动力学的论文：

Komar, PD 和 Miller, MC，1975。关于波浪和单向流下沉积物运动阈值的比较以及阈值实际评估的讨论：答复。沉积研究杂志，45（1）。

PROGRAM TSHOLD 
REAL LI, LO
G = 981.0 
PIE = 3.1416 
RHOW  = 1.00 
READ (6O,1) DIAM, RHOS
1 FORMAT (2X, F6.3,2X, F5.3) 
IF(DIAM .LT. 0.05) GO TO 5 
A = 0.463 * PIE 
B = 0.25 
GO TO 7 
5 A = 0.21 
B = 0.50 
7 PWR = 1.0 / (2.0 - B) 
FAC = (A * (RHOS - RHOW) * G/(RHOW * PIE**B))**PWR
FAC1 = FAC * DIAM**((1.0 - B) * PWR) 
T = 1.0 
15 J = 1.20 
LD = 156.13 * (T**2) 
UM = FAC1 * T**(B*PWR) 
WRITE(61,9) DIAM, T, UM
9 FORMAT(1H0, 10X, 17HGRAIN DIAMETER = ,F6.3,1X,2HCM //
1 11X, 14HWAVE PERIOD = ,F5.2, 1X, 3HSEC // 
2 11X, 22HORBITAL VELOCITY, UM = ,F6.2, 1X, 6HCM/SECl //
3 20X, 6HHEIGHT, 5X, 5HDEPTH, 8X, 3HH/L, 6X, 7HH/DEPTH //
4 22X, 2HCM, 8X, 2HCM /)
C INCREMENT WAVE HEIGHT, CALCULATE DEPTH 
H = 10.0 
DO 12 K = 1.60
SING = PIE * H / (UM * T) 
X = SING 
IF(X.LT.1.0) GO TO 30 
30 ASINH = X - 0.16666*X**3.0 + 0.07500* X ** 5.0 - 0.04464 * X ** 7.0
1 + 0.03038 * X ** 9.0 - 0.02237 * X ** 11.0 
32 LI = LD * (SINH(ASINH)/COSH(ASINH)) 
OPTH = ASINH * LI / 6.2832 
C CHECK WAVE STABILITY 
RATIO = H / DPTH 
IF(RATIO.GE.0.78) GO TO 11 
STEEP = H / LI 
TEST = 0.142 * (SINH(ASINH)/COSH(ASINH)) 
IF(STEEP.GE.TEST) GO TO 11 
WRITE(61,10) H, OPTH, STEEP, RATIO 
I0 FORMAT(IH0, 20X, F5.1, 4X, E9.3, 4X, F5.3, 4X, F4.2) 
11 H = H + 10.0 
12 CONTINUE 
T = T + 1.0 
15 CONTINUE 
END

Answer 1

问题更可能是旧的 Fortran 需要固定形式的代码格式，其中语句之前的空格数量非常重要。

以下是一些一般规则

• 正常语句从第 7 列及以后开始
• 行不能超过 72 列
• 第 6 列上的任何字符都表示该行是上一行的延续。我在上面的代码中看到以下几行9 FORMAT(..
• 位于第 1-5 列之间的数字表示标签，它可以是GO TO语句、DO语句或格式规范的目标。
• 第一列上的字符C（有时第一列上的任何字符）指示该行是注释行。

有关更多信息，请参阅https://people.cs.vt.edu/~asandu/Courses/MTU/CS2911/fortran_notes/node4.html。

根据上述规则，以下是如何使用正确的间距输入代码。我通过转换器运行F77代码，使其同时兼容F90和F77。下面的代码现在可以使用在线编译器进行编译。

      PROGRAM TSHOLD 
      REAL LI, LO 
      G = 981.0 
      PIE = 3.1416 
      RHOW  = 1.00 
      READ (60,1) DIAM, RHOS 
    1 FORMAT (2X, F6.3,2X, F5.3) 
      IF(DIAM .LT. 0.05) GO TO 5 
      A = 0.463 * PIE 
      B = 0.25 
      GO TO 7 
    5 A = 0.21 
      B = 0.50 
    7 PWR = 1.0 / (2.0 - B) 
      FAC = (A * (RHOS - RHOW) * G/(RHOW * PIE**B))**PWR 
      FAC1 = FAC * DIAM**((1.0 - B) * PWR) 
      T = 1.0 
     DO 15 J=1,20 
      LD = 156.13 * (T**2) 
      UM = FAC1 * T**(B*PWR) 
      WRITE(61,9) DIAM, T, UM 
    9 FORMAT(1H0, 10X, 17HGRAIN DIAMETER = ,F6.3,1X,2HCM //             &
     & 11X, 14HWAVE PERIOD = ,F5.2, 1X, 3HSEC //                        &
     & 11X, 22HORBITAL VELOCITY, UM = ,F6.2, 1X, 6HCM/SECl //           &
     & 20X, 6HHEIGHT, 5X, 5HDEPTH, 8X, 3HH/L, 6X, 7HH/DEPTH //          &
     & 22X, 2HCM, 8X, 2HCM /)                                           
! INCREMENT WAVE HEIGHT, CALCULATE DEPTH                                
      H = 10.0 
      DO 12 K = 1,60 
      SING = PIE * H / (UM * T) 
      X = SING 
      IF(X.LT.1.0) GO TO 30 
   30 ASINH = X - 0.16666*X**3.0 + 0.07500* X ** 5.0 - 0.04464 * X ** 7.&
     & + 0.03038 * X ** 9.0 - 0.02237 * X ** 11.0                       
   32 LI = LD * (SINH(ASINH)/COSH(ASINH)) 
      OPTH = ASINH * LI / 6.2832 
! CHECK WAVE STABILITY                                                  
      RATIO = H / DPTH 
      IF(RATIO.GE.0.78) GO TO 11 
      STEEP = H / LI 
      TEST = 0.142 * (SINH(ASINH)/COSH(ASINH)) 
      IF(STEEP.GE.TEST) GO TO 11 
      WRITE(61,10) H, OPTH, STEEP, RATIO 
   10 FORMAT(G14.4, 20X, F5.1, 4X, E9.3, 4X, F5.3, 4X, F4.2) 
   11 H = H + 10.0 
   12 CONTINUE 
      T = T + 1.0 
   15 CONTINUE 
      END                                           

我发现了几个转录错误，用点替换逗号，用字母 O 替换零，以及缺少DO语句。