mxdprb Subroutine

public pure subroutine mxdprb(n, a, x, job)

solution of a system of linear equations with a dense symmetric positive definite matrix a using the factorization a=trans(r)*r.

Method

back substitution

Arguments

Type IntentOptional Attributes Name
integer, intent(in) :: n

order of the matrix a.
real(kind=wp), intent(in) :: a(*)

a(n(n+1)/2) factorization a=trans(r)r.
real(kind=wp), intent(inout) :: x(*)

x(n) on input the right hand side of a system of linear equations. on output the solution of a system of linear equations.
integer, intent(in) :: job

option

  • if job=0 then x:=a*(-1)x.
  • if job>0 then x:=trans(r)*(-1)x.
  • if job<0 then x:=r*(-1)x.

Called by

proc~~mxdprb~~CalledByGraph proc~mxdprb mxdprb proc~dual_range_space_quad_prog psqp_class%dual_range_space_quad_prog proc~dual_range_space_quad_prog->proc~mxdprb proc~update_tri_decomp_general update_tri_decomp_general proc~dual_range_space_quad_prog->proc~update_tri_decomp_general proc~update_tri_decomp_general->proc~mxdprb proc~psqp psqp_class%psqp proc~psqp->proc~dual_range_space_quad_prog proc~psqpn psqp_class%psqpn proc~psqpn->proc~psqp

Source Code

      pure subroutine mxdprb(n,a,x,job)

      integer,intent(in) :: job  !! option
                                 !!
                                 !! * if job=0 then x:=a**(-1)*x.
                                 !! * if job>0 then x:=trans(r)**(-1)*x.
                                 !! * if job<0 then x:=r**(-1)*x.
      integer,intent(in) :: n  !! order of the matrix a.
      real(wp),intent(in) :: a(*)  !! a(n*(n+1)/2) factorization a=trans(r)*r.
      real(wp),intent(inout) :: x(*)  !! x(n)  on input the right hand side of a system of linear
                                      !! equations. on output the solution of a system of linear
                                      !! equations.

      integer :: i , ii , ij , j

      if ( job>=0 ) then
         ! phase 1 : x:=trans(r)**(-1)*x
         ij = 0
         do i = 1 , n
            do j = 1 , i - 1
               ij = ij + 1
               x(i) = x(i) - a(ij)*x(j)
            end do
            ij = ij + 1
            x(i) = x(i)/a(ij)
         end do
      endif
      if ( job<=0 ) then
         ! phase 2 : x:=r**(-1)*x
         ii = n*(n+1)/2
         do i = n , 1 , -1
            ij = ii
            do j = i + 1 , n
               ij = ij + j - 1
               x(i) = x(i) - a(ij)*x(j)
            end do
            x(i) = x(i)/a(ii)
            ii = ii - i
         end do
      endif
      end subroutine mxdprb