Jacobian Subroutine

    subroutine Jacobian ( n, vx, r, dr, x, rx )
    !! Calculate the Jacobian for transforming derivatives from the
    !! reference element to the physical element.
      integer(ip), intent(in) :: n
      !! The order of the element.
      real(wp), dimension(2), intent(in) :: vx
      !! A 1d real array of size 2 that contains the physical coordinates at
      !! the boundaries of the element.`
      real(wp), dimension(n+1), intent(in) :: r
      !! The node locations within the reference element, \(r_i\).
      real(wp), dimension(:,:), intent(in) :: dr
      !! The derivative matrix for the reference element, \(\mathcal{D}_r\).
      real(wp), dimension(n+1), intent(out) :: x
      !! On output contains the physical coordinates for this physical element,
      !! \(x_i\).
      real(wp), dimension(n+1), intent(out) :: rx
      !! On output contains \(\frac{d r_i}{d x_j}.
      real(wp), dimension(n+1) :: xr
      integer(ip) :: i,j
      
      do i = 1, n+1
        x(i) = vx(1)+0.5_wp*(1.0_wp+r(i))*(vx(2)-vx(1))
      end do
      xr = 0.0_wp
      do i = 1, n+1
        do j = 1, n+1
          xr(i) = xr(i) + dr(i,j) * x(j)
        end do
      end do
      rx = 1.0_wp / xr
    end subroutine Jacobian