vs报错0x00007FF63AE76A66 处(位于 ceshi1.exe 中)引发的异常: 0xC0000...

观星游月

vs 报错：0x00007FF63AE76A66 处(位于 ceshi1.exe 中)引发的异常: 0xC

[Fortran] 纯文本查看 复制代码

program gongcheng
    implicit none
    real(kind=8),parameter::G=1.3218278962236387e-42
    real(kind=8),parameter::m0=1.1158032783756145e60
    real(kind=8),parameter::R=20.0e18
    real(kind=8),parameter::pi=3.1415926535
    real(kind=8),parameter::h0=1.0e15
    real(kind=8)::p0(100)
    common p0
    real,external::f
    real,external::q
    integer::i
    real(kind=8)::rho(991),ene(991),pre(991)    
	real(kind=8)::x(20000),z(20000),y(20000),p(20000),m(20000),e(20000)
    double precision h,n,k1,k2,k3,k4,M1,M2,M3,M4
    
    do i=1,100
        p0(i)=i
    call rungekutta(x,y,z,f,q)    
    write(*,*)x(i+1),z(i+1)/m0
    end do
    end program
    
    function f(r,p,m)
    implicit none
    real(kind=8)::r(20000.d0),p(20000.d0),m(20000.d0),e(20000.d0)
    real(kind=8)::f(20000)
    real(kind=8)::rho(990.d0),ene(990.d0),pre(990.d0)
    integer::i
    real(kind=8),parameter::G=1.3218278962236387e-42
    real(kind=8),parameter::pi=3.1415926535
    open(unit=26,file="EOS.txt")
    do i=1.d0,990.d0
       read(26,*)rho(i),ene(i),pre(i)
    end do
    call spline(pre,ene,990,p,20000,e)
    f=-(e+p)*(G*m+4*pi*G*(r**3)*p)/(r*(r-2*G*m))
    return
    end
    
    function q(r,p,m)
    implicit none
    real(kind=8)::e(20000.d0),p(20000.d0),r(20000.d0),m(20000.d0)
    integer::i
    real(kind=8)::rho(990.d0),ene(990.d0),pre(990.d0)
    real::q(20000)
    real(kind=8),parameter::pi=3.1415926535 
    open(unit=26,file="EOS.txt")
    do i=1.d0,990.d0
       read(26,*)rho(i),ene(i),pre(i)
    end do
    call spline(pre,ene,990,p,20000,e)
    q=4*pi*e*r**2 
    return
    end
    
    
    subroutine rungekutta(x,y,z,f,q)
    real(kind=8)::x(20000),z(20000),y(20000),p(20000),m(20000)
    integer::i
    real,external::f
    real,external::q
    double precision h,n,k1,k2,k3,k4,M1,M2,M3,M4
    n=20.0e18/1.0e15
    h=1.0e15
    x(1)=0.0
    y(1)=p0
    z(1)=0.0
    do i=1,n-1
        y(i)=p(i)
        z(i)=m(i)
        x(i+1)=(i+1)*h
        k1=f(x(i),y(i))
        M1=q(x(i),y(i))
        k2=f(x(i)+0.5*h,y(i)+0.5*h*k1)
        M2=q(x(i)+0.5*h,y(i)+0.5*h*M1)
        k3=f(x(i)+0.5*h,y(i)+0.5*h*k2)
        M3=q(x(i)+0.5*h,y(i)+0.5*h*M2)
        k4=f(x(i)+h,y(i)+h*k3)
        M4=q(x(i)+h,y(i)+h*M3)
        y(i+1)=y(i)+h*(k1+k2+k3+k4)/6
        z(i+1)=z(i)+h*(M1+M2+M3+M4)/6 
        
        if (y(i+1)<=1.0*10**(-5))then
        exit
        end if
    end do
    end subroutine rungekutta

    
    
    subroutine spline(x,y,n,sx,m,f) !x,y:差值基点坐标，n：节点个数，sx：需要返回差值数值的x坐标，m：sx的容量，f：差值结果返回值                   
    implicit none
    integer:: i, j, k, m, n, n1, n2
    Real(kind=8)::x(n),y(n),sx(m),f(m)
    Real(kind=8):: s2(n), h(n), dy(n), s(n), e(n)     !直接都给数组定义N的空间大小，方便调用 
    Real(kind=8):: z, h1, h2, h3, h4
  
     n1=n-1!间距个数
     n2=n-2!内点个数
     
     do i = 1,n1
       h(i)  = x(i+1) - x(i)!计算点间距 
     dy(i) = (y(i+1) - y(i) ) / H(i)!计算斜率sj
    end do
    
     S2(1) = 0.d0; S2(N) = 0.d0!自然边界条件，两侧二阶导数为零
     
     DO i = 2, N1!求解内结点处的二阶导数
        S2(i) = 6.d0 * ( DY(i) - DY(i-1) )
     ENDDO
    
     Z = 0.5d0 / ( H(1) + H(2) )
     S(1) = -H(2) * Z
     E(1) = S2(2) * Z
     DO i = 2, N2
        k    = i - 1
        j    = i + 1
        Z    = 1.d0 / ( 2.d0*( H(i)+H(j) ) + H(i)*S(k) )
        S(i) = -H(j) * Z
        E(i) = ( S2(j)-H(i)*E(k) ) * Z
     ENDDO
    
     S2(N1) = E(N2)
     DO i = N2, 2, -1
        k     = i - 1
        S2(i) = S(k)*S2(i+1) + E(k)
     ENDDO
    
     DO i = 1, N1
        S(i) = ( S2(i+1) - S2(i) ) / H(i)
     ENDDO
    
     !具体求结果
     i = 2
     k = 1  
     DO j = 1, M
        DO 
           IF(sx(j)>x(i)) THEN     
                k = i
                i = i + 1
            ELSE
                EXIT
            ENDIF
        ENDDO
        H1    = SX(i) - X(k)
        H2    = SX (i)- X(i)
        H3    = H1 * H2
        H4    = S2(k) + H1*S(k)
        Z     = ( S2(i) + S2(k) + H4 ) / 6.d0
        F(j)  = Y(k) + H1*DY(k) + H3*Z
     ENDDO    
    end subroutine spline

000005: 写入位置 0x0000000000000000 时发生访问冲突。

fcode

你的问题太多了，我都不知道该怎么给你说起。

比如：
1. 数组定义时不要用浮点数做为数组的大小。
2. 函数的接口一定要明确，例如 function f(r,p,m) 是三个参数，但 rungekutta 函数里调用却写成了 k1=f(x(i),y(i)) 只有2个参数。并且 f q 函数的定义是返回大小为 20000 的数组，但 rungekutta 里面却认为返回的是一个单变量。
real,external::f 这样的写法并不会让编译器帮助你检查函数接口，所以 interface 是更好的选择。
3. 一定要写 implicit none，例如 rungekutta 函数里没有写，导致 p0 没有定义，也没有值。（昨天给你提出来，但你没有改）
昨天提的死循环的问题也没有改。
4. 大量相同大小的常数，例如 20000，990 可以写成 integer , parameter :: N = 2000 , S = 990 这样后期比较容易修改。

我给你的建议是，先从简单的代码写起。能得心应手的书写一二十行的代码之后，再来写几百行的代码。
每个程序单元，先做单元测试。然后再组织成稍大的程序。

例如，你先确定从 EOS.txt 中读取数组，单步调试，确保读取的数组是正确的。
然后，再调试插值函数，确保 spline 插值后的结果是符合自己预期的。（实际上，你的死循环就出现在插值函数里）
一般来说，从 EOS中读取数组，插值。是先做的，而你的代码逻辑是，每次对 f 和 q 求解，都去读 EOS文件、插值。显然是做了很多无用功。

最后，确保读取数据和插值正确了，再调试龙格库塔。

注意认真分析，每个函数应该设计几个参数？输入应该是什么，输出应该是什么？是单变量还是数组？数组的话，大小应该是多大？

观星游月

真的非常感谢，我打算从头开始再做一遍