refactor lmdif1 example

Author: jacobwilliams

examples/example_lmdif1.f90

@@ -1,61 +1,59 @@
-module testmod_dif1
+program example_lmdif1
+
+use minpack_module, only: wp, enorm, lmdif1
+use iso_fortran_env, only: nwrite => output_unit
+
 implicit none
-private
-public fcn, dp
 
-integer, parameter :: dp=kind(0d0)
+integer, parameter :: n = 3
+integer, parameter :: m = 15
+integer, parameter :: lwa = m*n+5*n+m
+
+integer :: info
+real(wp) :: tol, x(n), fvec(m), wa(lwa)
+integer :: iwa(n)
+
+! The following starting values provide a rough fit.
+x = [1.0_wp, 1.0_wp, 1.0_wp]
+
+! Set tol to the square root of the machine precision. Unless high precision
+! solutions are required, this is the recommended setting.
+tol = sqrt(epsilon(1.0_wp))
+
+call lmdif1(fcn, m, n, x, fvec, tol, info, iwa, wa, lwa)
+
+write(nwrite, '(5x,a,d15.7//,5x,a,16x,i10//,5x,a//(5x,3d15.7))') &
+        'FINAL L2 NORM OF THE RESIDUALS', enorm(m, fvec), &
+        'EXIT PARAMETER', info, &
+        'FINAL APPROXIMATE SOLUTION', x
 
 contains
 
 subroutine fcn(m, n, x, fvec, iflag)
-integer, intent(in) :: m, n
-integer, intent(inout) :: iflag
-real(dp), intent(in) :: x(n)
-real(dp), intent(out) :: fvec(m)
-
-integer :: i
-real(dp) :: tmp1, tmp2, tmp3, y(15)
-! Suppress compiler warning:
-y(1) = iflag
-y = [1.4D-1, 1.8D-1, 2.2D-1, 2.5D-1, 2.9D-1, 3.2D-1, 3.5D-1, 3.9D-1, &
-        3.7D-1, 5.8D-1, 7.3D-1, 9.6D-1, 1.34D0, 2.1D0, 4.39D0]
-
-do i = 1, 15
-    tmp1 = i
-    tmp2 = 16 - i
-    tmp3 = tmp1
-    if (i > 8) tmp3 = tmp2
-    fvec(i) = y(i) - (x(1) + tmp1/(x(2)*tmp2 + x(3)*tmp3))
-end do
-end subroutine
-
-end module
 
+    integer, intent(in) :: m
+    integer, intent(in) :: n
+    real(wp), intent(in) :: x(n)
+    real(wp), intent(out) :: fvec(m)
+    integer, intent(inout) :: iflag
 
-program example_lmdif1
-use minpack_module, only: enorm, lmdif1
-use testmod_dif1, only: dp, fcn
-implicit none
+    integer :: i
+    real(wp) :: tmp1, tmp2, tmp3
 
-integer :: info, m, n
-real(dp) :: tol, x(3), fvec(15)
-integer :: iwa(size(x))
-real(dp), allocatable :: wa(:)
+    real(wp),parameter :: y(15) = [1.4e-1_wp, 1.8e-1_wp, 2.2e-1_wp, 2.5e-1_wp, 2.9e-1_wp, &
+                                   3.2e-1_wp, 3.5e-1_wp, 3.9e-1_wp, 3.7e-1_wp, 5.8e-1_wp, &
+                                   7.3e-1_wp, 9.6e-1_wp, 1.34e0_wp, 2.1e0_wp, 4.39e0_wp]
 
-! The following starting values provide a rough fit.
-x = [1._dp, 1._dp, 1._dp]
+    if (iflag > 0) then ! just to avoid the compiler warning
+        do i = 1, 15
+            tmp1 = i
+            tmp2 = 16 - i
+            tmp3 = tmp1
+            if (i > 8) tmp3 = tmp2
+            fvec(i) = y(i) - (x(1) + tmp1/(x(2)*tmp2 + x(3)*tmp3))
+        end do
+    end if
 
-! Set tol to the square root of the machine precision. Unless high precision
-! solutions are required, this is the recommended setting.
-tol = sqrt(epsilon(1._dp))
-
-m = size(fvec)
-n = size(x)
-allocate(wa(m*n + 5*n + m))
-call lmdif1(fcn, size(fvec), size(x), x, fvec, tol, info, iwa, wa, size(wa))
-print 1000, enorm(size(fvec), fvec), info, x
-1000 format(5x, 'FINAL L2 NORM OF THE RESIDUALS', d15.7 // &
-            5x, 'EXIT PARAMETER', 16x, i10              // &
-            5x, 'FINAL APPROXIMATE SOLUTION'            // &
-            5x, 3d15.7)
-end program
+end subroutine fcn
+
+end program example_lmdif1