add support for interior face integrals

[samuelpmishLLNL]

This PR aims to implement a feature requested by @tupek2 to be able to integrate over interior faces. The main challenge is in handling DG/L2 spaces for this type of integral, as they have separate values on each side.

---

The the new kind of integral behaves similarly to the existing boundary integral, except that for L2 spaces it represents:

where the domain of integration $\Omega$ represents the specified collection of interior faces, $\phi_1, \phi_2$ are the basis functions from side 1 and side 2, and $u_1, u_2$ are the values of the DG field evaluated from side 1 and side 2. For example, consider a mesh made of two adjacent quadrilateral elements, with a L2<1, 1> field defined on it. The red disks represent the nodes of the DG field (with the positions of the nodes shifted inward slightly to exaggerate which side of the interior face they belong to), and the blue edge represents the domain $\Omega$.

In this example, the basis functions that appear in the integrand are:

where $N_0, N_1$ are the usual linear interpolation shape functions in 1D.

---

Here is an example of a residual calculation involving a vector-valued DG space and interior faces:

  ...

  constexpr int p = 2; // polynomial order
  using test_space  = L2<p, dim>;
  using trial_space = L2<p, dim>;

  constexpr int DERIVATIVE = 1;

  Domain interior_faces = InteriorFaces(*mesh);

  // Construct the new functional object using the specified test and trial spaces
  Functional<test_space(trial_space)> residual(&fespace, {&fespace});

  // register a new integral calculation with this residual object
  residual.AddInteriorFaceIntegral(
      Dimension<dim-1>{}, DependsOn<0>{},
      [=](double /*t*/, auto X, auto velocity) {

        // compute the (unit) surface normal
        auto dX_dxi = get<DERIVATIVE>(X);
        auto n = normalize(cross(dX_dxi));

        // extract the velocity values from each side of the interface
        // note: the orientation convention is such that the normal 
        //       computed as above will point from from side 1->2
        auto [u_1, u_2] = velocity; 

        // carry out whatever calculations are necessary to compute f_1, f_2
        auto a = dot(u_2 - u_1, n);
        auto f_1 = u_1 * a;
        auto f_2 = u_2 * a;

        // note: the values returned in the tuple are integrated 
        // against the basis functions from their respective sides
        return serac::tuple{f_1, f_2};

      }, interior_faces);

  ...

  // evaluate the residual with nodal values specified by U
  mfem::Vector U = ...; 
  mfem::Vector r = residual(t, U);

---

Here's a table showing which argument types are passed to qfunction and which weight functions are applied to the output of qfunctions, depending on the type of integral and trial/test space (respectively):

[samuelpmishLLNL]

/style

[samuelpmishLLNL]

I've temporarily disabled these serac "drivers" since:

• to my knowledge, no one uses them
• we don't really have an input-file interface for a lot of things (Domain included)
• it's caused confusion for new users (understandably, but incorrectly) thinking that this is the primary way to use serac

Keeping them enabled would add extra burden to this PR which I don't feel is justified in the near future.

[white238]

There is a user that is currently using it for testing the memory mapped file library.

[samuelpmishLLNL]

This is ready for review, so please take a look! My apologies to the reviewers for the number of files touched.