2016-11-18 4 views
0

저는 fortran을 사용하고 있으며 그 요소가 기능을하는 행렬의 파생물을 얻으려고합니다.요소가 기능을하는 행렬을 유도하는 방법은 무엇입니까?

program derivada_matrix 

    integer, parameter :: matrix_size = 5 
    integer :: i,j 
    real(8) :: time = 1.0 
    real(8),dimension (matrix_size, matrix_size) :: W 
    real(8),dimension (matrix_size, matrix_size) :: dW 

    call potent(time,W) 
    do i = 1, matrix_size 
     do j=1, matrix_size 
     call Derivada(time,W(i,j),dW(i,j)) 
     end do 
    end do 

    print*, 'matrix' 
    print*, W 
    print*, 'derivada', dW 


    end program 


    Subroutine Derivada (x1,F,D) 
    implicit none 
    Real*8 :: x1 
    Real*8 :: h= 1.0E-6 
    integer, parameter :: matrix_size = 5 
    real*8 :: D,F 
    external F 
    D = (1.0*F(x1-2*h) - 8.0*F(x1-h) + 8.0*F(x1+h) - 1.0*F(x1+2*h))/(12.0*h) 

return 
End subroutine Derivada 

subroutine potent(T,W) 
implicit none 
integer, parameter :: matrix_size = 5 
real(8),dimension(matrix_size,matrix_size) :: W 
Real(8):: T 
integer :: i,j 
do i = 1, matrix_size 
    do j=i,matrix_size 
     W(i,j) = 0.0 
     W(j,i) = W(i,j) 
    end do 
    W(i,i) = cos(T) 
end do 
RETURN 
END subroutine potent 

는 기본적으로 제 루틴은 다른 주 대각선에 제로 함수 (코사인)와이 유도체를하기로되어있어 상기 제 서브 루틴 검사 행렬을 생성한다. 이 내가

 call Derivada(time,W(i,j),dW(i,j)) 
         1 
    Warning: Expected a procedure for argument ‘f’ at (1) 

내가 두 번째 호출에 얻을 오류/경고 메시지

를 얻을 오류/경고 메시지입니다. W 매트릭스를 만들면 그 속성이 함수로 손실되고 각 요소를 파생시키기 위해 두 번째 호출에서 인수로 사용할 수 없기 때문에 추측합니다. 어떻게 개선 할 수 있습니까? 어떻게 프로그램/함수 서브 루틴을 만들어 그 입력이이 같은 행렬이고 출력이 미분이 될 수 있습니까 ??? 코드는하지 않습니다

감사

답변

0

당신은이하는 생각. W은 실제 배열입니다. 그냥 T 지점에서 계산 된 cos로 항목을 채우는 것입니다. 함수 배열이 아닙니다.

여기에 당신이하고 싶은 것에 대한보다 현대적인 접근 방식이 있습니다. extern 함수를 사용하는 대신 프로 시저 포인터를 사용하고 파생 형식으로 래핑하여 배열을 만들 수 있습니다. 그런 다음 파생물을 가져오고 자하는 함수를 가리 킵니다.

이 접근법의 한계는 인터페이스를 준수해야한다는 것입니다.이 경우 두 개의 인수를 허용하는 함수에서만 작동합니다. 선택적 인수를 사용하거나 파생 된 유형으로 래핑하거나 여러 인터페이스를 사용하는 등 여러 가지 방법으로이를 확장 할 수 있습니다.하지만 이는 귀하에게 달려 있습니다.

module deriv 
    use, intrinsic :: iso_fortran_env, only: dp => real64 
    implicit none 

    ! Wrapper for a procedure pointer so we can make an array of them 
    type ptr_wrapper 
    procedure(f), nopass, pointer :: func 
    end type ptr_wrapper 


    ! Use an interface for the function rather than "extern" functions 
    ! Change the interface to whatever you want. You could take multiple reals, integers, etc 
    abstract interface 
    pure function f(x1) 
     import 
     real(dp), intent(in) :: x1 
     real(dp) :: f 
    end function f 
    end interface 


contains 

    function derivada(x1, fx) result(d) 
    implicit none 
    real(dp), intent(in) :: x1 
    procedure(f), pointer :: fx 
    real(dp) :: d 

    real(dp) :: h = 1.0E-6 

    d = (1.0*fx(x1-2*h) - 8.0*fx(x1-h) + 8.0*fx(x1+h) - 1.0*fx(x1+2*h))/(12.0*h) 

    end function derivada 


    pure function test_func(x1) 
    real(dp), intent(in) :: x1 
    real(dp) test_func 

    test_func = 2d0 * x1 
    end function test_func 

    pure function test_func2(x1) 
    real(dp), intent(in) :: x1 
    real(dp) test_func2 

    test_func2 = x1**2 
    end function test_func2 


end module deriv 

program derivada_matrix 
    use, intrinsic :: iso_fortran_env, only: dp => real64 
    use deriv 
    implicit none 

    integer i, j 
    real(dp) :: dx 


    type(ptr_wrapper) :: w(2,2) 

    !Point to the function that you want each index to point to 
    w(1,1)%func => test_func 
    w(1,2)%func => test_func2 
    ! Pointing to the generic function wasn't working so I had to point to the specific double value trig functions 
    w(2,1)%func => dcos 
    w(2,2)%func => dsin 

    do i = 1, 2 
    do j = 1, 2 
     dx = derivada(1d0, w(i,j)%func) 
     print *, dx 
    end do 
    end do 


end program derivada_matrix 

그리고 그것은 이들 기능에 대한 미분을 출력 X = 1.0

2.00000000000000  
    2.00000000011102  
-0.841470984776962  
    0.540302305853706 
+0

내가 노력 할게요 내가 너희를 알려 드리겠습니다. 감사! – Daniel

+0

Pd : 제네릭 기능을 가리 키지 않는 이유에 대해 이해할 수 없습니다. 결과를 올바르게 인쇄 할 수 있었고 작동하지 않는 것은 무엇 이었습니까? – Daniel

+0

@Daniel ifort를 사용하여 제네릭 함수를 가리킬 때 10^-300과 같은 별난 값이 생겼습니다. 17 – Exascale