viscosity for DG

[jiajiaecho]

Added the viscosity to the multi-species DG.

---

This change is 

[adityakpandare]

Thank you, @jiajiaecho. I have some comments/questions below. Also, there are several instances of (a) trailing whitespaces, and (b) incorrect indentations in function-scope or for-loops. Could you fix these please (see https://quinoacomputing.github.io/contributing.html#styling for styling guidelines in Quinoa). I have yet to review your Surface and Boundary integrals, will do that soon. Thanks.

@adityakpandare reviewed 9 files and all commit messages, and made 8 comments.
Reviewable status: 9 of 12 files reviewed, 7 unresolved discussions (waiting on @jiajiaecho).

---

src/PDE/Integrate/Basis.cpp line 1064 at r2 (raw file):

                 const std::array< std::vector<tk::real>, 3 >& dBdx )
// *****************************************************************************
//  Compute the state variables for the tetrahedron element

Let's update the comments here reflecting that this computes gradient of the state vars, and the arguments too for dBdx.

---

src/PDE/Integrate/Basis.cpp line 1099 at r2 (raw file):

                       + U( e, mark+5 ) * dBdx[i][5]
                       + U( e, mark+6 ) * dBdx[i][6]
                       + U( e, mark+4 ) * dBdx[i][7] 

Should this be mark+7 and +8 and +9 on the following lines?

---

src/PDE/Integrate/Basis.cpp line 1101 at r2 (raw file):

                       + U( e, mark+4 ) * dBdx[i][7] 
                       + U( e, mark+5 ) * dBdx[i][8]
                       + U( e, mark+6 ) * dBdx[i][9];

These loops are quite inefficient- you have SoA data-layout for state_grad, but dBdx is AoS. As a result the loops above will cache-miss every time the second index of state_grad is accessed. Let's attempt at making state_grad also AoS.

---

src/PDE/Integrate/Volume.cpp line 141 at r2 (raw file):

         for (ncomp_t c=0; c<ncomp; ++c){
           grad_all.push_back(state_P_grad[c]);

This will seg-fault. Looping variable should go from 0 to nprim for P, not ncomp.

---

src/PDE/MultiSpecies/DGMultiSpecies.hpp line 1045 at r2 (raw file):

    
// 2. Compute viscous stress tensor

If this is 2., what's 1.?

---

src/PDE/MultiSpecies/DGMultiSpecies.hpp line 1093 at r2 (raw file):

  tau[1][0] = tau[0][1];
  tau[2][0] = tau[0][2];
  tau[2][1] = tau[1][2];

Let's extract this tau computation into a separate function in namespace tk that computes this tau based on input arg dudx, because 'modifiedGradientViscousFlux()' also uses it.

---

src/PDE/MultiSpecies/DGMultiSpecies.hpp line 1108 at r2 (raw file):

      for (std::size_t i=0; i<3; ++i) {
        for (std::size_t j=0; j<3; ++j) {
        fl[idx_2][i] += u[j] * tau[i][j]+conduct*dTdx[i];

This indexing is incorrect. The flux is for energy equation- should go into multispecies::energyIdx(nspec, 0).