program newton
c
c Use a Newton iteration to solve the equation
c
c John Mahaffy 2/14/95
c
c x**3+x-10=0
c
c xin - input initial guess to the solution
c x - current approximation to the solution
c f(x) - x**3+x-10
c dfdx - derivative of f with respect to x
c xo - previous guess for solution
c eps - convergence criterion
c dx - change in solution approximation
c it - number of iterations
c itmax - maximum number of iterations
c
c
implicit none
real x,f,dfdx,xo,eps,fx,dx,xin,abs
integer it,itmax
c c
c Use fortran statement functions to define f and its derivative
c Note that these statements must appear before other executable
c statements
c
c Now start executable fortran statements
c
eps=1.e-5
itmax=10
print *,'What is your initial guess for the solution?'
read *, xin
x=xin
c
c The following is a more traditional Fortran do loop. It is a more
c natural form for this particular task. It executes the instructions
c from the "do" through the instruction with label "100" for all
c values of "it" between 1 and itmax. Try setting itmax=0 and see
c what happens. If itmax=0 are the contents of the loop executed?
c
write(*,'(/,"Simple do loop",/)')
c
c Note the direct use of a format above rather than a reference to a
c format number. This is only useful for one time use of a given
c format structure
c
do 100 it=1,itmax
xo=x
fx=f(xo,dfdx)
c ^^^^^^^
c ^ This area between the parentheses is often refered
c ^ to as the calling sequence, and contains arguments
c ^ passed as input to the function and passed back as
c ^ output from the function
c
c
dx = -fx/dfdx
x=xo+dx
write (*,2000)it,x,fx,dx
if(abs(x-xo).lt.eps*abs(x)) go to 200
100 continue
print *, 'Iteration failed to converge'
200 x=xin
c
c I could loop over only even values of "it" with the following structure.
c It basically says loop over values of "it" from 0 to itmax in steps of 2.
c Check to see what happens when you set itmax=3 and input a starting guess
c of xo=6.0. Yes, this is a pointless way to do a Newton Iteration. It
c is just here as another example of the workings of DO loop indices.
c
write(*,'(/,"Do loop with a special increment on the index",/)')
do 250 it=0,itmax,2
xo=x
fx=f(xo,dfdx)
c ^^^^^^^
c ^ This area between the parentheses is often refered to as the
c ^ calling sequence, and contains arguments passed as input to
c ^ the function and passed back as output from the function
c
dx = -fx/dfdx
x=xo+dx
write (*,2000)it,x,fx,dx
if(abs(x-xo).lt.eps*abs(x)) go to 300
250 continue
c
c Note the value of "it" after the loop is completed. Because of the step
c of 2 "it" never directly gives the number of times you've been through
c the loop.
c
print *,'------------------------------------------------'
print *,' Value of loop index after exit from the loop'
300 print *, 'it = ', it
c
c
c
500 stop
2000 format(1x,'it=',i3,', x=',f8.3,', f(x)=',f8.3,', dx=',f8.3)
end
function f(x,deriv)
c ^^^^^^^
c ^ These are the "dummy" arguments of the function. They need
c ^ not be the same variable names as those appearing
c ^ in the calling sequence where the function is referenced.
c ^ However the number of arguments here should match the
c ^ number in the reference to the function (There are
c ^ exceptions to this rule but worry about that when you
c ^ become much better with Fortran).
c ^ It is also very important to match the data types between
c ^ corresponding arguments here and in the calling sequence
c ^ (There may be times later when you violate this rule too
c ^ but, again, you'll need to understand why you are tricking
c ^ the machine.
c
c Evaluate the function f(x)=x**3+x-10
c also return the derivative of the function
c
implicit none
real f,x,deriv
c
c John Mahaffy 2/14/95
c
c Function arguments can receive values from the current reference to
c the function, and can also pass data back to the program using the
c function. Documentation of arguments should clearly delimit those
c receiving input for the function from those transmitting output.
c
c INPUT
c
c x - independent variable
c
c OUTPUT
c
c f - value of the function
c deriv - derivative of the function with respect to x
c
deriv=3*x**2+1
f=x**3+x-10
return
end