not working

moving_domain/moving_domain_heat_equation_2.py

+import dolfin as df
+import numpy as np
+import matplotlib.pyplot as plt
+from matplotlib.widgets import Slider
+import progressbar
+
+# bother me only if there is an error
+df.set_log_level(df.LogLevel.ERROR)                             
+df.parameters['form_compiler']['optimize'] = True
+
+# parameters
+Nx = 2000
+L = 2*np.pi
+dt = 0.005
+T = 5
+D = 1.0
+k = 0.5
+times = np.arange(0, T, dt)
+
+# diffusion and advection
+def diffusion(func, testfunc, D):
+    return D * df.inner(func.dx(0), testfunc.dx(0))
+
+def advection(func, testfunc, vel):
+    return df.inner((vel * func).dx(0), testfunc)
+
+# mesh
+mesh = df.IntervalMesh(Nx, 0, L)
+x = mesh.coordinates()
+
+# function space
+conc_element = df.FiniteElement('P', mesh.ufl_cell(), 1)
+function_space = df.FunctionSpace(mesh, conc_element)
+
+# define velocity
+v = df.Constant(0.5)
+# v = df.interpolate(df.Expression('1.0', degree=1), function_space)
+
+# initial condition
+c0 = df.interpolate(df.Expression('1 + 0.2*cos(x[0])',degree=1), function_space)
+
+# define function, test function
+c = df.Function(function_space)
+tc = df.TestFunction(function_space)
+
+# weak form
+form = (df.inner((c - c0)/dt, tc) 
+        + diffusion(c, tc, D)
+        + advection(c, tc, v)
+        )*df.dx
+
+# define the arrays
+c_array = np.zeros((len(times), len(x)))
+c_array[0] = c0.compute_vertex_values(mesh)
+
+x_array = np.zeros((len(times), mesh.num_vertices()))
+x_array[0] = mesh.coordinates()[:,0]
+
+c_tot = np.zeros_like(times)
+c_tot[0] = df.assemble(c0 * df.dx(mesh))
+
+# time stepping
+for i in progressbar.progressbar(range(1,len(times))):
+    df.solve(form == 0, c)
+    c_tot[i] = df.assemble(c * df.dx(mesh))
+    c_array[i] = c.compute_vertex_values(mesh)
+    df.ALE.move(mesh, df.Expression('v*dt', v=v, dt=dt, degree=1))
+    x_array[i] = mesh.coordinates()[:,0]
+    c0.assign(c)
+
+
+# plotting
+
+# plot c(x,t) vs x for all t
+fig, axc = plt.subplots(1, 1, figsize=(8,6))
+axc.set_xlabel(r'$x$')
+axc.set_ylabel(r'$c(x,t)$')
+
+cplot, = axc.plot(x_array[0], c_array[0], 'r')
+axc.set_xlim(np.min(x_array)-2, np.max(x_array)+2)
+axc.set_ylim(np.min(c_array)-2, np.max(c_array)+2)
+
+def update(value):
+    ti = np.abs(times-value).argmin()
+    cplot.set_xdata(x_array[ti])
+
+    cplot.set_ydata(c_array[ti])
+    plt.draw()
+        
+sax = plt.axes([0.1, 0.92, 0.7, 0.02])
+slider = Slider(sax, r'$t/\tau$', min(times), max(times),
+        valinit=min(times), valfmt='%3.1f',
+        fc='#999999')
+slider.drawon = False
+slider.on_changed(update)
+
+# plot ctot vs time
+fig1, ax1 = plt.subplots(1, figsize=(8,6))
+fig1.subplots_adjust(left=0.1, bottom=0.1, top=0.9, right=0.9,
+        wspace=0.0, hspace=0.0)
+ax1.plot(times, c_tot)
+ax1.set_ylim([0,np.max(c_tot)+1])
+ax1.set_xlabel(r'$t$')
+ax1.set_title('Total concentration vs time')
+ax1.set_ylabel(r'$\int_{\Gamma} c$', color='#ff5b00')
+plt.show()
+

moving_domain_heat_equation_2.py

+import dolfin as df
+import numpy as np
+import matplotlib.pyplot as plt
+from matplotlib.widgets import Slider
+import progressbar
+
+# bother me only if there is an error
+df.set_log_level(df.LogLevel.ERROR)                             
+df.parameters['form_compiler']['optimize'] = True
+
+# parameters
+Nx = 2000
+L = 2*np.pi
+dt = 0.005
+T = 5
+D = 1.0
+k = 0.5
+times = np.arange(0, T, dt)
+
+# diffusion and advection
+def diffusion(func, testfunc, D):
+    return D * df.inner(func.dx(0), testfunc.dx(0))
+
+def advection(func, testfunc, vel):
+    return df.inner((vel * func).dx(0), testfunc)
+
+# mesh
+mesh = df.IntervalMesh(Nx, 0, L)
+x = mesh.coordinates()
+
+# function space
+conc_element = df.FiniteElement('P', mesh.ufl_cell(), 1)
+function_space = df.FunctionSpace(mesh, conc_element)
+
+# define velocity
+v = df.Constant(0.5)
+# v = df.interpolate(df.Expression('1.0', degree=1), function_space)
+
+# initial condition
+c0 = df.interpolate(df.Expression('1 + 0.2*cos(x[0])',degree=1), function_space)
+
+# define function, test function
+c = df.Function(function_space)
+tc = df.TestFunction(function_space)
+
+# weak form
+form = (df.inner((c - c0)/dt, tc) 
+        + diffusion(c, tc, D)
+        + advection(c, tc, v)
+        )*df.dx
+
+# define the arrays
+c_array = np.zeros((len(times), len(x)))
+c_array[0] = c0.compute_vertex_values(mesh)
+
+x_array = np.zeros((len(times), mesh.num_vertices()))
+x_array[0] = mesh.coordinates()[:,0]
+
+c_tot = np.zeros_like(times)
+c_tot[0] = df.assemble(c0 * df.dx(mesh))
+
+# time stepping
+for i in progressbar.progressbar(range(1,len(times))):
+    df.solve(form == 0, c)
+    c_tot[i] = df.assemble(c * df.dx(mesh))
+    c_array[i] = c.compute_vertex_values(mesh)
+    df.ALE.move(mesh, df.Expression('v*dt', v=v, dt=dt, degree=1))
+    x_array[i] = mesh.coordinates()[:,0]
+    c0.assign(c)
+
+
+# plotting
+
+# plot c(x,t) vs x for all t
+fig, axc = plt.subplots(1, 1, figsize=(8,6))
+axc.set_xlabel(r'$x$')
+axc.set_ylabel(r'$c(x,t)$')
+
+cplot, = axc.plot(x_array[0], c_array[0], 'r')
+axc.set_xlim(np.min(x_array)-2, np.max(x_array)+2)
+axc.set_ylim(np.min(c_array)-2, np.max(c_array)+2)
+
+def update(value):
+    ti = np.abs(times-value).argmin()
+    cplot.set_xdata(x_array[ti])
+
+    cplot.set_ydata(c_array[ti])
+    plt.draw()
+        
+sax = plt.axes([0.1, 0.92, 0.7, 0.02])
+slider = Slider(sax, r'$t/\tau$', min(times), max(times),
+        valinit=min(times), valfmt='%3.1f',
+        fc='#999999')
+slider.drawon = False
+slider.on_changed(update)
+
+# plot ctot vs time
+fig1, ax1 = plt.subplots(1, figsize=(8,6))
+fig1.subplots_adjust(left=0.1, bottom=0.1, top=0.9, right=0.9,
+        wspace=0.0, hspace=0.0)
+ax1.plot(times, c_tot)
+ax1.set_ylim([0,np.max(c_tot)+1])
+ax1.set_xlabel(r'$t$')
+ax1.set_title('Total concentration vs time')
+ax1.set_ylabel(r'$\int_{\Gamma} c$', color='#ff5b00')
+plt.show()
+