~ubuntu-branches/ubuntu/maverick/dolfin/maverick

« back to all changes in this revision

Viewing changes to demo/la/trilinos/python/demo.py

  • Committer: Bazaar Package Importer
  • Author(s): Johannes Ring
  • Date: 2009-10-12 14:13:18 UTC
  • mfrom: (1.1.1 upstream)
  • Revision ID: james.westby@ubuntu.com-20091012141318-hkbxl0sq555vqv7d
Tags: 0.9.4-1
* New upstream release. This version cleans up the design of the
  function class by adding a new abstraction for user-defined
  functions called Expression. A number of minor bugfixes and
  improvements have also been made.
* debian/watch: Update for new flat directory structure.
* Update debian/copyright.
* debian/rules: Use explicit paths to PETSc 3.0.0 and SLEPc 3.0.0.

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
1
""" This demo implements a Poisson equations solver
2
2
based on the demo "dolfin/demo/pde/poisson/python/demo.py"
3
 
in Dolfin using Epetra matrices, the AztecOO CG solver and ML 
4
 
AMG preconditioner 
 
3
in Dolfin using Epetra matrices, the AztecOO CG solver and ML
 
4
AMG preconditioner
5
5
"""
6
6
 
7
7
__author__ = "Kent-Andre Mardal (kent-and@simula.no)"
10
10
 
11
11
# Test for Trilinos:
12
12
try:
13
 
    from PyTrilinos import Epetra, AztecOO, TriUtils, ML 
 
13
    from PyTrilinos import Epetra, AztecOO, TriUtils, ML
14
14
except:
15
15
    print "You Need to have PyTrilinos with Epetra, AztecOO, TriUtils and ML installed for this demo to run",
16
16
    print "Exiting."
18
18
 
19
19
from dolfin import *
20
20
 
21
 
if not has_linear_algebra_backend("Epetra"):
 
21
if not has_la_backend("Epetra"):
22
22
    print "*** Warning: Dolfin is not compiled with Trilinos linear algebra backend"
23
23
    print "Exiting."
24
24
    exit()
25
25
 
26
 
dolfin_set("linear algebra backend", "Epetra")
 
26
parameters["linear_algebra_backend"] = "Epetra"
27
27
 
28
28
# Create mesh and finite element
29
29
mesh = UnitSquare(20,20)
37
37
# Define variational problem
38
38
v = TestFunction(V)
39
39
u = TrialFunction(V)
40
 
f = Function(V,"500.0 * exp(-(pow(x[0] - 0.5, 2) + pow(x[1] - 0.5, 2)) / 0.02)")
41
 
 
 
40
f = Expression("500.0 * exp(-(pow(x[0] - 0.5, 2) + pow(x[1] - 0.5, 2)) / 0.02)", V= V)
42
41
 
43
42
a = dot(grad(v), grad(u))*dx
44
 
L = v*f*dx 
 
43
L = v*f*dx
45
44
 
46
45
 
47
46
# Define boundary condition
49
48
bc = DirichletBC(V, u0, DirichletBoundary())
50
49
 
51
50
# Create linear system
52
 
A, b = assemble_system(a, L, bc) 
 
51
A, b = assemble_system(a, L, bc)
53
52
 
54
 
# Solution   
 
53
# Solution
55
54
U = Function(V)
56
55
 
57
 
# Fetch underlying epetra objects 
58
 
A_epetra = down_cast(A).mat() 
59
 
b_epetra = down_cast(b).vec() 
60
 
x_epetra = down_cast(U.vector()).vec() 
 
56
# Fetch underlying epetra objects
 
57
A_epetra = down_cast(A).mat()
 
58
b_epetra = down_cast(b).vec()
 
59
x_epetra = down_cast(U.vector()).vec()
61
60
 
62
61
# Sets up the parameters for ML using a python dictionary
63
 
MLList = {"max levels"        : 3, 
 
62
MLList = {"max levels"        : 3,
64
63
          "output"            : 10,
65
64
          "smoother: type"    : "ML symmetric Gauss-Seidel",
66
65
          "aggregation: type" : "Uncoupled",
67
66
          "ML validate parameter list" : False
68
67
}
69
68
 
70
 
# Create the preconditioner 
 
69
# Create the preconditioner
71
70
Prec = ML.MultiLevelPreconditioner(A_epetra, False)
72
71
Prec.SetParameterList(MLList)
73
72
Prec.ComputePreconditioner()
74
73
 
75
 
# Create solver and solve system 
 
74
# Create solver and solve system
76
75
Solver = AztecOO.AztecOO(A_epetra, x_epetra, b_epetra)
77
76
Solver.SetPrecOperator(Prec)
78
77
Solver.SetAztecOption(AztecOO.AZ_solver, AztecOO.AZ_cg)
79
78
Solver.SetAztecOption(AztecOO.AZ_output, 16)
80
79
Solver.Iterate(1550, 1e-5)
81
80
 
82
 
# Plot the solution 
 
81
# Plot the solution
83
82
plot(U)
84
83
interactive()
85
84