solve_nonlinear_system Subroutine

public subroutine solve_nonlinear_system(n, guesses, residual_fn, x_best, success)

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: n
real(kind=dp), intent(in) :: guesses(:,:)
procedure(nonlinear_residual) :: residual_fn
real(kind=dp), intent(out) :: x_best(n)
logical, intent(out) :: success

Calls

proc~~solve_nonlinear_system~~CallsGraph proc~solve_nonlinear_system solve_nonlinear_system proc~try_newton_solve try_newton_solve proc~solve_nonlinear_system->proc~try_newton_solve proc~numerical_jacobian numerical_jacobian proc~try_newton_solve->proc~numerical_jacobian proc~residual_norm residual_norm proc~try_newton_solve->proc~residual_norm proc~solve_small_linear_system solve_small_linear_system proc~try_newton_solve->proc~solve_small_linear_system

Called by

proc~~solve_nonlinear_system~~CalledByGraph proc~solve_nonlinear_system solve_nonlinear_system proc~try_solve_zhao_branch_a try_solve_zhao_branch_a proc~try_solve_zhao_branch_a->proc~solve_nonlinear_system proc~try_solve_zhao_branch_b try_solve_zhao_branch_b proc~try_solve_zhao_branch_b->proc~solve_nonlinear_system proc~try_solve_zhao_branch_c try_solve_zhao_branch_c proc~try_solve_zhao_branch_c->proc~solve_nonlinear_system proc~solve_zhao_branch_a solve_zhao_branch_a proc~solve_zhao_branch_a->proc~try_solve_zhao_branch_a proc~solve_zhao_branch_b solve_zhao_branch_b proc~solve_zhao_branch_b->proc~try_solve_zhao_branch_b proc~solve_zhao_branch_c solve_zhao_branch_c proc~solve_zhao_branch_c->proc~try_solve_zhao_branch_c proc~solve_zhao_unknowns solve_zhao_unknowns proc~solve_zhao_unknowns->proc~try_solve_zhao_branch_a proc~solve_zhao_unknowns->proc~try_solve_zhao_branch_b proc~solve_zhao_unknowns->proc~try_solve_zhao_branch_c proc~solve_zhao_unknowns->proc~solve_zhao_branch_a proc~solve_zhao_unknowns->proc~solve_zhao_branch_b proc~solve_zhao_unknowns->proc~solve_zhao_branch_c proc~resolve_sheath_injection_context resolve_sheath_injection_context proc~resolve_sheath_injection_context->proc~solve_zhao_unknowns proc~init_particle_batch_from_config init_particle_batch_from_config proc~init_particle_batch_from_config->proc~resolve_sheath_injection_context

Source Code

  subroutine solve_nonlinear_system(n, guesses, residual_fn, x_best, success)
    integer, intent(in) :: n
    real(dp), intent(in) :: guesses(:, :)
    procedure(nonlinear_residual) :: residual_fn
    real(dp), intent(out) :: x_best(n)
    logical, intent(out) :: success

    integer :: guess_idx
    real(dp) :: x_trial(n), best_norm, trial_norm
    logical :: trial_success

    success = .false.
    if (size(guesses, 1) /= n) error stop 'solve_nonlinear_system guess dimension mismatch.'
    x_best = guesses(:, 1)
    best_norm = huge(1.0d0)
    do guess_idx = 1, size(guesses, 2)
      call try_newton_solve(n, guesses(:, guess_idx), residual_fn, x_trial, trial_norm, trial_success)
      if (trial_norm < best_norm) then
        best_norm = trial_norm
        x_best = x_trial
      end if
      if (trial_success .and. trial_norm < nonlinear_tol) then
        success = .true.
        x_best = x_trial
        return
      end if
    end do

    success = best_norm < nonlinear_tol
  end subroutine solve_nonlinear_system