title: Data type mismatch in fortran
body: <p>I've written a rudimentary algorithm in Fortran 95 to calculate the gradient of a function (an example of which is prescribed in the code) using central differences augmented with a procedure known as Richardson extrapolation.  </p>
  
  <pre><code>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)) &lt; 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
  </code></pre>
  
  <p>In single precision it runs fine and gives me decent results.  But when I try to change to double precision as shown in the code, I get an error when trying to compile claiming that the assignment statement</p>
  
  <pre><code>d(j) = ( f(n,xhup) - f(n,xhdown) ) / (2.d0*h)
  </code></pre>
  
  <p>is producing a type mismatch <code>real(4)/real(8)</code>.  I have tried several different declarations of double precision, appended every appropriate double precision constant in the code with <code>d0</code>, and I get the same error every time.  I'm a little stumped as to how the function <code>f</code> is possibly producing a single precision number.</p>