~njansson/dolfin/hpc

« back to all changes in this revision

Viewing changes to sandbox/la/test/demo.py

  • Committer: Garth N. Wells
  • Date: 2008-03-29 09:34:25 UTC
  • Revision ID: gnw20@cam.ac.uk-20080329093425-3ea2vhjvccq56zvi
Add some basic build & install instructions.

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
 
2
 
from dolfin import *
3
 
from sys import path 
4
 
path.append("../poisson")
5
 
 
6
 
 
7
 
from BlockLinearAlgebra import *
8
 
 
9
 
# Create mesh and finite element
10
 
N = 4 
11
 
mesh = UnitSquare(N,N)
12
 
mesh.disp()
13
 
CG = FiniteElement("Lagrange", "triangle", 1)
14
 
DG = FiniteElement("DG", "triangle", 0)
15
 
 
16
 
# Source term
17
 
class Source1(Function):
18
 
    def __init__(self, element, mesh):
19
 
        Function.__init__(self, element, mesh)
20
 
    def eval(self, values, x):
21
 
        values[0] = exp(x[0]*x[1]) 
22
 
 
23
 
# Source term
24
 
class Source2(Function):
25
 
    def __init__(self, element, mesh):
26
 
        Function.__init__(self, element, mesh)
27
 
    def eval(self, values, x):
28
 
        values[0] = cos(x[0]*x[1]) 
29
 
 
30
 
 
31
 
v = TestFunction(CG)
32
 
u = TrialFunction(DG)
33
 
w = TestFunction(DG)
34
 
f1 = Source1(CG, mesh)
35
 
f2 = Source2(DG, mesh)
36
 
f3 = Source2(CG, mesh)
37
 
 
38
 
a = u*v*dx + dot(grad(v), grad(u))*dx 
39
 
L1 = v*f1*dx 
40
 
L2 = w*f2*dx 
41
 
L3 = v*f2*dx 
42
 
 
43
 
# Assemble matrices
44
 
A = assemble(a, mesh)
45
 
b1 = assemble(L1, mesh)
46
 
b2 = assemble(L2, mesh)
47
 
b3 = assemble(L3, mesh)
48
 
 
49
 
file = File("A.m"); file << A; 
50
 
file = File("b1.m"); file << b1;  
51
 
file = File("b2.m"); file << b2;  
52
 
file = File("b3.m"); file << b3;  
53
 
 
54
 
print ""
55
 
print "testing Matrix Vector Product"
56
 
print " testing Vector (transposed)", 
57
 
x = b2.copy()
58
 
x.zero()
59
 
A.mult(b1, x, True)
60
 
print x.inner(x)
61
 
 
62
 
x.zero()
63
 
print " testing BlockVector (transposed)", 
64
 
AA = BlockMatrix(1,1); AA[0,0] = A; 
65
 
bb1 = BlockVector(1) ; bb1[0] = b1;  
66
 
xx = BlockVector(1); xx[0] = x; 
67
 
 
68
 
AA.prod(bb1, xx, True);
69
 
print xx.inner(xx); 
70
 
 
71
 
 
72
 
print ""
73
 
print "testing Matrix Vector Product"
74
 
print " testing Matrix::mult ", 
75
 
x = b1.copy()
76
 
x.zero()
77
 
A.mult(b2, x, False)
78
 
print x.inner(x)
79
 
print " testing Matrix::__mul__ ", 
80
 
x = A*b2
81
 
print x.inner(x)
82
 
 
83
 
 
84
 
print " testing BlockVector", 
85
 
 
86
 
x.zero()
87
 
bb2 = BlockVector(1) ; bb2[0] = b2;  
88
 
xx = BlockVector(1); xx[0] = x; 
89
 
 
90
 
AA.prod(bb2, xx, False);
91
 
print xx.inner(xx); 
92
 
 
93
 
 
94
 
print ""
95
 
print "testing Vector Addition "
96
 
tmp = b1.copy()
97
 
print " testing Vector", 
98
 
tmp += b3 
99
 
print tmp.inner(tmp)
100
 
 
101
 
tmp = b1.copy()
102
 
bb1 = BlockVector(1) ; bb1[0] = tmp;  
103
 
bb3 = BlockVector(1) ; bb3[0] = b3;  
104
 
bb1 += bb3
105
 
print " testing BlockVector", bb1.inner(bb1)
106
 
 
107
 
print ""
108
 
print "testing Vector Subtraction "
109
 
tmp = b1.copy()
110
 
print " testing Vector", 
111
 
tmp -= b3 
112
 
print tmp.inner(tmp)
113
 
 
114
 
tmp = b1.copy()
115
 
bb1 = BlockVector(1) ; bb1[0] = tmp;  
116
 
bb3 = BlockVector(1) ; bb3[0] = b3;  
117
 
bb1 -= bb3
118
 
print " testing BlockVector", bb1.inner(bb1)
119
 
 
120
 
 
121
 
print ""
122
 
print "testing Vector Multiplication "
123
 
tmp = b1.copy()
124
 
print " testing Vector", 
125
 
tmp *=3.7 
126
 
print tmp.inner(tmp)
127
 
 
128
 
tmp = b1.copy()
129
 
bb1 = BlockVector(1) ; bb1[0] = tmp;  
130
 
bb3 = BlockVector(1) ; bb3[0] = b3;  
131
 
bb1 *= 3.7
132
 
print " testing BlockVector", bb1.inner(bb1)
133
 
 
134
 
 
135
 
 
136
 
 
137
 
 
138
 
 
139
 
 
140
 
 
141
 
 
142
 
 
143
 
 
144