如何使用 Fortran 以更聪明的方式计算复杂函数？

Question

我想使用fortran：

其中ψ（ psi ）是一个复杂的fortran变量。 现在，我通过定义两个新的复杂变量来解决此问题：
ir =（1.0,0.0）和 ii =（0.0,1.0）。
我使用这些仅选择方程的真实或虚构部分。通过这种方式，我将方程式分别求解为真实和虚构的部分。代码在这里：

do i = 1,nn    
  mod2 = (abs(psi(i)))**2
  psi(i) = ir*(-beta*imag(der2(i)) + alpha*mod2*imag(psi(i))) + ii*(beta*real(der2(i)) - alpha*mod2*real(psi(i)))
end do

其中 psi 和 der2 是复杂的数组， nn elements。
我想以更好的方式求解此方程，而不会将其分为两个方程式。我试图以这种方式解决它：

mod2 = abs(psi)**2
psi = -ii*(-beta*der2+alpha*mod2*psi)

但是它不起作用，因为相对于我使用的第一种方法，我获得了完全不同的值。对我来说，这是不起作用的，因为在第二种方法中，我没有评估实际部分。这是对的吗？
例如，我的PSI阵列的10°元素变为：\

1. （-6.39094355774850412E-003，-6.04041029332164168E-003）（使用1°方法） 003
2. ）（2 °方法）\

有什么建议？ 谢谢你！

Answer 1

问题是 mod2 = abs（psi）** 2 应该是 mod2 = abs（dat）** 2

---

，但我认为 mod2的正确计算是 sum（real（dat*conjg（dat）））

我通过两种方法获得了相同的结果，并带有一些任意数据：

eq1=
 (-595.000000000000,-120.200000000000)
 (-713.800000000000,-1.40000000000000)
 (-832.600000000000,117.400000000000)
 (-951.400000000000,236.200000000000)
 (-1070.20000000000,355.000000000000)
 (-1189.00000000000,473.800000000000)
 (-1307.80000000000,592.600000000000)
 (-1426.60000000000,711.400000000000)
 (-1545.40000000000,830.200000000000)
 (-1664.20000000000,949.000000000000)
 eq2=
 (-595.000000000000,-120.200000000000)
 (-713.800000000000,-1.40000000000000)
 (-832.600000000000,117.400000000000)
 (-951.400000000000,236.200000000000)
 (-1070.20000000000,355.000000000000)
 (-1189.00000000000,473.800000000000)
 (-1307.80000000000,592.600000000000)
 (-1426.60000000000,711.400000000000)
 (-1545.40000000000,830.200000000000)
 (-1664.20000000000,949.000000000000)

。

program Console1
use,intrinsic :: iso_fortran_env
implicit none
! Variables
integer, parameter :: sp=real32, wp=real64
complex(wp), parameter :: ir = (1d0,0d0), ii = (0d0,1d0)
real(wp), parameter :: alpha = 0.1d0, beta = 0.2d0
integer, parameter :: n=10
complex(wp) :: dat(n), psi(n), der2(n)
integer :: i

! Body of Console1
    
    dat = [ ( (2d0-i)*ir - (4d0+i)*ii, i=1,n ) ]
    der2 = [ ( -(5d0+i)*ir + (1d0-i)*ii, i=1,n ) ]

    psi = eq1(dat, der2)
    
    print *, "eq1="
    do i=1,n
    print *, psi(i)
    end do

    psi = eq2(dat, der2)
    
    print *, "eq2="
    do i=1,n
    print *, psi(i)
    end do
    

contains

function eq1(dat, der2) result(psi)
complex(wp), intent(in) :: dat(:), der2(size(psi))
complex(wp) :: psi(size(dat))
real(wp) :: mod2        
integer :: i
    mod2 = sum( real( dat*conjg(dat) ) )
    do i=1, size(psi)
        psi(i) = ir*(-beta*imag(der2(i)) + alpha*mod2*imag(dat(i))) + ii*(beta*real(der2(i)) - alpha*mod2*real(dat(i)))
    end do
end function

function eq2(dat, der2) result(psi)
complex(wp), intent(in) :: dat(:), der2(size(dat))
complex(wp) :: psi(size(dat))
real(wp) :: mod2        
    mod2 = sum( real( dat*conjg(dat) ) )
    psi = -ii*(-beta*der2+alpha*mod2*dat)    
end function

end program Console1

还带有您的定义 |ψ|的定义 我们也有一致的结果。

eq1=
 (-13.0000000000000,-3.80000000000000)
 (-21.4000000000000,-1.40000000000000)
 (-34.6000000000000,3.40000000000000)
 (-53.8000000000000,11.8000000000000)
 (-80.2000000000000,25.0000000000000)
 (-115.000000000000,44.2000000000000)
 (-159.400000000000,70.6000000000000)
 (-214.600000000000,105.400000000000)
 (-281.800000000000,149.800000000000)
 (-362.200000000000,205.000000000000)
 eq2=
 (-13.0000000000000,-3.80000000000000)
 (-21.4000000000000,-1.40000000000000)
 (-34.6000000000000,3.40000000000000)
 (-53.8000000000000,11.8000000000000)
 (-80.2000000000000,25.0000000000000)
 (-115.000000000000,44.2000000000000)
 (-159.400000000000,70.6000000000000)
 (-214.600000000000,105.400000000000)
 (-281.800000000000,149.800000000000)
 (-362.200000000000,205.000000000000)

随着代码

function eq1(dat, der2) result(psi)
complex(wp), intent(in) :: dat(:), der2(size(psi))
complex(wp) :: psi(size(dat))
real(wp) :: mod2(size(dat))
integer :: i
    mod2 = abs(dat)**2 
    do i=1, size(psi)
        psi(i) = ir*(-beta*imag(der2(i)) + alpha*mod2(i)*imag(dat(i))) + ii*(beta*real(der2(i)) - alpha*mod2(i)*real(dat(i)))
    end do
end function

function eq2(dat, der2) result(psi)
complex(wp), intent(in) :: dat(:), der2(size(dat))
complex(wp) :: psi(size(dat))
real(wp) :: mod2(size(dat))     
    mod2 = abs(dat)**2 
    psi = -ii*(-beta*der2+alpha*mod2*dat)    
end function

---

返回原始问题，我可以将问题复制

 eq1=
 (-13.0000000000000,-3.80000000000000)
 (-21.4000000000000,-1.40000000000000)
 (-34.6000000000000,3.40000000000000)
 (-53.8000000000000,11.8000000000000)
 (-80.2000000000000,25.0000000000000)
 (-115.000000000000,44.2000000000000)
 (-159.400000000000,70.6000000000000)
 (-214.600000000000,105.400000000000)
 (-281.800000000000,149.800000000000)
 (-362.200000000000,205.000000000000)
 eq2=
 (0.000000000000000E+000,-1.20000000000000)
 (0.200000000000000,-1.40000000000000)
 (0.400000000000000,-1.60000000000000)
 (0.600000000000000,-1.80000000000000)
 (0.800000000000000,-2.00000000000000)
 (-1952.64000000000,779.256000000000)
 (NaN,NaN)
 (1.40000000000000,-2.60000000000000)
 (NaN,NaN)
 (1.80000000000000,-3.00000000000000)

到 mod2 = abs（psi）** 2 而不是 mod2 = abs（dat）** 2