res Subroutine

public pure subroutine res(t, c, cprime, cj, delta, ires, rpar, ipar)

Uses

  • proc~~res~4~~UsesGraph proc~res~4 res module~web_par web_par proc~res~4->module~web_par module~daskr_kinds daskr_kinds module~web_par->module~daskr_kinds iso_fortran_env iso_fortran_env module~daskr_kinds->iso_fortran_env

This routine computes the residuals vector, using subroutine f for the right-hand sides.

Arguments

Type IntentOptional Attributes Name
real(kind=rk), intent(in) :: t
real(kind=rk), intent(in) :: c(*)
real(kind=rk), intent(in) :: cprime(*)
real(kind=rk), intent(in) :: cj
real(kind=rk), intent(out) :: delta(*)
integer, intent(in) :: ires
real(kind=rk), intent(inout) :: rpar(*)
integer, intent(in) :: ipar(*)

Calls

proc~~res~4~~CallsGraph proc~res~4 res proc~f~2 f proc~res~4->proc~f~2 proc~rates rates proc~f~2->proc~rates

Variables

Type Visibility Attributes Name Initial
integer, public :: i
integer, public :: ic0
integer, public :: ici
integer, public :: iyoff
integer, public :: jx
integer, public :: jy

Source Code

   pure subroutine res(t, c, cprime, cj, delta, ires, rpar, ipar)
   !! This routine computes the residuals vector, using subroutine [[f]] for the right-hand
   !! sides.
      use web_par, only: ns, np, mx, my, mxns
      real(rk), intent(in) :: t
      real(rk), intent(in) :: c(*)
      real(rk), intent(in) :: cprime(*)
      real(rk), intent(in) :: cj
      real(rk), intent(out) :: delta(*)
      integer, intent(in) :: ires
      real(rk), intent(inout) :: rpar(*)
      integer, intent(in) :: ipar(*)

      integer :: i, ic0, ici, iyoff, jx, jy

      call f(t, c, delta, rpar)

      do jy = 1, my
         iyoff = mxns*(jy - 1)
         do jx = 1, mx
            ic0 = iyoff + ns*(jx - 1)
            do i = 1, ns
               ici = ic0 + i
               if (i > np) then
                  delta(ici) = -delta(ici)
               else
                  delta(ici) = cprime(ici) - delta(ici)
               end if
            end do
         end do
      end do

   end subroutine res