Program c05pdae ! C05PDA Example Program Text ! Mark 24 Release. NAG Copyright 2012. ! .. Use Statements .. Use nag_library, Only: c05pda, dnrm2, nag_wp, x02ajf ! .. Implicit None Statement .. Implicit None ! .. Parameters .. Integer, Parameter :: n = 9, nout = 6 Integer, Parameter :: ldfjac = n ! .. Local Scalars .. Real (Kind=nag_wp) :: factor, fnorm, xtol Integer :: icount, ifail, irevcm, j, k, lr, mode ! .. Local Arrays .. Real (Kind=nag_wp), Allocatable :: diag(:), fjac(:,:), fvec(:), qtf(:), & r(:), w(:,:), x(:) Real (Kind=nag_wp) :: rwsav(10) Integer :: iwsav(15) Logical :: lwsav(2) ! .. Intrinsic Procedures .. Intrinsic :: sqrt ! .. Executable Statements .. Write (nout,*) 'C05PDA Example Program Results' Allocate (diag(n),fjac(ldfjac,n),fvec(n),qtf(n),r(n*(n+ & 1)/2),w(n,4),x(n)) ! The following starting values provide a rough solution. x(1:n) = -1.0E0_nag_wp xtol = sqrt(x02ajf()) diag(1:n) = 1.0E0_nag_wp mode = 2 factor = 100.0E0_nag_wp icount = 0 irevcm = 0 ifail = -1 revcomm: Do Call c05pda(irevcm,n,x,fvec,fjac,ldfjac,xtol,diag,mode,factor,r,lr, & qtf,w,lwsav,iwsav,rwsav,ifail) Select Case (irevcm) Case (1) icount = icount + 1 ! Insert print statements here to monitor progess if desired Cycle revcomm Case (2) ! Evaluate functions at current point fvec(1:n) = (3.0E0_nag_wp-2.0E0_nag_wp*x(1:n))*x(1:n) + 1.0E0_nag_wp fvec(2:n) = fvec(2:n) - x(1:(n-1)) fvec(1:(n-1)) = fvec(1:(n-1)) - 2.0E0_nag_wp*x(2:n) Cycle revcomm Case (3) ! Evaluate Jacobian at current point fjac(1:n,1:n) = 0.0E0_nag_wp Do k = 1, n fjac(k,k) = 3.0E0_nag_wp - 4.0E0_nag_wp*x(k) If (k/=1) Then fjac(k,k-1) = -1.0E0_nag_wp End If If (k/=n) Then fjac(k,k+1) = -2.0E0_nag_wp End If End Do Cycle revcomm Case Default Exit revcomm End Select End Do revcomm Select Case (ifail) Case (0) ! The NAG name equivalent of dnrm2 is f06ejf fnorm = dnrm2(n,fvec,1) Write (nout,*) Write (nout,99999) 'Final 2 norm of the residuals after', icount, & ' iterations is ', fnorm Write (nout,*) Write (nout,*) 'Final approximate solution' Write (nout,99998)(x(j),j=1,n) Case (3:) Write (nout,*) Write (nout,*) 'Approximate solution' Write (nout,99998)(x(j),j=1,n) End Select 99999 Format (1X,A,I4,A,E12.4) 99998 Format (5X,3F12.4) End Program c05pdae