modified_quasi_cs_newton_method_T1 Subroutine

private impure subroutine modified_quasi_cs_newton_method_T1(this, F, x0, x_sol)

Arguments

Type IntentOptional Attributes Name
class(nlsolver), intent(inout) :: this
procedure(Fun17) :: F
complex(kind=rk), intent(in), dimension(:) :: x0
complex(kind=rk), intent(out), dimension(size(x0)) :: x_sol

Calls

proc~~modified_quasi_cs_newton_method_t1~~CallsGraph proc~modified_quasi_cs_newton_method_t1 forsolver::modified_quasi_cs_newton_method_T1 derivative derivative proc~modified_quasi_cs_newton_method_t1->derivative interface~solve forsolver::solve proc~modified_quasi_cs_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~~modified_quasi_cs_newton_method_t1~~CalledByGraph proc~modified_quasi_cs_newton_method_t1 forsolver::modified_quasi_cs_newton_method_T1 proc~newton_complex_step_rel_t1 forsolver::nlsolver%newton_complex_step_rel_T1 proc~newton_complex_step_rel_t1->proc~modified_quasi_cs_newton_method_t1 none~solve forsolver::nlsolver%solve none~solve->proc~newton_complex_step_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 modified_quasi_cs_newton_method_T1(this, F, x0,  x_sol)
      interface
         impure function Fun17(x) result(res)
            import rk
            complex(rk), dimension(:), intent(in)  :: x
            complex(rk), dimension(:), allocatable :: res
         end function Fun17
      end interface

      procedure(Fun17)  :: F

      class(nlsolver),                  intent(inout) :: this
      complex(rk), dimension(:),        intent(in)    :: x0
      complex(rk), dimension(size(x0)), intent(out)   :: x_sol
      complex(rk), dimension(size(x0))                :: xk
      complex(rk), dimension(:),   allocatable        :: F_val
      real(rk),    dimension(:,:), allocatable        :: dFdx_val
      real(rk)                                        :: criteriaFun
      integer                                         :: k
      logical                                         :: convergenz
      real(rk), dimension(size(x0))                   :: dk
      real(rk)                                        :: alphak

      k          = 0
      xk         = x0
      convergenz = .false.

      do while (.not. convergenz .and. k < this%maxit)
         F_val    = F(xk)
         if (mod(k, this%nmp) == 0) dFdx_val = derivative(f=F, x=real(xk,kind=rk), h=this%cs_tol)

         criteriaFun = norm2(real(F_val, kind=rk))

         if (this%verbosity == 1) then
            print '(g0, e12.4)', k, criteriaFun
         end if

         if (criteriaFun <= this%TolFun) then
            convergenz = .true.
            x_sol      = xk
            return
         else
            dk = - solve(dFdx_val, real(F_val, kind=rk), this%lin_method)
            alphak = 1.0_rk
            xk = xk + alphak*dk
            k = k + 1
         end if
      end do
   end subroutine modified_quasi_cs_newton_method_T1