From ef873b44d88f95d62822c0c50d3733753264bdea Mon Sep 17 00:00:00 2001 From: Seyed Ali Ghasemi Date: Sun, 9 Jul 2023 18:02:57 +0200 Subject: [PATCH] Fix bugs. Update tests. --- ford.yml | 2 +- fpm.toml | 106 ++++++++++++++++++++-------------------- src/forsolver.f90 | 36 +++++++------- test/test10.f90 | 8 +-- test/test11.f90 | 10 ++-- test/test12.f90 | 12 ++--- test/test13.f90 | 8 +-- test/test14.f90 | 10 ++-- test/test5.f90 | 2 +- test/test6.f90 | 2 +- test/test9.f90 | 6 +-- test/test_functions.f90 | 44 +++++++---------- 12 files changed, 119 insertions(+), 127 deletions(-) diff --git a/ford.yml b/ford.yml index d63601c..410475a 100644 --- a/ford.yml +++ b/ford.yml @@ -1,5 +1,5 @@ project: forsolver -version: v0.3.0 +version: v0.3.1 year: 2023 project_github: https://github.com/gha3mi/forsolver author: Seyed Ali Ghasemi diff --git a/fpm.toml b/fpm.toml index 70b91ef..07d79f7 100644 --- a/fpm.toml +++ b/fpm.toml @@ -1,5 +1,5 @@ name = "forsolver" -version = "0.3.0" +version = "0.3.1" author = "Seyed Ali Ghasemi" maintainer = "info@gha3mi.com" copyright = "Copyright (c) 2023, Seyed Ali Ghasemi" @@ -20,70 +20,70 @@ source-form = "free" kinds = {git="https://github.com/gha3mi/kinds.git"} fordiff = {git="https://github.com/gha3mi/fordiff.git"} -[[test]] -name = "test1" -source-dir = "test" -main = "test1.f90" +# [[test]] +# name = "test1" +# source-dir = "test" +# main = "test1.f90" -[[test]] -name = "test2" -source-dir = "test" -main = "test2.f90" +# [[test]] +# name = "test2" +# source-dir = "test" +# main = "test2.f90" -[[test]] -name = "test3" -source-dir = "test" -main = "test3.f90" +# [[test]] +# name = "test3" +# source-dir = "test" +# main = "test3.f90" -[[test]] -name = "test4" -source-dir = "test" -main = "test4.f90" +# [[test]] +# name = "test4" +# source-dir = "test" +# main = "test4.f90" -[[test]] -name = "test5" -source-dir = "test" -main = "test5.f90" +# [[test]] +# name = "test5" +# source-dir = "test" +# main = "test5.f90" -[[test]] -name = "test6" -source-dir = "test" -main = "test6.f90" +# [[test]] +# name = "test6" +# source-dir = "test" +# main = "test6.f90" -[[test]] -name = "test7" -source-dir = "test" -main = "test7.f90" +# [[test]] +# name = "test7" +# source-dir = "test" +# main = "test7.f90" -[[test]] -name = "test8" -source-dir = "test" -main = "test8.f90" +# [[test]] +# name = "test8" +# source-dir = "test" +# main = "test8.f90" -[[test]] -name = "test9" -source-dir = "test" -main = "test9.f90" +# [[test]] +# name = "test9" +# source-dir = "test" +# main = "test9.f90" -[[test]] -name = "test10" -source-dir = "test" -main = "test10.f90" +# [[test]] +# name = "test10" +# source-dir = "test" +# main = "test10.f90" -[[test]] -name = "test11" -source-dir = "test" -main = "test11.f90" +# [[test]] +# name = "test11" +# source-dir = "test" +# main = "test11.f90" -[[test]] -name = "test12" -source-dir = "test" -main = "test12.f90" +# [[test]] +# name = "test12" +# source-dir = "test" +# main = "test12.f90" -[[test]] -name = "test13" -source-dir = "test" -main = "test13.f90" +# [[test]] +# name = "test13" +# source-dir = "test" +# main = "test13.f90" [[test]] name = "test14" diff --git a/src/forsolver.f90 b/src/forsolver.f90 index 4f62447..1527bd8 100644 --- a/src/forsolver.f90 +++ b/src/forsolver.f90 @@ -210,7 +210,7 @@ end function dFun if (this%verbosity == 1) then print '(a)', '-----------------------------------------------' print '(a)', 'maxit x0 tol' - print '(i3, 10x, f12.8, 10x, e12.4)', this%maxit, x0, this%TolFun + print '(i5, 10x, f12.8, 10x, e12.4)', this%maxit, x0, this%TolFun print '(a)', '-----------------------------------------------' print '(a)', 'start newton' print '(a)', '-----------------------------------------------' @@ -261,8 +261,8 @@ impure function dFun(x) end function dFun end interface - procedure(Fun) :: F - procedure(dFun) :: dFdx + procedure(Fun) :: F + procedure(dFun), optional :: dFdx class(nlsolver), intent(inout) :: this real(rk), dimension(:), intent(in) :: x0 @@ -272,7 +272,7 @@ end function dFun if (this%verbosity == 1) then print '(a)', '-----------------------------------------------' print '(a)', 'maxit tol' - print '(i3, 10x, e12.4)', this%maxit, this%TolFun + print '(i5, 10x, e12.4)', this%maxit, this%TolFun print '(a)', '-----------------------------------------------' print '(a)', 'start newton' print '(a)', '-----------------------------------------------' @@ -327,7 +327,7 @@ end function Fun if (this%verbosity == 1) then print '(a)', '-----------------------------------------------' print '(a)', 'maxit x0 tol' - print '(i3, 10x, f12.8, 10x, e12.4)', this%maxit, real(x0, kind=rk), this%TolFun + print '(i5, 10x, f12.8, 10x, e12.4)', this%maxit, real(x0, kind=rk), this%TolFun print '(a)', '-----------------------------------------------' print '(a)', 'start newton' print '(a)', '-----------------------------------------------' @@ -373,7 +373,7 @@ end function Fun if (this%verbosity == 1) then print '(a)', '-----------------------------------------------' print '(a)', 'maxit tol' - print '(i3, 10x, f12.8, e12.4)', this%maxit, this%TolFun + print '(i5, 10x, f12.8, e12.4)', this%maxit, this%TolFun print '(a)', '-----------------------------------------------' print '(a)', 'start newton' print '(a)', '-----------------------------------------------' @@ -537,7 +537,7 @@ end function dFun criteriaFun = abs(F_val) if (this%verbosity == 1) then - print '(i3, f12.4, 4x, e12.4, 4x, e12.4)', k, xk, F_val, dFdx_val + print '(i5, f12.4, 4x, e12.4, 4x, e12.4)', k, xk, F_val, dFdx_val end if if (criteriaFun <= this%TolFun) then @@ -598,7 +598,7 @@ end function dFun criteriaFun = abs(F_val) if (this%verbosity == 1) then - print '(i3, f12.4, 4x, e12.4, 4x, e12.4)', k, xk, F_val, dFdx_val + print '(i5, f12.4, 4x, e12.4, 4x, e12.4)', k, xk, F_val, dFdx_val end if if (criteriaFun <= this%TolFun) then @@ -654,7 +654,7 @@ end function Fun criteriaFun = abs(F_val) if (this%verbosity == 1) then - print '(i3, f12.4, 4x, e12.4, 4x, e12.4)', k, xk, F_val, dFdx_val + print '(i5, f12.4, 4x, e12.4, 4x, e12.4)', k, xk, F_val, dFdx_val end if if (criteriaFun <= this%TolFun) then @@ -710,7 +710,7 @@ end function Fun criteriaFun = abs(F_val) if (this%verbosity == 1) then - print '(i3, f12.4, 4x, e12.4, 4x, e12.4)', k, xk, F_val, dFdx_val + print '(i5, f12.4, 4x, e12.4, 4x, e12.4)', k, xk, F_val, dFdx_val end if if (criteriaFun <= this%TolFun) then @@ -771,7 +771,7 @@ end function dFun criteriaFun = norm2(F_val) if (this%verbosity == 1) then - print '(i3, e12.4)', k, criteriaFun + print '(i5, e12.4)', k, criteriaFun end if if (criteriaFun <= this%TolFun) then @@ -832,7 +832,7 @@ end function dFun criteriaFun = norm2(F_val) if (this%verbosity == 1) then - print '(i3, e12.4)', k, criteriaFun + print '(i5, e12.4)', k, criteriaFun end if if (criteriaFun <= this%TolFun) then @@ -886,7 +886,7 @@ end function Fun criteriaFun = norm2(F_val) if (this%verbosity == 1) then - print '(i3, e12.4)', k, criteriaFun + print '(i5, e12.4)', k, criteriaFun end if if (criteriaFun <= this%TolFun) then @@ -941,7 +941,7 @@ end function Fun criteriaFun = norm2(F_val) if (this%verbosity == 1) then - print '(i3, e12.4)', k, criteriaFun + print '(i5, e12.4)', k, criteriaFun end if if (criteriaFun <= this%TolFun) then @@ -996,7 +996,7 @@ end function Fun criteriaFun = abs(F_val) if (this%verbosity == 1) then - print '(i3, f12.4, 4x, e12.4, 4x, e12.4)', k, real(xk, kind=rk), real(F_val, kind=rk), dFdx_val + print '(i5, f12.4, 4x, e12.4, 4x, e12.4)', k, real(xk, kind=rk), real(F_val, kind=rk), dFdx_val end if if (criteriaFun <= this%TolFun) then @@ -1050,7 +1050,7 @@ end function Fun criteriaFun = abs(F_val) if (this%verbosity == 1) then - print '(i3, f12.4, 4x, e12.4, 4x, e12.4)', k, real(xk, kind=rk), real(F_val, kind=rk), dFdx_val + print '(i5, f12.4, 4x, e12.4, 4x, e12.4)', k, real(xk, kind=rk), real(F_val, kind=rk), dFdx_val end if if (criteriaFun <= this%TolFun) then @@ -1104,7 +1104,7 @@ end function Fun criteriaFun = norm2(real(F_val, kind=rk)) if (this%verbosity == 1) then - print '(i3, e12.4)', k, criteriaFun + print '(i5, e12.4)', k, criteriaFun end if if (criteriaFun <= this%TolFun) then @@ -1158,7 +1158,7 @@ end function Fun criteriaFun = norm2(real(F_val, kind=rk)) if (this%verbosity == 1) then - print '(i3, e12.4)', k, criteriaFun + print '(i5, e12.4)', k, criteriaFun end if if (criteriaFun <= this%TolFun) then diff --git a/test/test10.f90 b/test/test10.f90 index 7bfa543..46adccb 100644 --- a/test/test10.f90 +++ b/test/test10.f90 @@ -7,15 +7,15 @@ program test10 implicit none type(nlsolver) :: nls - real(rk), dimension(3) :: x_sol + real(rk), dimension(2) :: x_sol call nls%set_options(& nl_method = 'newton-modified',& - nmp = 2,& + nmp = 1,& maxit = 100,& - TolFun = 1e-12_rk,& + TolFun = 1e-15_rk,& verbosity = 1) - call nls%solve(F=F3, dFdx=dF3dx, x0=[5.0_rk,3.0_rk,1.0_rk], x_sol=x_sol) + call nls%solve(F=F3, dFdx=dF3dx, x0=[-1.0_rk,-1.0_rk], x_sol=x_sol) end program test10 \ No newline at end of file diff --git a/test/test11.f90 b/test/test11.f90 index 9acbb04..a147add 100644 --- a/test/test11.f90 +++ b/test/test11.f90 @@ -7,16 +7,16 @@ program test11 implicit none type(nlsolver) :: nls - real(rk), dimension(3) :: x_sol + real(rk), dimension(2) :: x_sol call nls%set_options(& nl_method = 'newton-quasi-fd',& fdm_method = 'central',& - fdm_tol = 1e-8_rk,& - maxit = 100,& - TolFun = 1e-12_rk,& + fdm_tol = 1e-6_rk,& + maxit = 1000000,& + TolFun = 1e-13_rk,& verbosity = 1) - call nls%solve(F=F3, dFdx=dF3dx, x0=[5.0_rk,3.0_rk,1.0_rk], x_sol=x_sol) + call nls%solve(F=F3, x0=[-1.0_rk,-1.0_rk], x_sol=x_sol) end program test11 \ No newline at end of file diff --git a/test/test12.f90 b/test/test12.f90 index 00e8e29..6a289b5 100644 --- a/test/test12.f90 +++ b/test/test12.f90 @@ -7,17 +7,17 @@ program test12 implicit none type(nlsolver) :: nls - real(rk), dimension(3) :: x_sol + real(rk), dimension(2) :: x_sol call nls%set_options(& nl_method = 'newton-quasi-fd-modified',& fdm_method = 'central',& - fdm_tol = 1e-8_rk,& - nmp = 2,& - maxit = 100,& - TolFun = 1e-12_rk,& + fdm_tol = 1e-6_rk,& + nmp = 1,& + maxit = 1000000,& + TolFun = 1e-13_rk,& verbosity = 1) - call nls%solve(F=F3, dFdx=dF3dx, x0=[5.0_rk,3.0_rk,1.0_rk], x_sol=x_sol) + call nls%solve(F=F3, x0=[-1.0_rk,-1.0_rk], x_sol=x_sol) end program test12 \ No newline at end of file diff --git a/test/test13.f90 b/test/test13.f90 index fe20dea..944d36a 100644 --- a/test/test13.f90 +++ b/test/test13.f90 @@ -7,15 +7,15 @@ program test13 implicit none type(nlsolver) :: nls - complex(rk), dimension(3) :: x_sol + complex(rk), dimension(2) :: x_sol call nls%set_options(& nl_method = 'newton-quasi-cs',& cs_tol = 1e-200_rk,& - maxit = 100,& - TolFun = 1e-15_rk,& + maxit = 10800,& + TolFun = 1e-13_rk,& verbosity = 1) - call nls%solve(F=F4, x0=[(5.0_rk,0.0_rk) ,(3.0_rk,0.0_rk), (1.0_rk,0.0_rk)], x_sol=x_sol) + call nls%solve(F=F4, x0=[(-1.0_rk,0.0_rk) ,(-1.0_rk,0.0_rk)], x_sol=x_sol) end program test13 diff --git a/test/test14.f90 b/test/test14.f90 index 4fff404..6e14922 100644 --- a/test/test14.f90 +++ b/test/test14.f90 @@ -7,16 +7,16 @@ program test14 implicit none type(nlsolver) :: nls - complex(rk), dimension(3) :: x_sol + complex(rk), dimension(2) :: x_sol call nls%set_options(& nl_method = 'newton-quasi-cs-modified',& cs_tol = 1e-200_rk,& - nmp = 2,& - maxit = 100,& - TolFun = 1e-15_rk,& + nmp = 1,& + maxit = 10800,& + TolFun = 1e-13_rk,& verbosity = 1) - call nls%solve(F=F4, x0=[(5.0_rk,0.0_rk) ,(3.0_rk,0.0_rk), (1.0_rk,0.0_rk)], x_sol=x_sol) + call nls%solve(F=F4, x0=[(-1.0_rk,0.0_rk) ,(-1.0_rk,0.0_rk)], x_sol=x_sol) end program test14 diff --git a/test/test5.f90 b/test/test5.f90 index 3be941e..c433def 100644 --- a/test/test5.f90 +++ b/test/test5.f90 @@ -17,6 +17,6 @@ program test5 TolFun = 1e-15_rk,& verbosity = 1) - call nls%solve(F=F1, dFdx=dF1dx, x0=10.0_rk, x_sol=x_sol) + call nls%solve(F=F1, x0=10.0_rk, x_sol=x_sol) end program test5 \ No newline at end of file diff --git a/test/test6.f90 b/test/test6.f90 index 3e60f1c..b6cf513 100644 --- a/test/test6.f90 +++ b/test/test6.f90 @@ -18,6 +18,6 @@ program test6 TolFun = 1e-15_rk,& verbosity = 1) - call nls%solve(F=F1, dFdx=dF1dx, x0=10.0_rk, x_sol=x_sol) + call nls%solve(F=F1, x0=10.0_rk, x_sol=x_sol) end program test6 \ No newline at end of file diff --git a/test/test9.f90 b/test/test9.f90 index 23ae5b0..c14de40 100644 --- a/test/test9.f90 +++ b/test/test9.f90 @@ -7,14 +7,14 @@ program test9 implicit none type(nlsolver) :: nls - real(rk), dimension(3) :: x_sol + real(rk), dimension(2) :: x_sol call nls%set_options(& nl_method = 'newton',& maxit = 100,& - TolFun = 1e-12_rk,& + TolFun = 1e-15_rk,& verbosity = 1) - call nls%solve(F=F3, dFdx=dF3dx, x0=[5.0_rk,3.0_rk,1.0_rk], x_sol=x_sol) + call nls%solve(F=F3, dFdx=dF3dx, x0=[-1.0_rk,-1.0_rk], x_sol=x_sol) end program test9 \ No newline at end of file diff --git a/test/test_functions.f90 b/test/test_functions.f90 index 4c0bf3d..fbe1e8b 100644 --- a/test/test_functions.f90 +++ b/test/test_functions.f90 @@ -26,33 +26,25 @@ end function F2 function F3(x) result(F_val) real(rk), dimension(:), intent(in) :: x real(rk), dimension(:), allocatable :: F_val - integer :: i - allocate(F_val(3)) - F_val = 0.0_rk - do i = 1,3 - F_val(i) = x(i)**3 - end do + allocate(F_val(2)) + + F_val(1) = 2.0_rk*x(1) - 400.0_rk*x(1) * (x(2) - x(1)**2) - 2 + F_val(2) = 200.0_rk*x(2) - 200.0_rk*x(1)**2 + end function F3 function dF3dx(x) result(dFdx_val) real(rk), dimension(:), intent(in) :: x real(rk), dimension(:,:), allocatable :: dFdx_val - integer :: i, j - real(rk), dimension(3,3) :: Idt - - Idt = 0.0_rk - Idt(1,1) = 1.0_rk - Idt(2,2) = 1.0_rk - Idt(3,3) = 1.0_rk - - allocate(dFdx_val(3,3)) - dFdx_val = 0.0_rk - do i = 1,3 - do j = 1,3 - dFdx_val(i,j) = dFdx_val(i,j) + 3.0_rk*Idt(i,j)*x(i)**2 - end do - end do + + allocate(dFdx_val(2,2)) + + dFdx_val(1,1) = 1200.0_rk*x(1)**2 - 400.0_rk*x(2) + 2 + dFdx_val(1,2) = - 400.0_rk*x(1) + + dFdx_val(2,1) = - 400.0_rk*x(1) + dFdx_val(2,2) = 200.0_rk end function dF3dx @@ -61,11 +53,11 @@ function F4(x) result(F_val) complex(rk), dimension(:), allocatable :: F_val integer :: i - allocate(F_val(3)) - F_val = 0.0_rk - do i = 1,3 - F_val(i) = x(i)**3 - end do + allocate(F_val(2)) + + F_val(1) = 2.0_rk*x(1) - 400.0_rk*x(1) * (x(2) - x(1)**2) - 2 + F_val(2) = 200.0_rk*x(2) - 200.0_rk*x(1)**2 + end function F4 end module functions_module