tlinalg_dolfin.py - pism - [fork] customized build of PISM, the parallel ice sheet model (tillflux branch)
(HTM) git clone git://src.adamsgaard.dk/pism
(DIR) Log
(DIR) Files
(DIR) Refs
(DIR) LICENSE
---
tlinalg_dolfin.py (3468B)
---
1 ############################################################################
2 #
3 # This file is a part of siple.
4 #
5 # Copyright 2010, 2014 David Maxwell
6 #
7 # siple is free software: you can redistribute it and/or modify
8 # it under the terms of the GNU General Public License as published by
9 # the Free Software Foundation, either version 2 of the License, or
10 # (at your option) any later version.
11 #
12 ############################################################################
13
14 from .linalg_abstract import AbstractVector
15 import dolfin
16
17 class DolfinGenericVector(AbstractVector):
18 """Implements the siple.linalg.AbstractVector interface for a dolfin generic vector."""
19 def __init__(self,u):
20 if isinstance(u,dolfin.Function):
21 self._core = u.vector()
22 elif isinstance(u,dolfin.GenericVector):
23 self._core = u
24 else:
25 raise ValueError("An DolfinGenericVector only be constructed from a dolfin Function or a dolfin GenericVector: found a %s" % u)
26
27 def _set_from_abstract(self,rhs):
28 self._core[:] = rhs._core[:]
29
30 def _set_from_array(self,rhs):
31 self._core[:] = rhs[:]
32
33 def acc(self,rhs):
34 self._core += rhs._core
35
36 def scale(self,t):
37 self._core *= t
38
39 def axpy(self,t,v):
40 self._core.axpy(t,v.core())
41
42 def copy(self):
43 return DolfinGenericVector(self._core.copy())
44
45 def vector_like(self):
46 w = self._core.factory().create_vector();
47 w.resize(self._core.size())
48 return DolfinGenericVector(w)
49
50 def zero_like(self):
51 z = self.vector_like()
52 z._core.zero()
53 return z
54
55 def dim(self):
56 return self._core.size()
57
58 def core(self):
59 return self._core
60
61 def norm(self,name):
62 return self._core.norm(name)
63
64 class DolfinFunctionVector(AbstractVector):
65 """Implements the siple.linalg.AbstractVector interface for a dolfin Function."""
66
67 def __init__(self,u):
68 if isinstance(u,dolfin.Function):
69 self._core = u
70 elif isinstance(u,dolfin.FunctionSpace):
71 self._core = dolfin.Function(u)
72 else:
73 raise ValueError("An DolfinFunctionSpaceVector wraps a dolfin Function or FunctionSpace: found a %s" % u)
74
75 def set(self,rhs):
76 self._core.vector()[:] = rhs._core.vector()[:]
77
78 def acc(self,rhs):
79 self._core.vector()[:] += rhs._core.vector()[:]
80
81 def scale(self,t):
82 self._core.vector()[:] *= t
83
84 def axpy(self,t,v):
85 self._core.vector().axpy(t,v.core().vector())
86
87 def copy(self):
88 rv = dolfin.Function(self._core.function_space())
89 rv.vector()[:] = self._core.vector()[:]
90 return DolfinFunctionSpaceVector(rv)
91
92 def vector_like(self):
93 w = dolfin.Function(self._core.function_space())
94 return DolfinFunctionSpaceVector(w)
95
96 def zero_like(self):
97 z = self.vector_like()
98 z._core.vector().zero()
99 return z
100
101 def dim(self):
102 return self._core.vector().size()
103
104 def core(self):
105 return self._core
106
107 def norm(self,name):
108 return self._core.vector().norm(name)
109
110
111
112 def apply_form(F,x,y,scratch=None):
113 vx = x
114 if (isinstance(x,dolfin.Function)):
115 vx = x.vector();
116
117 vy = y
118 if (isinstance(y,dolfin.Function)):
119 vy = y.vector();
120
121 if scratch is None:
122 Fy = F*vy
123 else:
124 F.mult(vy,scratch)
125 Fy = scratch
126 return vx.inner(Fy)
127
128 def lhsVectorFor(A,out=None):
129 outv=None
130 if not out is None:
131 outv = out.core()
132 if (outv is None) or (outv.size() != A.size(0)):
133 outv = A.factory().create_vector()
134 outv.resize(A.size(0))
135 out = DolfinGenericVector(outv)
136 return out