From aeeee3e4edcc0db2a7376034536d2df102d91fce Mon Sep 17 00:00:00 2001 From: Daniel Shapero Date: Mon, 8 Jan 2024 12:23:38 -0800 Subject: [PATCH] Remove dualform library --- .github/workflows/compile.yml | 2 - demos/dualform/pyproject.toml | 12 ----- demos/dualform/src/dualform/__init__.py | 1 - demos/dualform/src/dualform/ice_shelf.py | 42 ------------------ demos/dualform/src/dualform/ice_stream.py | 53 ----------------------- demos/dualform/src/dualform/viscosity.py | 21 --------- 6 files changed, 131 deletions(-) delete mode 100644 demos/dualform/pyproject.toml delete mode 100644 demos/dualform/src/dualform/__init__.py delete mode 100644 demos/dualform/src/dualform/ice_shelf.py delete mode 100644 demos/dualform/src/dualform/ice_stream.py delete mode 100644 demos/dualform/src/dualform/viscosity.py diff --git a/.github/workflows/compile.yml b/.github/workflows/compile.yml index 6394819..c3c03d3 100644 --- a/.github/workflows/compile.yml +++ b/.github/workflows/compile.yml @@ -23,8 +23,6 @@ jobs: pip install git+https://github.com/icepack/icepack2.git - name: Check out git repository uses: actions/checkout@v3 - - name: Install Python package - run: pip install --editable demos/dualform - name: Compile TeX document env: EARTHDATA_USERNAME: ${{ secrets.EARTHDATA_USERNAME }} diff --git a/demos/dualform/pyproject.toml b/demos/dualform/pyproject.toml deleted file mode 100644 index 4c3ad07..0000000 --- a/demos/dualform/pyproject.toml +++ /dev/null @@ -1,12 +0,0 @@ -[project] -name = "dualform" -version = "0.0.1" -authors = [ - { name="Daniel Shapero", email="shapero@uw.edu" }, - { name="Gonzalo Gonzalez de Diego", email="Gonzalo.GonzalezDeDiego@maths.ox.ac.uk" }, -] -description = "Dual forms of momentum conservation equations for ice flow" -requires-python = ">=3.8" - -[project.urls] -"Homepage" = "https://github.com/icepack/dual-problems" diff --git a/demos/dualform/src/dualform/__init__.py b/demos/dualform/src/dualform/__init__.py deleted file mode 100644 index 6690f5c..0000000 --- a/demos/dualform/src/dualform/__init__.py +++ /dev/null @@ -1 +0,0 @@ -from . import ice_shelf, ice_stream diff --git a/demos/dualform/src/dualform/ice_shelf.py b/demos/dualform/src/dualform/ice_shelf.py deleted file mode 100644 index 0193dc3..0000000 --- a/demos/dualform/src/dualform/ice_shelf.py +++ /dev/null @@ -1,42 +0,0 @@ -import ufl -import firedrake -from firedrake import Constant, inner, sym, grad, dx, ds, dS -import icepack -from .viscosity import power as viscous_power - - -__all__ = ["viscous_power", "boundary", "constraint"] - - -# Physical constants -ρ_I = Constant(icepack.constants.ice_density) -ρ_W = Constant(icepack.constants.water_density) -ρ = ρ_I * (1 - ρ_I / ρ_W) -g = Constant(icepack.constants.gravity) - - -def boundary(**kwargs): - # Get all the dynamical fields and boundary conditions - u, M, h = map(kwargs.get, ("velocity", "membrane_stress", "thickness")) - outflow_ids = kwargs["outflow_ids"] - - # Get the unit outward normal vector to the terminus - mesh = u.ufl_domain() - ν = firedrake.FacetNormal(mesh) - - return 0.5 * ρ * g * h**2 * inner(u, ν) * ds(outflow_ids) - - -def constraint(**kwargs): - u, M, h = map(kwargs.get, ("velocity", "membrane_stress", "thickness")) - ε = sym(grad(u)) - return (-h * inner(M, ε) - inner(0.5 * ρ * g * grad(h**2), u)) * dx - - -def constraint_edges(**kwargs): - u, h = map(kwargs.get, ("velocity", "thickness")) - mesh = ufl.domain.extract_unique_domain(u) - ν = firedrake.FacetNormal(mesh) - u_ν = inner(u, ν) - return 0.5 * ρ * g * (h("+")**2 * u_ν("+") + h("-")**2 * u_ν("-")) * dS - diff --git a/demos/dualform/src/dualform/ice_stream.py b/demos/dualform/src/dualform/ice_stream.py deleted file mode 100644 index ecf4f03..0000000 --- a/demos/dualform/src/dualform/ice_stream.py +++ /dev/null @@ -1,53 +0,0 @@ -import firedrake -from firedrake import Constant, inner, sym, grad, dx, ds -import icepack -from .viscosity import power as viscous_power - - -__all__ = ["viscous_power", "friction_power", "boundary", "constraint"] - - -# Physical constants -ρ_I = Constant(icepack.constants.ice_density) -ρ_W = Constant(icepack.constants.water_density) -g = Constant(icepack.constants.gravity) - - -def boundary(**kwargs): - # Get all the dynamical fields and boundary conditions - u, h, s = map(kwargs.get, ("velocity", "thickness", "surface")) - outflow_ids = kwargs["outflow_ids"] - - # Get the unit outward normal vector to the terminus - mesh = u.ufl_domain() - ν = firedrake.FacetNormal(mesh) - - # Compute the forces per unit length at the terminus from the glacier - # and from the ocean (assuming that sea level is at z = 0) - f_I = 0.5 * ρ_I * g * h**2 - d = firedrake.min_value(0, s - h) - f_W = 0.5 * ρ_W * g * d**2 - - return (f_I - f_W) * inner(u, ν) * ds(outflow_ids) - - -def friction_power(**kwargs): - τ = kwargs["basal_stress"] - parameter_names = ( - "friction_yield_speed", "friction_yield_stress", "sliding_law_exponent" - ) - u_c, τ_c, m = map(kwargs.get, parameter_names) - K = u_c / τ_c ** m - τ_2 = inner(τ, τ) - τ_m = τ_2 if float(m) == 1 else τ_2 ** ((m + 1) / 2) - return K / (m + 1) * τ_m * dx - - -def constraint(**kwargs): - field_names = ( - "velocity", "membrane_stress", "basal_stress", "thickness", "surface" - ) - u, M, τ, h, s = map(kwargs.get, field_names) - f = kwargs.get("floating", Constant(1.0)) - ε = sym(grad(u)) - return (-h * inner(M, ε) + inner(f * τ - ρ_I * g * h * grad(s), u)) * dx diff --git a/demos/dualform/src/dualform/viscosity.py b/demos/dualform/src/dualform/viscosity.py deleted file mode 100644 index 8d0371a..0000000 --- a/demos/dualform/src/dualform/viscosity.py +++ /dev/null @@ -1,21 +0,0 @@ -from firedrake import inner, tr, dx -import icepack - - -def power(**kwargs): - # Get all the dynamical fields - u, M, h = map(kwargs.get, ("velocity", "membrane_stress", "thickness")) - - # Get the parameters for the constitutive relation - parameter_names = ( - "viscous_yield_strain", "viscous_yield_stress", "flow_law_exponent" - ) - ε, τ, n = map(kwargs.get, parameter_names) - A = ε / τ**n - - mesh = u.ufl_domain() - d = mesh.geometric_dimension() - - M_2 = (inner(M, M) - tr(M) ** 2 / (d + 1)) / 2 - M_n = M_2 if float(n) == 1 else M_2 ** ((n + 1) / 2) - return 2 * h * A / (n + 1) * M_n * dx