gels_rel Subroutine

private pure subroutine gels_rel(A, b, x, info)

Uses

  • proc~~gels_rel~~UsesGraph proc~gels_rel forsolver::gels_rel module~external_interfaces_solver external_interfaces_solver proc~gels_rel->module~external_interfaces_solver kinds kinds module~external_interfaces_solver->kinds

solves an overdetermined or underdetermined linear system using gels.

Arguments

Type IntentOptional Attributes Name
real(kind=rk), intent(in), dimension(:, :), contiguous :: A
real(kind=rk), intent(in), dimension(:), contiguous :: b
real(kind=rk), intent(out), dimension(max(1, size(A, 2))) :: x
integer, intent(out) :: info

Calls

proc~~gels_rel~~CallsGraph proc~gels_rel forsolver::gels_rel interface~gels external_interfaces_solver::gels proc~gels_rel->interface~gels

Called by

proc~~gels_rel~~CalledByGraph proc~gels_rel forsolver::gels_rel proc~solver_lin forsolver::solver_lin proc~solver_lin->proc~gels_rel interface~solve forsolver::solve interface~solve->proc~solver_lin proc~modified_newton_method_t1 forsolver::modified_newton_method_T1 proc~modified_newton_method_t1->interface~solve proc~modified_quasi_cs_newton_method_t1 forsolver::modified_quasi_cs_newton_method_T1 proc~modified_quasi_cs_newton_method_t1->interface~solve proc~modified_quasi_fd_newton_method_t1 forsolver::modified_quasi_fd_newton_method_T1 proc~modified_quasi_fd_newton_method_t1->interface~solve proc~newton_method_t1 forsolver::newton_method_T1 proc~newton_method_t1->interface~solve proc~quasi_cs_newton_method_t1 forsolver::quasi_cs_newton_method_T1 proc~quasi_cs_newton_method_t1->interface~solve proc~quasi_fd_newton_method_t1 forsolver::quasi_fd_newton_method_T1 proc~quasi_fd_newton_method_t1->interface~solve program~test_solver1 test_solver1 program~test_solver1->interface~solve program~test_solver15 test_solver15 program~test_solver15->interface~solve program~test_solver16 test_solver16 program~test_solver16->interface~solve program~test_solver17 test_solver17 program~test_solver17->interface~solve program~test_solver2 test_solver2 program~test_solver2->interface~solve 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 proc~newton_complex_step_rel_t1->proc~quasi_cs_newton_method_t1 proc~newton_rel_t1 forsolver::nlsolver%newton_rel_T1 proc~newton_rel_t1->proc~modified_newton_method_t1 proc~newton_rel_t1->proc~modified_quasi_fd_newton_method_t1 proc~newton_rel_t1->proc~newton_method_t1 proc~newton_rel_t1->proc~quasi_fd_newton_method_t1 none~solve forsolver::nlsolver%solve none~solve->proc~newton_complex_step_rel_t1 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

   pure subroutine gels_rel(A, b, x, info)

      use external_interfaces_solver

      ! inputs:
      real(rk), dimension(:, :), contiguous, intent(in) :: A    ! input matrix A
      real(rk), dimension(:),    contiguous, intent(in) :: b    ! right-hand side matrix b

      ! outputs:
      real(rk), dimension(max(1, size(A, 2))), intent(out) :: x    ! solution matrix x
      integer,                                 intent(out) :: info ! result info

      ! local variables
      character(1)                          :: trans
      integer                               :: m, n, lda, ldb, lwork, nrhs
      real(rk), allocatable                 :: work(:)
      real(rk)                              :: work1(1)
      real(rk), dimension(:,:), allocatable :: a_copy
      real(rk), dimension(:,:), allocatable :: b_copy

      ! 
      trans = 'n'

      ! get dimensions
      nrhs = 1 ! size(b, 2)
      m    = size(A, 1)
      n    = size(A, 2)
      lda  = max(1, m)
      ldb  = max(1, max(m, n))

      ! copy the input matrices
      a_copy = a
      allocate(b_copy(ldb, nrhs))
      b_copy(:, 1) = b

      ! calculate the optimal size of the work array
      call gels(trans, m, n, nrhs, a_copy, lda, b_copy, ldb, work1, -1, info)

      ! allocate work array
      lwork = nint(work1(1))
      allocate(work(lwork))

      ! call gels subroutine
      call gels(trans, m, n, nrhs, a_copy, lda, b_copy, ldb, work, lwork, info)

      ! copy the solution matrix
      if (info == 0) then
         if (trans == 'n') x = b_copy(1:n, 1) ! nrhs = 1
         if (trans == 't') x = b_copy(1:m, 1) ! nrhs = 1
      else
         error stop 'gels failed'
      end if

      ! deallocate workspace
      deallocate(work)
   end subroutine gels_rel