numerical_jacobian Subroutine

public subroutine numerical_jacobian(n, x, f0, residual_fn, jac)

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: n
real(kind=dp), intent(in) :: x(n)
real(kind=dp), intent(in) :: f0(n)
procedure(nonlinear_residual) :: residual_fn
real(kind=dp), intent(out) :: jac(n,n)

Called by

proc~~numerical_jacobian~~CalledByGraph proc~numerical_jacobian numerical_jacobian proc~try_newton_solve try_newton_solve proc~try_newton_solve->proc~numerical_jacobian proc~solve_nonlinear_system solve_nonlinear_system proc~solve_nonlinear_system->proc~try_newton_solve 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 numerical_jacobian(n, x, f0, residual_fn, jac)
    integer, intent(in) :: n
    real(dp), intent(in) :: x(n), f0(n)
    procedure(nonlinear_residual) :: residual_fn
    real(dp), intent(out) :: jac(n, n)

    integer :: j
    real(dp) :: h, xh(n), fh(n)

    do j = 1, n
      h = 1.0d-6*max(1.0d0, abs(x(j)))
      xh = x
      xh(j) = xh(j) + h
      call residual_fn(xh, fh)
      jac(:, j) = (fh - f0)/h
    end do
  end subroutine numerical_jacobian