pjac Subroutine

public subroutine pjac(res, ierr_res, neq, t, c, cdot, rewt, savr, wk, h, cj, rwp, iwp, ierr_pjac, rpar, ipar)

Uses

  • proc~~pjac~~UsesGraph proc~pjac pjac module~daskr daskr proc~pjac->module~daskr module~daskr_rbdpre daskr_rbdpre proc~pjac->module~daskr_rbdpre module~daskr_rbgpre daskr_rbgpre proc~pjac->module~daskr_rbgpre module~daskr_kinds daskr_kinds module~daskr->module~daskr_kinds module~daskr_rbdpre->module~daskr_kinds module~daskr_rbgpre->module~daskr_kinds iso_fortran_env iso_fortran_env module~daskr_kinds->iso_fortran_env

This routine interfaces to subroutines pjac_rbdpre or pjac_rbgpre, depending on the flag jbg = ipar(2), to generate and preprocess the block-diagonal Jacobian corresponding to the reaction term.

  • If jbg == 0, we call pjac_rbdpre, without block-grouping.
  • If jbg == 1, we call pjac_rbgpre, with block-grouping.

Array rpar, containing the current vector, is passed to pjac_rbdpre and pjac_rbgpre as argument r0, consistent with the loading of rpar in subroutine f. The procedure name rates is passed as the name of the routine which computes the individual blocks of .

Arguments

Type IntentOptional Attributes Name
procedure(res_t) :: res
integer, intent(out) :: ierr_res
integer, intent(in) :: neq
real(kind=rk), intent(in) :: t
real(kind=rk), intent(inout) :: c(*)
real(kind=rk), intent(in) :: cdot(*)
real(kind=rk), intent(in) :: rewt(*)
real(kind=rk), intent(inout) :: savr(*)
real(kind=rk), intent(inout) :: wk(*)
real(kind=rk), intent(in) :: h
real(kind=rk), intent(in) :: cj
real(kind=rk), intent(inout) :: rwp(*)
integer, intent(inout) :: iwp(*)
integer, intent(inout) :: ierr_pjac
real(kind=rk), intent(inout) :: rpar(*)
integer, intent(inout) :: ipar(*)

Calls

proc~~pjac~~CallsGraph proc~pjac pjac proc~pjac_rbdpre pjac_rbdpre proc~pjac->proc~pjac_rbdpre proc~pjac_rbgpre pjac_rbgpre proc~pjac->proc~pjac_rbgpre proc~dgefa dgefa proc~pjac_rbdpre->proc~dgefa proc~pjac_rbgpre->proc~dgefa interface~idamax idamax proc~dgefa->interface~idamax

Source Code

   subroutine pjac( &
      res, ierr_res, neq, t, c, cdot, rewt, savr, wk, h, cj, rwp, iwp, ierr_pjac, rpar, ipar)
   !! This routine interfaces to subroutines [[pjac_rbdpre]] or [[pjac_rbgpre]], depending
   !! on the flag `jbg = ipar(2)`, to generate and preprocess the block-diagonal Jacobian
   !! corresponding to the reaction term.
   !!
   !! * If `jbg == 0`, we call [[pjac_rbdpre]], without block-grouping.
   !! * If `jbg == 1`, we call [[pjac_rbgpre]], with block-grouping.
   !!
   !! Array `rpar`, containing the current \(R\) vector, is passed to [[pjac_rbdpre]] and
   !! [[pjac_rbgpre]] as argument `r0`, consistent with the loading of `rpar` in subroutine
   !! `f`. The procedure name `rates` is passed as the name of the routine which computes
   !! the individual blocks of \(R\).

      use daskr_rbdpre, only: pjac_rbdpre
      use daskr_rbgpre, only: pjac_rbgpre
      use daskr, only: res_t

      procedure(res_t) :: res
      integer, intent(out) :: ierr_res
      integer, intent(in) :: neq
      real(rk), intent(in) :: t
      real(rk), intent(inout) :: c(*)
      real(rk), intent(in) :: cdot(*)
      real(rk), intent(in) :: rewt(*)
      real(rk), intent(inout) :: savr(*)
      real(rk), intent(inout) :: wk(*)
      real(rk), intent(in) :: h
      real(rk), intent(in) :: cj
      real(rk), intent(inout) :: rwp(*)
      integer, intent(inout) :: iwp(*)
      integer, intent(inout) :: ierr_pjac
      real(rk), intent(inout) :: rpar(*)
      integer, intent(inout) :: ipar(*)

      integer :: jbg

      jbg = ipar(2)
      if (jbg == 0) then
         call pjac_rbdpre(t, c, rpar, rates, wk, rewt, cj, rwp, iwp, ierr_pjac)
      else
         call pjac_rbgpre(t, c, rpar, rates, wk, rewt, cj, rwp, iwp, ierr_pjac)
      end if

   end subroutine pjac