changeset 98:46598cedb076

New example based on the one avilable with bim pkg. *
author gedeone-octave <marcovass89@hotmail.it>
date Mon, 12 Aug 2013 20:42:39 +0200
parents dc9987325fea
children 9818e6302e7b
files example/Advection-Diffusion/3d_fem.m example/Advection-Diffusion/Advection_Diffusion.ufl example/Advection-Diffusion/Reaction.ufl
diffstat 3 files changed, 79 insertions(+), 0 deletions(-) [+]
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/example/Advection-Diffusion/3d_fem.m	Mon Aug 12 20:42:39 2013 +0200
@@ -0,0 +1,59 @@
+# the example is based on the one described here for the 
+# bim pkg:
+# http://wiki.octave.org/Bim_package#3D_Time_dependent_problem
+
+pkg load msh
+pkg load fpl
+pkg load fem-fenics
+
+fem_init_env ();
+problem = 'Advection_Diffusion';
+fem_create_rhs (problem);
+
+problem = 'Reaction';
+fem_create_rhs (problem);
+fem_create_fs (problem);
+
+x = linspace (0, 1, 40);
+y = z = linspace (0, 1, 20);
+msho = msh3m_structured_mesh (x, y, z, 1, 1:6);
+mshd = Mesh (msho);
+
+x = msho.p(1, :).';
+y = msho.p(2, :).';
+z = msho.p(3, :).';
+ 
+x0 = .2; sx = .1;
+y0 = .2; sy = .1;
+z0 = .8; sz = .1;
+ 
+u = exp (- ((x-x0)/(2*sx)) .^2 - ((y-y0)/(2*sy)) .^2 - ((z-z0)/(2*sz)) .^2);
+
+V  = Reaction_FunctionSpace (mshd);
+
+f = Expression ('phi', @(x,y,z) x + y -z);
+g = Expression ('mu', @(x,y,z) 0.01);
+
+a = Advection_Diffusion_BilinearForm (V, f, g);
+m = Reaction_BilinearForm (V);
+A = assemble (a);
+M = assemble (m);
+
+# Lumped reaction matrix
+ML = diag (sum (M, 1));
+
+function du = f (u, t, A, M)
+  du = - M \ (A * u);
+endfunction 
+ 
+time = linspace (0, 1, 30);
+lsode_options ("integration method", "adams");
+U = lsode (@(u, t) f(u, t, A, ML), u, time);
+ 
+for ii = 1:1:numel (time)
+  name = sprintf ("u_%3.3d", ii);
+  delete ([name ".vtu"]);
+  fpl_vtk_write_field (name, msho, {U(ii,:)', 'u'}, {}, 1);
+endfor
+
+
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/example/Advection-Diffusion/Advection_Diffusion.ufl	Mon Aug 12 20:42:39 2013 +0200
@@ -0,0 +1,12 @@
+# Copyright (C) 2005-2009 Anders Logg
+element = FiniteElement("Lagrange", tetrahedron, 1)
+
+
+u = TrialFunction(element)
+v = TestFunction(element)
+
+mu = Coefficient(element)
+phi = Coefficient(element)
+
+a = mu*inner(grad(u), grad(v))*dx - u * inner (grad(phi), grad (v))*dx
+L = 0
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/example/Advection-Diffusion/Reaction.ufl	Mon Aug 12 20:42:39 2013 +0200
@@ -0,0 +1,8 @@
+# Copyright (C) 2005-2009 Anders Logg
+element = FiniteElement("Lagrange", tetrahedron, 1)
+
+u = TrialFunction(element)
+v = TestFunction(element)
+
+a = u*v*dx
+L = 0