Взятие производных полиномов Лежандра в Фортране - PullRequest
Новая
       87

Взятие производных полиномов Лежандра в Фортране

0 голосов
/ 07 декабря 2018

Мой код компилируется, но он не работает должным образом.Я не уверен, является ли его неспособность работать из-за математической ошибки и / или есть ли проблема в моем синтаксисе кодирования.Код показан ниже.Если кто-нибудь может мне понять, что не так, пожалуйста, дайте мне знать. 

Program received signal SIGSEGV: Segmentation fault - invalid memory reference.
Backtrace for this error:
#0  0x7F1E52D74E08
#1  0x7F1E52D73F90
#2  0x7F1E526BB4AF
#3  0x400BBC in legendrepoly.3381 at project.f90:?
#4  0x40138C in MAIN__ at project.f90:?
Segmentation fault (core dumped)



program Project

use iso_fortran_env
implicit none

integer(int32)               :: Nmax   ! maximum order of Legendre polynomial
integer(int32)               :: step   ! step size in angle from 0 to 180 degrees
real(real64), parameter      :: pi=3.1415926
real(real64)                 :: mu
real(real64)                 :: theta=0
real(real64), dimension(:,:), allocatable :: P, dP      ! Legendre polynomial
integer(int32)               :: k, h, i, n, l

! User Input
print *, "Please enter maximum order of Legendre polynomials:"
read (*,*) Nmax
print *, "Please enter integer step size (in degrees between 0 and 180):"
read (*,*) step

! Output
mu = cos(theta*(pi/180.0_real64))
k = 180.0/step
h = 2*cos(real(step)*(pi/180.0_real64))
allocate (  P(Nmax+1,180)  )
allocate ( dP(Nmax+1,180) )

call legendrepoly(l,mu,Nmax,step)
print *, dP

contains
subroutine legendrepoly(l,mu,Nmax,step)
real(real64), dimension(:,:), allocatable :: P, dP
real(real64)                 :: a
real(real64), intent(in)     :: mu
integer(int32)               :: i, Nmax, step
integer(int32), intent(in)   :: l

a = real(l)
do i = 1, 181, step
   P(1,i-1) = 1.0
   P(2,i-1) = mu
end do

do i = 1, 181, step
   do n = 2, Nmax
      P(l+1,i-1) = ((2.0*a+1.0) * mu * P(l,i-1) - a*P(l-1,i-1)) / (a+1.0)
      dP(l+1,i-1) = (P(l+1,i)-P(l+1,i-2))/2
   end do
end do

end subroutine legendrepoly

end program Project
...