积分无法计算
本帖最后由 hcj 于 2018-11-9 09:09 编辑PROGRAM MAIN
EXTERNAL D
DOUBLE PRECISION D,e,D1
DO e=0.1,1.0,0.01
D1=D(e)
WRITE(*,*)D1
END DO
END
FUNCTION D(e)
EXTERNAL F5
DOUBLE PRECISION D,F5,G5,e
CALL FLAGS2(F5,G5,e)
C WRITE(*,*)G5
D=G5
END
FUNCTION F5(b,e)
EXTERNAL IS,IC
DOUBLE PRECISION F5,IS,IC,b,e
F5=(1/3.1415926)*(COS(b*(e-IC(b))))*EXP(-1*b*IS(b))
RETURN
END
FUNCTION IS(b)
EXTERNAL F4
DOUBLE PRECISION IS,F4,G4,b
CALL FLAGS2(F4,G4,b)
C WRITE(*,*)G4
IS=G4
END
FUNCTION F4(V1,b)
EXTERNAL N
DOUBLE PRECISION F4,N,V1,b
F4=(SIN(b*V1))*N(V1)
RETURN
END
FUNCTION IC(b)
EXTERNAL F3
DOUBLE PRECISION IC,F3,G3,b
CALL FLAGS2(F3,G3,b)
C WRITE(*,*)G3
IC=G3
END
FUNCTION F3(V1,b)
EXTERNAL N
DOUBLE PRECISION F3,b,V1
F3=(COS(b*V1))*N(V1)
RETURN
END
FUNCTION N(V1)
EXTERNAL F1,F2
DOUBLE PRECISION N,V1,F1,F2,V
N=F1(V)-F2(V1)
RETURN
END
FUNCTION F2(V1)
EXTERNAL F
DOUBLE PRECISION F2,V1,F,A,B,G1
A=0.0
B=V1
EPS=0.005
CALL FLAGS1(A,B,F,EPS,G1)
C WRITE(*,*)G1
F2=G1
END
FUNCTION F1(V)
EXTERNAL F
DOUBLE PRECISION F1,F,G,V
CALL FLAGS(F,G)
C WRITE(*,*)G
F1=G
END
FUNCTION F(V)
DOUBLE PRECISION F,V
F=(3/(2*3.1415926*V))*LOG((COSH(SQRT(1/(2*V)))**2-
* COS(SQRT(1/(2*V)))**2)/(2*V))
RETURN
END
SUBROUTINE FLAGS(F,G)
DIMENSION T(5),C(5)
DOUBLE PRECISION F,G,T,C,X
DATA C/0.6790941054,1.638487956,2.769426772,
* 4.31594400,7.104896230/
DATA T/0.26355990,1.41340290,3.59642600,
* 7.08580990,12.64080000/
G=0.0D0
DO 10 I=1,5
X=T(I)
G=G+F(X)*C(I)
10 CONTINUE
END
SUBROUTINE FLAGS1(A,B,F,EPS,G)
DIMENSION T(5),C(5)
DOUBLE PRECISION A,B,F,G,T,C,S,P,H,AA,BB,W,X,Q
DATA T/-0.061798459,-0.5384693101,0.0,
* 0.5384693101,0.9061798459/
DATA C/0.2369268851,0.4786286705,0.5688888889,
* 0.4786286705,0.2369268851/
M=1
S=(B-A)*0.001
P=0.0
10 H=(B-A)/M
G=0.0
DO 30 I=1,M
AA=A+(I-1)*H
BB=A+I*H
W=0.0
DO 20 J=1,5
X=((BB-AA)*T(J)+(BB+AA))/2.0
W=W+F(X)*C(J)
20 CONTINUE
G=G+W
30 CONTINUE
G=G*H/2.0
Q=ABS(G-P)/(1.0+ABS(G))
IF ((Q.GE.EPS).AND.(ABS(H).GT.ABS(S)))THEN
P=G
M=M+1
GOTO 10
END IF
RETURN
END
SUBROUTINE FLAGS2(F,G,b)
DIMENSION T(5),C(5)
DOUBLE PRECISION F,G,T,C,X,b
DATA C/0.6790941054,1.638487956,2.769426772,
* 4.31594400,7.104896230/
DATA T/0.26355990,1.41340290,3.59642600,
* 7.08580990,12.64080000/
G=0.0D0
DO 10 I=1,5
X=T(I)
G=G+F(X,b)*C(I)
10 CONTINUE
END
带着未知量积分，把未知量一直顺延下去，在最后的一个积分里积掉，剩余一个未知量。然后做一个do循环。
这个方法为什么无法计算。与已知数据无法拟合
1. 我建议你换一本教科书。
都什么年代了，还用这么古老的代码风格，真是痛心。
放弃固定格式吧，放弃全部大写吧，放弃DATA语句吧，放弃 DO 数字 CONTINUE 吧
老代码尚且可以忍受，自己书写代码，为什么还要用这种古老的风格呢？？
2. 一定要书写 Implicit None，否则你出错了自己都不知道。比如 F3 函数里，你忘了定义 N 是 double，编译器误以为是 integer，结果就出现不可预料的后果。
3. 你的问题结果非常大，double 已经不足以承受。
特别是 F5 函数
F5=(1/3.1415926)*(COS(b*(e-IC(b))))*EXP(-1*b*IS(b))
当 b 到 7.8 以后，IS 返回就是 -190 多。EXP 已经 10的500次方了。 楚香饭 发表于 2018-11-9 09:53
1. 我建议你换一本教科书。
都什么年代了，还用这么古老的代码风格，真是痛心。
放弃固定格式吧，放弃全部 ...
谢谢。。。。。。。。。。。。
页:
[1]