# HG changeset patch # User gedeone-octave # Date 1376332959 -7200 # Node ID 46598cedb07661f20ba57af4816b378c99942f90 # Parent dc9987325fea149b128cbb6af7bf0b2b0ca8197d New example based on the one avilable with bim pkg. * diff -r dc9987325fea -r 46598cedb076 example/Advection-Diffusion/3d_fem.m --- /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 + + diff -r dc9987325fea -r 46598cedb076 example/Advection-Diffusion/Advection_Diffusion.ufl --- /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 diff -r dc9987325fea -r 46598cedb076 example/Advection-Diffusion/Reaction.ufl --- /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