Fortran Coder

标题: 求助大佬，动态数组如何作子程序形参 [打印本页]

---

作者: Proteus    时间: 2022-7-30 10:21
标题: 求助大佬，动态数组如何作子程序形参


遇到一个问题，我在主程序中定义了动态数组a，它要传递到子程序中运行，但是在子程序中就不是动态数组了，求助各位大佬，如何解决这个问题。我需要动态数组作为形参的目的是，这个数组太大了，如果定义静态数组会超过堆栈内存，我需要定义动态数组来计算a。感谢感谢！！！

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

program main

    use mod_assembly

    use mod_solver

    implicit none

    integer n

    real(kind=8),allocatable::a(:,:),b(:)

    real(kind=8),allocatable::x(:)

    real(kind=8) m,s

    n=4000

    allocate(a(n,n),b(n),x(n))

    call assembly(n,a,b)

    call lu_crout(a,b,x,n)

    end program main

---

作者: kyra    时间: 2022-7-30 11:00
1. 主程序中的动态数组，分配之后传递给子程序。无需额外操作，不会造成堆栈溢出。
2. 如果不需要在子程序中完成分配，释放。对子程序来说，虚参就是普通数组（实际是对应的实参可以是动态数组，也可以是普通数组）。
3. 如果要在子程序中完成分配，释放。则需要把虚参定义为可分配数组。
并且该子程序需要书写 interface 才能在主程序使用。（也可把子程序包含在module中，避免手动书写interface)
例如

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

subroutine assembly(a)

  real(8),allocatable::a(:,:)

end subroutine assembly



call assembly(a)

---

作者: Proteus    时间: 2022-7-30 15:41

kyra 发表于 2022-7-30 11:00
1. 主程序中的动态数组，分配之后传递给子程序。无需额外操作，不会造成堆栈溢出。
2. 如果不需要在子程序 ...

感谢感谢！！
我现在出现的问题是，在主程序中定义了动态数组a，并allocate了，再传递给子程序，在子程序中对a进行赋值等运算操作，还是会出现堆栈溢出的问题，请问这个要怎么解决呢？
请问我需要在子程序中书写interface吗？

---

作者: kyra    时间: 2022-7-30 16:51
原因可能不是因为主程序中分配的数组。

那你看看是否子程序中定义了其他（非虚参）的大数组？

如果没有，需要你给出足够多的，能让别人复现问题的代码。

---

作者: Proteus    时间: 2022-7-30 19:29

kyra 发表于 2022-7-30 16:51
原因可能不是因为主程序中分配的数组。

那你看看是否子程序中定义了其他（非虚参）的大数组？

您看看我的代码：
就是用LU分解求一个线性方程组。

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

program main

    use mod_assembly

    use mod_solver

    implicit none

    integer n

    real(kind=8),allocatable::a(:,:),b(:)

    real(kind=8),allocatable::x(:)

    n=40000

    allocate(a(n,n),b(n),x(n))

    call assembly(n,a,b)

    call lu_crout(a,b,x,n)

    end program main



module mod_assembly

    contains

    

    subroutine assembly(n,a,b)

    implicit none

    integer n

    real(kind=8) a(n,n),b(n)

    integer i,j

    a=0.0_8

    do i=1,n

        do j=1,n

            if (i.eq.j) then

                a(i,j)=3.0_8

            elseif (abs(i-j).eq.2) then

                a(i,j)=1.0_8

            elseif (abs(i-j).eq.1) then

                a(i,j)=2.0_8

            endif

        enddo

    enddo

    b=1.0_8

    return

    end subroutine assembly

    end module mod_assembly



module mod_solver

    contains

    subroutine lu_crout(A,b,xx,N)

    implicit none

    integer N

    real(kind=8) A(n,n),b(n),xx(n)

    real(kind=8) L(N,N),U(N,N)

    real(kind=8) y(N)

    call crout(A,L,U,N)

    call downtri(L,b,y,N)

    call uptri(U,y,xx,N)

    end subroutine lu_crout

    

    subroutine crout(A,L,U,N)

    implicit none

    integer N,i,k,r,j

    real(kind=8) s

    real(kind=8) A(n,n)

    real(kind=8) L(N,N),U(N,N)

    l=0.0_8

    u=0.0_8

    L(:,1)=a(:,1)

    U(1,:)=a(1,:)/L(1,1)

    do k=2,N

        do i=k,n

            s=0

            do r=1,k-1

                s=s+l(i,r)*u(r,k)

            end do

            l(i,k)=a(i,k)-s

        end do

        do j=k+1,n

            s=0

            do r=1,k-1

                s=s+l(k,r)*u(r,j)

            end do

            u(k,j)=(a(k,j)-s)/l(k,k)

        end do

        u(k,k)=1

    end do

    return

    end subroutine crout

    

    subroutine uptri(A,b,xx,N)

implicit none

integer::i,j,N

real(kind=8) A(n,n),b(N),xx(N)

xx(N)=b(N)/A(N,N)

do i=n-1,1,-1

    xx(i)=b(i)

   do j=i+1,N

    xx(i)=xx(i)-a(i,j)*xx(j)

   end do

    xx(i)=xx(i)/A(i,i)

end do

return

    end subroutine uptri

    

    subroutine downtri(A,b,xx,N)

    implicit none

    integer N,k,i

    real(kind=8) A(n,n),b(N),xx(N)

    xx(1)=b(1)/a(1,1)

    do k=2,N

        xx(k)=b(k)

        do i=1,k-1

            xx(k)=xx(k)-a(k,i)*xx(i)

        end do

        xx(k)=xx(k)/a(k,k)

    end do

    end subroutine downtri

    end module mod_solver

---

作者: kyra    时间: 2022-7-30 21:36
lu_crout 函数中     real(kind=8) L(N,N),U(N,N) 这两个（非虚参）数组占用了太多的堆栈。改成动态数组。

---

作者: Proteus    时间: 2022-7-31 09:52

kyra 发表于 2022-7-30 21:36
lu_crout 函数中     real(kind=8) L(N,N),U(N,N) 这两个（非虚参）数组占用了太多的堆栈。改成动态数组。 ...

感谢大佬！！

---

欢迎光临 Fortran Coder (http://bbs.fcode.cn/)

Powered by Discuz! X3.2