C x=t^2 と置換すると
C  ∞ sin x      ∞ sin t^2
C S   -----dx = S   -------- 2t dt
C  0  x^2+1      0  t^{4}+1
C
C 部分積分で
C        ∞   x^2-1     
C 1 - 2 S   ---------sin xdx
C        0  (x^2+1)^3  
C
C x=e^t と置換すると
C  ∞ sin x      ∞  sin e^t
C S   -----dx = S    -------- e^t dt
C  0  x^2+1      -∞ e^{2t}+1
C
C x=e^{e^t}-1 と置換すると
C  ∞ sin x      ∞  sin (e^{e^t}-1)
C S   -----dx = S    --------------- e^{e^t} e^t dt
C  0  x^2+1      -∞ (e^{e^t}-1)^2+1
C
      PROGRAM OSCINT
      IMPLICIT REAL*8 (A-H,O-Z)
 1000 WRITE(*,*)'Menu : '
      WRITE(*,*)'1...Direct calculation by Trapezoidal rule '
     +           //'on [0,10000]'
      WRITE(*,*)'2...Change of variable x=t**2 and Trapezoidal on '
     +           //'[0,10000]'
      WRITE(*,*)'3...Twice integration by parts and Trapezoidal on '
     +           //'[0,1000]'
      WRITE(*,*)'4...Change of variable x = e**t and Trapezoidal on '
     +           //'[-20,20]'
      WRITE(*,*)'5...Finish'
      READ(*,*) K
      GO TO (1,2,3,4,5),K  ! 計算型 GO TO 文 (C の switch 文)
      GO TO 1000
C 直接計算
    1 A=0.0D0
      B=10000.0D0
      N=100000
      DO 20 I=1,4
          H=(B-A)/N
          S=F(A)/2*H
          X=A
          DO 10 J=1,N-1
              X=X+H
              S=S+F(X)*H
   10     CONTINUE     
          S=S+F(B)/2*H
          WRITE(*,300)H,' & ',S,' \\\\'
          N=N*10
   20 CONTINUE
      GO TO 1000
C      N=10
C      A=-5.0D0
C      B=5.0D0
C x=t^2 と置換  (0 < t < ∞)
    2 A=0.0D0
      B=10000
      N=100000
      DO 250 I=1,4
          H=(B-A)/N
          S=F2(A)/2*H
          X=A
          DO 150 J=1,N-1
              X=X+H
              S=S+F2(X)*H
  150     CONTINUE     
          S=S+F2(B)/2*H
          WRITE(*,300)H,' & ',S,' \\\\'
          N=N*10
  250 CONTINUE
      GO TO 1000
C 2 回部分積分し，分母の次数を上げる
    3 A=0.0D0
      B=1000.0D0
      N=10000
      DO 270 I=1,5
          H=(B-A)/N
          S=F5(A)/2*H
          X=A
          DO 170 J=1,N-1
              X=X+H
              S=S+F5(X)*H
  170     CONTINUE     
          S=S+F5(B)/2*H
          WRITE(*,300)H,' & ',1+S+S,' \\\\'
          N=N*10
  270 CONTINUE
      GO TO 1000
C x = e^t と置換 (-∞ < t < ∞)
    4 A=-20.0D0
      B=20.0D0
      N=400
      DO 200 I=1,5
          H=(B-A)/N
          S=F1(A)/2*H
          X=A
          DO 100 J=1,N-1
              X=X+H
              S=S+F1(X)*H
  100     CONTINUE     
          S=S+F1(B)/2*H
          WRITE(*,300)H,' & ',S,' \\\\'
          N=N*10
  200 CONTINUE
      GO TO 1000
  300 FORMAT(1H ,F18.15,A3,F18.15,A3)
  400 FORMAT(1H ,A13,F18.15)
    5 STOP
      END
C
      REAL*8 FUNCTION F(X)
      IMPLICIT REAL*8 (A-H,O-Z)
      F=DSIN(X)/(X*X+1)
      RETURN
      END
C
      REAL*8 FUNCTION F1(X)
      IMPLICIT REAL*8 (A-H,O-Z)
      Y=DEXP(X)
      F1=Y*DSIN(Y)/(Y*Y+1)
      RETURN
      END
C
      REAL*8 FUNCTION F2(X)
      IMPLICIT REAL*8 (A-H,O-Z)
C      Y=DEXP(X)
C      Z=DEXP(Y)
C      ZI=1/Z
C      F2=Y*DSIN(Z-1)/((1-ZI)*(Z-1)+ZI)
      Y=X*X
      F2=DSIN(Y)/(Y*Y+1)*X*2
      RETURN
      END
C
      REAL*8 FUNCTION F5(X)
      IMPLICIT REAL*8 (A-H,O-Z)
      Y=X*X
      F5=DSIN(X)*(1-Y*3)/(Y+1)**3
      RETURN
      END