Tipo de dato no coincidente en fortran

Question

May 08, 2013, 07:07 PM

gradient double-precision fortran95 floating-point

Tipo de dato no coincidente en fortran

He escrito un algoritmo rudimentario en Fortran 95 para calcular el gradiente de una función (un ejemplo de lo cual se prescribe en el código) utilizando diferencias centrales aumentadas con un procedimiento conocido como extrapolación de Richardson.

function f(n,x)
! The scalar multivariable function to be differentiated

integer :: n
real(kind = kind(1d0)) :: x(n), f

f = x(1)**5.d0 + cos(x(2)) + log(x(3)) - sqrt(x(4))

end function f
!=====!
!=====!
!=====!

program gradient
!==============================================================================!
! Calculates the gradient of the scalar function f at x=0using a finite        !
! difference approximation, with a low order Richardson extrapolation.         !
!==============================================================================!

parameter (n = 4, M = 25)
real(kind = kind(1d0)) :: x(n), xhup(n), xhdown(n), d(M), r(M), dfdxi, h0, h, gradf(n)

h0 = 1.d0
x  = 3.d0

! Loop through each component of the vector x and calculate the appropriate
! derivative
do i = 1,n
! Reset step size
h = h0

    ! Carry out M successive central difference approximations of the derivative
do j = 1,M
        xhup = x
        xhdown = x
        xhup(i) = xhup(i) + h
        xhdown(i) = xhdown(i) - h
        d(j) = ( f(n,xhup) - f(n,xhdown) ) / (2.d0*h)
    h = h / 2.d0
    end do

    r = 0.d0
    do k = 3,M      r(k) = ( 64.d0*d(k) - 20.d0*d(k-1) + d(k-2) ) / 45.d0
        if ( abs(r(k) - r(k-1)) < 0.0001d0 ) then
        dfdxi = r(k)
            exit
        end if
    end do

    gradf(i) = dfdxi
end do

! Print out the gradient
write(*,*) " "
write(*,*) " Grad(f(x)) = "
write(*,*) " "
do i = 1,n
    write(*,*) gradf(i)
end do

end program gradient

En precisión simple funciona bien y me da resultados decentes. Pero cuando intento cambiar a doble precisión como se muestra en el código, recibo un error al intentar compilar afirmando que la declaración de asignación

d(j) = ( f(n,xhup) - f(n,xhdown) ) / (2.d0*h)

está produciendo un desajuste de tiporeal(4)/real(8). He intentado varias declaraciones diferentes de precisión doble, agregando cada constante de precisión doble apropiada en el código cond0, y me sale el mismo error cada vez. Estoy un poco perplejo en cuanto a cómo funciona la funciónf posiblemente esté produciendo un solo número de precisión.

Deja tu comentario