psol_banpre Subroutine

public pure subroutine psol_banpre(neq, t, y, ydot, savr, wk, cj, wght, rwp, iwp, b, epslin, ierr, rpar, ipar)

Uses

  • proc~~psol_banpre~~UsesGraph proc~psol_banpre psol_banpre module~linpack linpack proc~psol_banpre->module~linpack module~daskr_kinds daskr_kinds module~linpack->module~daskr_kinds iso_fortran_env iso_fortran_env module~daskr_kinds->iso_fortran_env

This routine solves the linear system for the banded preconditioner , given a vector , using the LU decomposition produced by jac_banpre. The solution is carried out by dgbsl from LINPACK.

Arguments

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

Problem size.

real(kind=rk), intent(in) :: t

Current independent variable (not used).

real(kind=rk), intent(in) :: y(*)

Current dependent variables (not used).

real(kind=rk), intent(in) :: ydot(*)

Current derivatives of dependent variables (not used).

real(kind=rk), intent(in) :: savr(*)

Current residual evaluated at (t, y, ydot) (not used).

real(kind=rk), intent(in) :: wk(*)

Real work space available to this subroutine (not used).

real(kind=rk), intent(in) :: cj

Scalar used in forming the system Jacobian (not used).

real(kind=rk), intent(in) :: wght(*)

Error weights for computing norms (not used).

real(kind=rk), intent(inout) :: rwp(*)

Real work array containing the LU decomposition of .

integer, intent(inout) :: iwp(*)

Integer array containing matrix pivot information.

real(kind=rk), intent(inout) :: b(*)

Right-hand side vector on input; solution on output.

real(kind=rk), intent(in) :: epslin

Tolerance for linear system (not used).

integer, intent(inout) :: ierr

Error flag.

real(kind=rk), intent(inout) :: rpar(*)

Real array used for communication between the calling program and user routines (not used).

integer, intent(inout) :: ipar(*)

Integer array used for communication between the calling program and user routines (not used).


Calls

proc~~psol_banpre~~CallsGraph proc~psol_banpre psol_banpre proc~dgbsl dgbsl proc~psol_banpre->proc~dgbsl interface~daxpy daxpy proc~dgbsl->interface~daxpy interface~ddot ddot proc~dgbsl->interface~ddot

Source Code

   pure subroutine psol_banpre( &
      neq, t, y, ydot, savr, wk, cj, wght, rwp, iwp, b, epslin, ierr, rpar, ipar)
   !! This routine solves the linear system \(P x = b\) for the banded preconditioner \(P\),
   !! given a vector \(b\), using the LU decomposition produced by [[jac_banpre]]. The solution
   !! is carried out by [[dgbsl]] from LINPACK.

      use linpack, only: dgbsl

      integer, intent(in) :: neq
        !! Problem size.
      real(rk), intent(in) :: t
        !! Current independent variable (not used).
      real(rk), intent(in) :: y(*)
        !! Current dependent variables (not used).
      real(rk), intent(in) :: ydot(*)
        !! Current derivatives of dependent variables (not used).
      real(rk), intent(in) :: savr(*)
        !! Current residual evaluated at `(t, y, ydot)` (not used).
      real(rk), intent(in) :: wk(*)
        !! Real work space available to this subroutine (not used).
      real(rk), intent(in) :: cj
        !! Scalar used in forming the system Jacobian (not used).
      real(rk), intent(in) :: wght(*)
        !! Error weights for computing norms (not used).
      real(rk), intent(inout) :: rwp(*)
        !! Real work array containing the LU decomposition of \(P\).
      integer, intent(inout) :: iwp(*)
        !! Integer array containing matrix pivot information.
      real(rk), intent(inout) :: b(*)
        !! Right-hand side vector on input; solution on output.
      real(rk), intent(in) :: epslin
        !! Tolerance for linear system (not used).
      integer, intent(inout) :: ierr
        !! Error flag.
      real(rk), intent(inout) :: rpar(*)
        !! Real array used for communication between the calling program and user routines
        !! (not used).
      integer, intent(inout) :: ipar(*)
        !! Integer array used for communication between the calling program and user routines 
        !! (not used).
      
      integer :: meband, ml, mu

      ml = ipar(1)
      mu = ipar(2)
      meband = 2*ml + mu + 1
      
      ierr = 0
      call dgbsl(rwp, meband, neq, ml, mu, iwp, b, 0)

   end subroutine psol_banpre