uinit Subroutine

   pure subroutine uinit(u, udot, rpar, ipar)
   !! This routine computes and loads the vector of initial values.
   !! The initial \(u\) values are given by the polynomial:
   !!
   !! $$ u = 16 x (1-x) y (1-y) $$
   !!
   !! The initial \(\dot{u}\) values are set to zero; ([[daskr]] corrects these during the  
   !! first time step).

      real(rk), intent(out) :: u(:)
      real(rk), intent(out) :: udot(:)
      real(rk), intent(in) :: rpar(:)
      integer, intent(in) :: ipar(:)

      integer :: i, ioff, j, k, m
      real(rk) :: dx, xj, yk

      m = ipar(4)
      dx = rpar(1)

      do k = 0, m + 1
         yk = k*dx
         ioff = (m + 2)*k
         do j = 0, m + 1
            xj = j*dx
            i = ioff + j + 1
            u(i) = 16*xj*(one - xj)*yk*(one - yk)
         end do
      end do

      udot = zero

   end subroutine uinit