Solve 1D transoport with moving boundary

README.md

1D Transport with Moving Boundary

Solve 1D transport with moving boundary, but with a fixed frame of reference to demonstrate how the level set method works

Install with Nix

Install Nix using these instuctions and this. Once you have a working installation of Nix use,

$ nix-shell --pure
[nix-shell]$ python main.py

to run it.

Install with Conda

Get Conda and then

$ conda env create -f environment.yml
$ conda activate transport-1D
(transport-1D)$ python main.py

DocumentforMoffat.pdf

environment.yml

name: transport-1D
channels:
  - conda-forge
dependencies:
  - python=3.7
  - fipy=3.3

main.py

from fipy import TransientTerm, DiffusionTerm, ImplicitSourceTerm
from fipy import CellVariable, Grid1D, FaceVariable, Viewer
import numpy as np


nx = 100
dx = 0.01
lx = 1.
k_value = 0.001
delta = 1.
coeff = 1.
conc0 = 1.
dt = 0.1

def level_set_function(position, time):
    return np.array(position) - time * delta * k_value


mesh = Grid1D(nx=nx, Lx=lx)

level_set_var = CellVariable(mesh=mesh)
level_set_var[:] = level_set_function(mesh.x, 0.)

conc = CellVariable(mesh=mesh, hasOld=True)
conc[:] = conc0
conc.constrain(1, mesh.facesRight)

mask = (level_set_var >= 0) & (level_set_var < dx)
x_center = mesh.x[mask][0]

diff_coeff = coeff * (mesh.faceCenters[0] > x_center)

implicit_source = mask * delta / dx / conc

eqn = TransientTerm() == DiffusionTerm(diff_coeff) - ImplicitSourceTerm(implicit_source)

time = 0.

viewer = Viewer(vars=conc)
viewer.plot()

for step in range(1000):
    conc.updateOld()
    sweep = 0
    res = 1.
    while sweep < 4:
        res = eqn.sweep(var=conc, dt=dt)
        sweep += 1
        print('sweep:', sweep, 'res:', res)
    print()
    print(step)
    print()
    if step % 10 == 0:
        viewer.plot()
    time += dt
    level_set_var[:] = level_set_function(mesh.x, time)

shell.nix

let
  pkgs = import (builtins.fetchGit {
    url = "https://github.com/NixOS/nixpkgs.git";
    rev = "92e1376cc3e37ee72f1417df05032556d39853c1";
    ref = "master";
  }) { };
  pypkgs = pkgs.python3Packages;
in
  pkgs.mkShell {
    nativeBuildInputs = with pypkgs; [ fipy jupyter pip ] ++ fipy.propagatedBuildInputs;
    postShellHook = ''
      jupyter nbextension install --py widgetsnbextension --user
      jupyter nbextension enable widgetsnbextension --user --py

      SOURCE_DATE_EPOCH=$(date +%s)
      export PYTHONUSERBASE=$PWD/.local
      export USER_SITE=`python -c "import site; print(site.USER_SITE)"`
      export PYTHONPATH=$PYTHONPATH:$USER_SITE
      export PATH=$PATH:$PYTHONUSERBASE/bin
    '';
  }