# Path: gusto/configuration.py
# def set_log_handler(comm):
#     def __init__(self, **kwargs):
#     def __setattr__(self, name, value):
#     def name(self):
# class Configuration(object):
# class TimesteppingParameters(Configuration):
# class OutputParameters(Configuration):
# class CompressibleParameters(Configuration):
# class ShallowWaterParameters(Configuration):
# class EadyParameters(Configuration):
# class CompressibleEadyParameters(CompressibleParameters, EadyParameters):
# class AdvectionOptions(Configuration, metaclass=ABCMeta):
# class EmbeddedDGOptions(AdvectionOptions):
# class RecoveredOptions(AdvectionOptions):
# class SUPGOptions(AdvectionOptions):
#     N = 0.01  # Brunt-Vaisala frequency (1/s)
#     T_0 = 273.15  # ref. temperature
#     H = None  # mean depth
#     H = None
#     L = None
#     N = sqrt(EadyParameters.Nsq)


from firedrake import (split, LinearVariationalProblem, Constant,
                       LinearVariationalSolver, TestFunctions, TrialFunctions,
                       TestFunction, TrialFunction, lhs, rhs, FacetNormal,
                       div, dx, jump, avg, dS_v, dS_h, ds_v, ds_t, ds_b, ds_tb, inner,
                       dot, grad, Function, VectorSpaceBasis, BrokenElement,
                       FunctionSpace, MixedFunctionSpace, DirichletBC)
from firedrake.petsc import flatten_parameters
from firedrake.parloops import par_loop, READ, INC
from pyop2.profiling import timed_function, timed_region
from gusto.configuration import logger, DEBUG
from gusto import thermodynamics
from abc import ABCMeta, abstractmethod, abstractproperty
        import numpy as np
        if mu is not None:
            eqn += dt*mu*inner(w, k)*inner(u, k)*dx

        aeqn = lhs(eqn)
        Leqn = rhs(eqn)

        # Function for the hybridized solutions
        self.urhol0 = Function(M)

        hybridized_prb = LinearVariationalProblem(aeqn, Leqn, self.urhol0)
        hybridized_solver = LinearVariationalSolver(hybridized_prb,
                                                    solver_parameters=self.solver_parameters,
                                                    options_prefix='ImplicitSolver')
        self.hybridized_solver = hybridized_solver

        # Project broken u into the HDiv space using facet averaging.
        # Weight function counting the dofs of the HDiv element:
        shapes = {"i": Vu.finat_element.space_dimension(),
                  "j": np.prod(Vu.shape, dtype=int)}
        weight_kernel = """
        for (int i=0; i<{i}; ++i)
            for (int j=0; j<{j}; ++j)
                w[i*{j} + j] += 1.0;
        """.format(**shapes)

        self._weight = Function(Vu)
        par_loop(weight_kernel, dx, {"w": (self._weight, INC)})

        # Averaging kernel
        self._average_kernel = """

# Path: gusto/configuration.py
# def set_log_handler(comm):
#     def __init__(self, **kwargs):
#     def __setattr__(self, name, value):
#     def name(self):
# class Configuration(object):
# class TimesteppingParameters(Configuration):
# class OutputParameters(Configuration):
# class CompressibleParameters(Configuration):
# class ShallowWaterParameters(Configuration):
# class EadyParameters(Configuration):
# class CompressibleEadyParameters(CompressibleParameters, EadyParameters):
# class AdvectionOptions(Configuration, metaclass=ABCMeta):
# class EmbeddedDGOptions(AdvectionOptions):
# class RecoveredOptions(AdvectionOptions):
# class SUPGOptions(AdvectionOptions):
#     N = 0.01  # Brunt-Vaisala frequency (1/s)
#     T_0 = 273.15  # ref. temperature
#     H = None  # mean depth
#     H = None
#     L = None
#     N = sqrt(EadyParameters.Nsq)


# Path: gusto/linear_solvers.py
from firedrake import (split, LinearVariationalProblem, Constant,
                       LinearVariationalSolver, TestFunctions, TrialFunctions,
                       TestFunction, TrialFunction, lhs, rhs, FacetNormal,
                       div, dx, jump, avg, dS_v, dS_h, ds_v, ds_t, ds_b, ds_tb, inner,
                       dot, grad, Function, VectorSpaceBasis, BrokenElement,
                       FunctionSpace, MixedFunctionSpace, DirichletBC)
from firedrake.petsc import flatten_parameters
from firedrake.parloops import par_loop, READ, INC
from pyop2.profiling import timed_function, timed_region
from gusto.configuration import logger, DEBUG
from gusto import thermodynamics
from abc import ABCMeta, abstractmethod, abstractproperty
        import numpy as np

        if mu is not None:
            eqn += dt*mu*inner(w, k)*inner(u, k)*dx

        aeqn = lhs(eqn)
        Leqn = rhs(eqn)

        # Function for the hybridized solutions
        self.urhol0 = Function(M)

        hybridized_prb = LinearVariationalProblem(aeqn, Leqn, self.urhol0)
        hybridized_solver = LinearVariationalSolver(hybridized_prb,
                                                    solver_parameters=self.solver_parameters,
                                                    options_prefix='ImplicitSolver')
        self.hybridized_solver = hybridized_solver

        # Project broken u into the HDiv space using facet averaging.
        # Weight function counting the dofs of the HDiv element:
        shapes = {"i": Vu.finat_element.space_dimension(),
                  "j": np.prod(Vu.shape, dtype=int)}
        weight_kernel = """
        for (int i=0; i<{i}; ++i)
            for (int j=0; j<{j}; ++j)
                w[i*{j} + j] += 1.0;
        """.format(**shapes)

        self._weight = Function(Vu)
        par_loop(weight_kernel, dx, {"w": (self._weight, INC)})

        # Averaging kernel
        self._average_kernel = """

        for (int i=0; i<{i}; ++i)