Mercurial > fem-fenics-eugenio
view example/Poisson/Poisson.m @ 158:17358a4eb648
Move the example to their own folder.
author | gedeone-octave <marcovass89@hotmail.it> |
---|---|
date | Thu, 12 Sep 2013 15:19:31 +0200 |
parents | example/Poisson.m@93a4ee13aa75 |
children |
line wrap: on
line source
pkg load fem-fenics msh import_ufl_Problem ('Poisson') # Create mesh and define function space x = y = linspace (0, 1, 33); mesh = Mesh(msh2m_structured_mesh (x, y, 1, 1:4)); V = FunctionSpace('Poisson', mesh); # Define boundary condition bc = DirichletBC(V, @(x, y) 0.0, [2;4]); f = Expression ('f', @(x,y) 10*exp(-((x - 0.5)^2 + (y - 0.5)^2) / 0.02)); g = Expression ('g', @(x,y) sin (5.0 * x)); a = BilinearForm ('Poisson', V); L = LinearForm ('Poisson', V, f, g); # Compute solution [A, b] = assemble_system (a, L, bc); sol = A \ b; u = Function ('u', V, sol); # Save solution in VTK format save (u, 'poisson') # Plot solution plot (u);