!*******************************************************************************
!* This program finds a machine zero.
!*
!*------------------------------------------------------------------------------
!* Input: 'given_number' = a number of your choice
!*
!* Output: Machine zero and the exponent of 2^i representing the machine zero
!*------------------------------------------------------------------------------
!*
!* Definition: "Machine zero" is defined as the largest number that does not
!* change a given number if added to it.
!*
!* Note: In a formal definition of machine zero or machine epsilon,
!* 'given number' is 1.0. Here, you can choose it.
!*
!* Note: Suppose that you iterate by u_{new} = u_{old} + du.
!* You say "It's converged" when du is small enough, i,e,
!* if du does not change u_{old}. That happens when the machine
!* cannot recognize du (too small), so that u_{old} doesn't change at all.
!* And we may say "It's converged to machine zero".
!* But then the machine zero depends on the value of u_{old}.
!* You can experience it with this code.
!*
!* Note: We make a loop over i=1,2,3,4,....., and at each i,
!* add 2^(-i) to a given number : given number + 2^(-i)
!* compare the result with the given_number.
!* If the given_number doesn't change at all, then,
!* the machine cannot recognize 2^(-i), i.e., machine zero.
!*
!* Examples:
!*
!* (1)Type in the number for which the machine zero is sought -> 1
!* Machine zero = 1.110223024625157E-016
!*
!* (2)Type in the number for which the machine zero is sought -> 10000000
!* Machine zero = 9.313225746154785E-010
!*
!* (3)Type in the number for which the machine zero is sought -> 1.0e+10
!* Machine zero = 9.536743164062500E-007
!*
!*
!* Version 02
!* 11-29-12: Definition of "machine zero" is given.
!* Modified the code to find "machine zero" for a given number.
!*
!* Katate Masatsuka, February 2009. http://www.cfdbooks.com
!*******************************************************************************
program machine_zero
implicit none
integer, parameter :: sp = kind(1.0)
integer, parameter :: dp = selected_real_kind(2*precision(1.0_sp))
! Note: Replace all 'dp' below by 'sp' to find machine zero in single precision.
real(dp) :: given_number !<- Input number (for which machine zero is sought.)
real(dp) :: temp, two = 2.0_dp
integer :: i
write(*,*)
write(*,'(a)',advance="no") "Type in the number for which the machine zero is sought -> "
read(*,*) given_number
write(*,*)
!Find the machine zero by adding 2^i (i=1,2,3,...) until given_number is not changed.
i = 0
do
i = i + 1
temp = given_number + two**(-i)
if (temp == given_number) exit !<- This means that 2**(-i) is considered as zero.
end do
write(*,*) "Machine zero = ", two**(-i)
write(*,*) "which is equal to 2^i: ", " i =", i
write(*,*)
end program machine_zero
!*******************************************************************************