newton_rel_T1 Subroutine

private impure subroutine newton_rel_T1(this, F, dFdx, x0, x_sol)

Type Bound

nlsolver

Arguments

Type IntentOptional Attributes Name
class(nlsolver), intent(inout) :: this
procedure(Fun2) :: F
procedure(dFun2), optional :: dFdx
real(kind=rk), intent(in), dimension(:) :: x0
real(kind=rk), intent(out), dimension(size(x0)) :: x_sol

Calls

proc~~newton_rel_t1~~CallsGraph proc~newton_rel_t1 forsolver::nlsolver%newton_rel_T1 proc~modified_newton_method_t1 forsolver::modified_newton_method_T1 proc~newton_rel_t1->proc~modified_newton_method_t1 proc~modified_quasi_fd_newton_method_t1 forsolver::modified_quasi_fd_newton_method_T1 proc~newton_rel_t1->proc~modified_quasi_fd_newton_method_t1 proc~newton_method_t1 forsolver::newton_method_T1 proc~newton_rel_t1->proc~newton_method_t1 proc~quasi_fd_newton_method_t1 forsolver::quasi_fd_newton_method_T1 proc~newton_rel_t1->proc~quasi_fd_newton_method_t1 interface~solve forsolver::solve proc~modified_newton_method_t1->interface~solve derivative derivative proc~modified_quasi_fd_newton_method_t1->derivative proc~modified_quasi_fd_newton_method_t1->interface~solve proc~newton_method_t1->interface~solve proc~quasi_fd_newton_method_t1->derivative proc~quasi_fd_newton_method_t1->interface~solve proc~solver_lin forsolver::solver_lin interface~solve->proc~solver_lin proc~gels_rel forsolver::gels_rel proc~solver_lin->proc~gels_rel proc~gesv_rel forsolver::gesv_rel proc~solver_lin->proc~gesv_rel interface~gels external_interfaces_solver::gels proc~gels_rel->interface~gels interface~gesv external_interfaces_solver::gesv proc~gesv_rel->interface~gesv

Called by

proc~~newton_rel_t1~~CalledByGraph proc~newton_rel_t1 forsolver::nlsolver%newton_rel_T1 none~solve forsolver::nlsolver%solve none~solve->proc~newton_rel_t1 program~test_solver10 test_solver10 program~test_solver10->none~solve program~test_solver11 test_solver11 program~test_solver11->none~solve program~test_solver12 test_solver12 program~test_solver12->none~solve program~test_solver13 test_solver13 program~test_solver13->none~solve program~test_solver14 test_solver14 program~test_solver14->none~solve program~test_solver3 test_solver3 program~test_solver3->none~solve program~test_solver4 test_solver4 program~test_solver4->none~solve program~test_solver5 test_solver5 program~test_solver5->none~solve program~test_solver6 test_solver6 program~test_solver6->none~solve program~test_solver7 test_solver7 program~test_solver7->none~solve program~test_solver8 test_solver8 program~test_solver8->none~solve program~test_solver9 test_solver9 program~test_solver9->none~solve

Source Code

   impure subroutine newton_rel_T1(this, F, dFdx, x0,  x_sol)
      interface
         impure function Fun2(x) result(res)
            import rk
            real(rk), dimension(:), intent(in)  :: x
            real(rk), dimension(:), allocatable :: res
         end function Fun2

         impure function dFun2(x) result(res)
            import rk
            real(rk), dimension(:), intent(in)    :: x
            real(rk), dimension(:,:), allocatable :: res
         end function dFun2
      end interface

      procedure(Fun2)            :: F
      procedure(dFun2), optional :: dFdx

      class(nlsolver),               intent(inout) :: this
      real(rk), dimension(:),        intent(in)    :: x0
      real(rk), dimension(size(x0)), intent(out)   :: x_sol
      integer                                      :: i

      if (this%verbosity == 1) then
         print '(a)', '-----------------------------------------------'
         print '(a)', 'maxit             tol'
         print '(g0, 10x, e12.4)', this%maxit, this%TolFun
         print '(a)', '-----------------------------------------------'
         print '(a)', 'start newton'
         print '(a)', '-----------------------------------------------'
         print '(a)', 'it     ||F||'
      end if

      select case (this%nl_method)
       case ('newton')
         call newton_method_T1(this, F, dFdx, x0,  x_sol)
       case ('newton-modified')
         call modified_newton_method_T1(this, F, dFdx, x0,  x_sol)
       case ('newton-quasi-fd')
         call quasi_fd_newton_method_T1(this, F, x0,  x_sol)
       case ('newton-quasi-fd-modified')
         call modified_quasi_fd_newton_method_T1(this, F, x0,  x_sol)
         ! case ('newton-quasi-bfgs')
         !    call quasi_bfgs_newton_method_T1(this, F, x0,  x_sol)
         ! case ('newton-quasi-bfgs-modified')
         !    call modified_quasi_bfgs_newton_method_T1(this, F, x0,  x_sol)
      end select

      if (this%verbosity == 1) then
         print '(a)', '-----------------------------------------------'
         print '(a)', 'end newton'
         print '(a)', '-----------------------------------------------'
         do i = 1,size(x_sol)
            print '(a, g0)', 'x_sol = ', x_sol(i)
         end do
      end if

   end subroutine newton_rel_T1