DiffusiveFluxMethod_ESAtmo2D_t Subroutine

  subroutine DiffusiveFluxMethod_ESAtmo2D_t(this)
    !! Fill diffFlux%interior with the constant-coefficient Laplacian
    !! flux at every interior node:
    !!
    !!   F_d(rho)      = 0                              (no mass diffusion)
    !!   F_d(rho*v_i)  = -nu    * d(rho*v_i)/dx_d
    !!   F_d(rho*theta)= -kappa * d(rho*theta)/dx_d
    !!   F_d(Phi)      = 0
    !!
    !! solutionGradient%interior(i,j,iel,iVar,d) holds d(s_iVar)/dx_d.
    implicit none
    class(ESAtmo2D_t),intent(inout) :: this
    ! Local
    integer :: i,j,iel,d

    do concurrent(i=1:this%solution%N+1,j=1:this%solution%N+1, &
                  iel=1:this%mesh%nElem,d=1:2)
      this%diffFlux%interior(i,j,iel,1,d) = 0.0_prec
      this%diffFlux%interior(i,j,iel,2,d) = &
        -this%nu*this%solutionGradient%interior(i,j,iel,2,d)
      this%diffFlux%interior(i,j,iel,3,d) = &
        -this%nu*this%solutionGradient%interior(i,j,iel,3,d)
      this%diffFlux%interior(i,j,iel,4,d) = &
        -this%kappa*this%solutionGradient%interior(i,j,iel,4,d)
      ! Geopotential carries no diffusion.
      this%diffFlux%interior(i,j,iel,5,d) = 0.0_prec
    enddo

  endsubroutine DiffusiveFluxMethod_ESAtmo2D_t