Mercurial > fem-fenics-eugenio
annotate inst/example/Evolution/Evolution.m @ 164:7190852cf57f
Update the examples giving two input parameters to the BilinearForm.
author | gedeone-octave <marcovass89@hotmail.it> |
---|---|
date | Sat, 05 Oct 2013 10:40:06 +0100 |
parents | d03d998e64af |
children | e85390ed620d |
rev | line source |
---|---|
157
08dc5547f4d6
Update the examples to the new naming conventions.
gedeone-octave <marcovass89@hotmail.it>
parents:
95
diff
changeset
|
1 pkg load fem-fenics msh |
73 | 2 |
3 problem = 'Evolution'; | |
157
08dc5547f4d6
Update the examples to the new naming conventions.
gedeone-octave <marcovass89@hotmail.it>
parents:
95
diff
changeset
|
4 import_ufl_Problem (problem); |
73 | 5 |
157
08dc5547f4d6
Update the examples to the new naming conventions.
gedeone-octave <marcovass89@hotmail.it>
parents:
95
diff
changeset
|
6 mesh = Mesh (msh2m_structured_mesh (0:0.05:1, 0:0.05:1, 1, 1:4)); |
73 | 7 |
157
08dc5547f4d6
Update the examples to the new naming conventions.
gedeone-octave <marcovass89@hotmail.it>
parents:
95
diff
changeset
|
8 V = FunctionSpace(problem, mesh); |
73 | 9 |
89
d06b423169fa
The examples now follow the new naming convention.
gedeone-octave <marcovass89@hotmail.it>
parents:
73
diff
changeset
|
10 bc = DirichletBC (V, @(x,y) 1, 1:4); |
73 | 11 |
12 t = 0; | |
13 dt = 0.01; | |
14 T = 0.3; | |
15 | |
89
d06b423169fa
The examples now follow the new naming convention.
gedeone-octave <marcovass89@hotmail.it>
parents:
73
diff
changeset
|
16 k = Expression ('dt', @(x,y) dt); |
73 | 17 |
89
d06b423169fa
The examples now follow the new naming convention.
gedeone-octave <marcovass89@hotmail.it>
parents:
73
diff
changeset
|
18 u0 = Expression ('u0', @(x,y) 10*exp(-((x - 0.5)^2 + (y - 0.5)^2) / 0.02)); |
73 | 19 |
164
7190852cf57f
Update the examples giving two input parameters to the BilinearForm.
gedeone-octave <marcovass89@hotmail.it>
parents:
159
diff
changeset
|
20 a = BilinearForm (problem, V, V, k); |
157
08dc5547f4d6
Update the examples to the new naming conventions.
gedeone-octave <marcovass89@hotmail.it>
parents:
95
diff
changeset
|
21 L = LinearForm (problem, V, k, u0); |
94
4f688918f76a
Used new available functions assemble and assemble-systems.
gedeone-octave <marcovass89@hotmail.it>
parents:
89
diff
changeset
|
22 |
4f688918f76a
Used new available functions assemble and assemble-systems.
gedeone-octave <marcovass89@hotmail.it>
parents:
89
diff
changeset
|
23 A = assemble (a, bc); |
73 | 24 |
25 # solve the problem for each time step | |
26 # We need to update only the lhs | |
27 while t < T | |
28 t += dt; | |
29 | |
30 # we can pass u0 to the lhs indifferently as a fem_coeff or | |
31 # as a fem_func | |
157
08dc5547f4d6
Update the examples to the new naming conventions.
gedeone-octave <marcovass89@hotmail.it>
parents:
95
diff
changeset
|
32 L = LinearForm (problem, V, k, u0); |
94
4f688918f76a
Used new available functions assemble and assemble-systems.
gedeone-octave <marcovass89@hotmail.it>
parents:
89
diff
changeset
|
33 b = assemble (L, bc); |
73 | 34 |
35 u = A \ b; | |
36 | |
89
d06b423169fa
The examples now follow the new naming convention.
gedeone-octave <marcovass89@hotmail.it>
parents:
73
diff
changeset
|
37 u0 = Function ('u0', V, u); |
73 | 38 |
39 #press Q to make the plot continue | |
157
08dc5547f4d6
Update the examples to the new naming conventions.
gedeone-octave <marcovass89@hotmail.it>
parents:
95
diff
changeset
|
40 plot (u0); |
08dc5547f4d6
Update the examples to the new naming conventions.
gedeone-octave <marcovass89@hotmail.it>
parents:
95
diff
changeset
|
41 end |