Remove remove_initial_w (#576)

examples/compressible_euler/mountain_hydrostatic.py

@@ -18,7 +18,7 @@
     TrapeziumRule, SUPGOptions, ZComponent, Perturbation,
     CompressibleParameters, HydrostaticCompressibleEulerEquations,
     CompressibleSolver, compressible_hydrostatic_balance, HydrostaticImbalance,
-    SpongeLayerParameters, MinKernel, MaxKernel, remove_initial_w, logger
+    SpongeLayerParameters, MinKernel, MaxKernel, logger
 )
 
 mountain_hydrostatic_defaults = {
@@ -243,8 +243,7 @@ def mountain_hydrostatic(
 
     theta0.assign(theta_b)
     rho0.assign(rho_b)
-    u0.project(as_vector([initial_wind, 0.0]))
-    remove_initial_w(u0)
+    u0.project(as_vector([initial_wind, 0.0]), bcs=eqns.bcs['u'])
 
     stepper.set_reference_profiles([('rho', rho_b), ('theta', theta_b)])
 

gusto/initialisation/hydrostatic_initialisation.py

@@ -12,9 +12,10 @@
 from gusto.recovery import Recoverer, BoundaryMethod
 
 
-__all__ = ["boussinesq_hydrostatic_balance",
-           "compressible_hydrostatic_balance", "remove_initial_w",
-           "saturated_hydrostatic_balance", "unsaturated_hydrostatic_balance"]
+__all__ = [
+    "boussinesq_hydrostatic_balance", "compressible_hydrostatic_balance",
+    "saturated_hydrostatic_balance", "unsaturated_hydrostatic_balance"
+]
 
 
 def boussinesq_hydrostatic_balance(equation, b0, p0, top=False, params=None):
@@ -219,23 +220,6 @@ def compressible_hydrostatic_balance(equation, theta0, rho0, exner0=None,
         rho0.interpolate(thermodynamics.rho(parameters, theta0, exner))
 
 
-def remove_initial_w(u):
-    """
-    Removes the vertical component of a velocity field.
-
-    Args:
-        u (:class:`Function`): the velocity field to be altered.
-    """
-    Vu = u.function_space()
-    Vv = FunctionSpace(Vu._ufl_domain, Vu.ufl_element()._elements[-1])
-    bc = DirichletBC(Vu[0], 0.0, "bottom")
-    bc.apply(u)
-    uv = Function(Vv).project(u)
-    ustar = Function(u.function_space()).project(uv)
-    uin = Function(u.function_space()).assign(u - ustar)
-    u.assign(uin)
-
-
 def saturated_hydrostatic_balance(equation, state_fields, theta_e, mr_t,
                                   exner0=None, top=False,
                                   exner_boundary=Constant(1.0),