最近用fortran编写最小二乘曲线模拟程序，一直出现问题
program Console1implicit none
integer::i,j,n,status,m,f,k
real*8 Xmin,Ymin,px,px1,average,yy,dt1,dt2,dt3,z,d1,p,c,d2,g,q,dt,err
real*8,allocatable::X(:),Y(:),d(:),h(:),a(:)
real s(20),t(20),b(20)
open(10,file="E:\Fortran\task\2017.8.3\Console1\Console1\11.dat",status="old")
k=0
do while(.true.)
read(10,*,iostat=status)
if(status/=0)exit
k=k+1
end do
n=k
allocate(X(n),Y(n))
rewind(10)
do j=1,n
read(10,*)X(j),Y(j)
end do
Xmin=minval(X)
Ymin=minval(Y)
close(10)
if(Ymin<=0.0)then
Y=log10(Y-Ymin+1.0)
else
Y=log10(Y)
end if
if(Xmin<=0.0)then
X=log10(X-Xmin+1.0)
else
X=log10(X)
end if
read (*,*)m
px=(X(1)+X(n))/2.0
X=X-px
f=m-1
px1=((real(n)+1.0)/sum(x**((2.0)*f)))**((1.0)/(2.0*f))
x=px1*x
d=y
loop1:do
call hpir1(x,y,a,n,m,dt1,dt2,dt3)
write(9,*)(i,a(i),i=1,m)
yy=0
average=sum(x)/real(n)
outter:do j=1,n
h(0)=0
inner:do i=1,m
h(j)=h(j-1)+a(i)*((x(j)-average)**(i-1))
end do inner
yy=max(yy,abs(y(j)-h(j)))
end do outter
if(y(j)>h(j)+yy/real(2.0))then
y(j)=h(j)+yy/real(2.0)
elseif(y(j)<h(j)-yy/real(2.0))then
y(j)=h(j)-yy/real(2.0)
end if
call RMS(i,y,d,n,err)
if(err<1.0e-6.or.j==999)then
exit loop1
else
j=j+1
end if
end do loop1
open(9,file="result01")
write(9,400)j,m,n,err,dt1,yy
400 format(3(i3.3,2x),es12.4,2x,es12.4,2x,f12.4)
end program
subroutine RMS(i,y,d,n,err)
implicit none
integer,intent(in)::n
real*8,intent(in)::d(i),y(i)
real*8 err
integer::i
err=0
do i=1,n
err=err+(d(i)-y(i))**2
end do
err=sqrt(err/real(n))
return
end
subroutine hpir1(x,y,a,n,m,dt1,dt2,dt3)
real*8 x(n),y(n),a(m),s(20),t(20),b(20)
double precision dt1,dt2,dt3,z,d1,p,c,d2,g,q,dt
do 5 i=1,m
5 a(i)=0.0
if (m>n) m=n
if (m>20) m=20
z=0.0
do 10 i=1,n
10 z=z+x(i)/n
b(1)=1.0
d1=n
p=0.0
c=0.0
do 20 i=1,n
p=p+(x(i)-z)
c=c+y(i)
20 continue
c=c/d1
p=p/d1
a(1)=c*b(1)
if (m>1) then
t(2)=1.0
t(1)=-p
d2=0.0
c=0.0
g=0.0
do 30 i=1,n
q=x(i)-z-p
d2=d2+q*q
c=y(i)*q+c
g=(x(i)-z)*q*q+g
30 continue
c=c/d2
p=g/d2
q=d2/d1
d1=d2
a(2)=c*t(2)
a(1)=c*t(1)+a(1)
endif
do 100 j=3,m
s(j)=t(j-1)
s(j-1)=-p*t(j-1)+t(j-2)
if (j>=4) then
do 40k=j-2,2,-1
40 s(k)=-p*t(k)+t(k-1)-q*b(k)
end if
s(1)=-p*T(1)-q*b(1)
d2=0.0
c=0.0
g=0.0
do 70 i=1,n
q=s(j)
do 60 k=j-1,1,-1
60q=q*(x(i)-z)+s(k)
d2=d2+q*q
c=y(i)*q+c
g=(x(i)-z)*q*q+g
70 continue
c=c/d2
p=g/d2
q=d2/d1
d1=d2
a(j)=c*s(j)
t(j)=s(j)
do 80 k=j-1,1,-1
a(k)=c*s(k)+a(k)
b(k)=t(k)
t(k)=s(k)
80 continue
100 continue
dt1=0.0
dt2=0.0
dt3=0.0
do 120 i=1,n
q=a(m)
do 110 k=m-1,1,-1
110 q=q*(x(i)-z)+a(k)
dt=q-y(i)
if (abs(dt)>dt3) dt3=abs(dt)
dt1=dt1+dt*dt
dt2=dt2+abs(dt)
120 continue
return
end
楼主是fortran新手，参考的徐士良的算法集，可是不知道为什么会报错 对了，再补充一下，我用的IVF，请大神帮帮忙，程序有点大，谢谢 Array a in Main Program has NOT allocated !!! fcode 发表于 2017-8-5 04:45
Array a in Main Program has NOT allocated !!!
可是程序第五行已经定义了它是可变数组啊 第五行，只是定义它是“可以”分配的，然而并没有分配。
应该像 15 行，对 X Y 分配那样，同样的对 a 进行分配。 fcode 发表于 2017-8-5 04:45
Array a in Main Program has NOT allocated !!!
可是我添加了allocate a(m)进去了，还是报错 燕雨蔷薇 发表于 2017-8-6 17:19
第五行，只是定义它是“可以”分配的，然而并没有分配。
应该像 15 行，对 X Y 分配那样，同样的对 a 进行 ...
请问添在什么地方呢，我试了好几个地方，包括添加在读取m后面也报错
错误 3 Compilation Aborted (code 1) E:\Fortran\task\2017.8.3\Console1\Console1\Console1.f90 1
错误 1 error #5082: Syntax error, found IDENTIFIER 'A' when expecting one of: ( E:\Fortran\task\2017.8.3\Console1\Console1\Console1.f90 55
错误 2 error #6724: An allocate/deallocate object must have the ALLOCATABLE or POINTER attribute. E:\Fortran\task\2017.8.3\Console1\Console1\Console1.f90 55
我直接把m赋值了，这样就不会出现问题了，可是出现了数组越界问题:-(
页:
[1]