~dinko-metalac/calculus-app2/trunk

« back to all changes in this revision

Viewing changes to lib/py/sympy/galgebra/tests/test_GA.py

  • Committer: dinko.metalac at gmail
  • Date: 2015-04-14 13:28:14 UTC
  • Revision ID: dinko.metalac@gmail.com-20150414132814-j25k3qd7sq3warup
new sympy

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
#!/usr/bin/python
2
 
#test_GA.py
3
 
 
4
 
"""
5
 
The reference D&L is "Geometric Algebra for Physicists" by Doran and Lasenby
6
 
"""
7
 
 
8
 
import sys
9
 
 
10
 
from sympy.external import import_module
11
 
numpy = import_module('numpy')
12
 
if not numpy:
13
 
    disabled = True
14
 
else:
15
 
    sys.path.append('../')
16
 
    from sympy.galgebra.GA import MV, ZERO, HALF, S
17
 
    import sympy
18
 
 
19
 
 
20
 
def F(x, n, nbar):
21
 
    """
22
 
    Conformal Mapping Function from 3D Euclidean space to 5D conformal space
23
 
    where the images of all maps are null vectors.
24
 
    """
25
 
    Fx = HALF*((x*x)*n + 2*x - nbar)
26
 
    #print 'F(x) =',Fx
27
 
    return(Fx)
28
 
 
29
 
 
30
 
def test_rmul():
31
 
    """
32
 
    Test for commutative scalar multiplication.  Leftover from when sympy and
33
 
    numpy were not working together and __mul__ and __rmul__ would not give the
34
 
    same answer.
35
 
    """
36
 
    x, y, z = MV.setup('x y z')
37
 
    a, b, c = sympy.symbols('a b c')
38
 
    assert 5*x == x*5
39
 
    assert HALF*x == x*HALF
40
 
    assert a*x == x*a
41
 
 
42
 
 
43
 
def test_contraction():
44
 
    """
45
 
    Test for inner product and left and right contraction
46
 
    """
47
 
 
48
 
    e_1, e_2, e_3 = MV.setup('e_1 e_2 e_3', '1 0 0, 0 1 0, 0 0 1', offset=1)
49
 
 
50
 
    assert ((e_1 ^ e_3) | e_1) == -e_3
51
 
    assert ((e_1 ^ e_3) > e_1) == -e_3
52
 
    assert (e_1 | (e_1 ^ e_3)) == e_3
53
 
    assert (e_1 < (e_1 ^ e_3)) == e_3
54
 
    assert ((e_1 ^ e_3) < e_1) == 0
55
 
    assert (e_1 > (e_1 ^ e_3)) == 0
56
 
 
57
 
 
58
 
def test_substitution():
59
 
 
60
 
    e_x, e_y, e_z = MV.setup('e_x e_y e_z', '1 0 0, 0 1 0, 0 0 1', offset=1)
61
 
    x, y, z = sympy.symbols('x y z')
62
 
 
63
 
    X = x*e_x + y*e_y + z*e_z
64
 
    Y = X.subs([(x, 2), (y, 3), (z, 4)])
65
 
    assert Y == 2*e_x + 3*e_y + 4*e_z
66
 
 
67
 
 
68
 
def test_vector_extraction():
69
 
    """
70
 
    Show that conformal bivector encodes two points. See D&L Section 10.4.1
71
 
    """
72
 
    metric = ' 0 -1 #,' + \
73
 
             '-1  0 #,' + \
74
 
             ' #  # #,'
75
 
 
76
 
    P1, P2, a = MV.setup('P1 P2 a', metric)
77
 
    """
78
 
    P1 and P2 are null vectors and hence encode points in conformal space.
79
 
    Show that P1 and P2 can be extracted from the bivector B = P1^P2. a is a
80
 
    third vector in the conformal space with a.B not 0.
81
 
    """
82
 
    B = P1 ^ P2
83
 
    Bsq = B*B
84
 
    ap = a - (a ^ B)*B
85
 
    Ap = ap + ap*B
86
 
    Am = ap - ap*B
87
 
    P1dota = sympy.Symbol('(P1.a)')
88
 
    P2dota = sympy.Symbol('(P2.a)')
89
 
    Ap_test = (-2*P2dota)*P1
90
 
    Am_test = (-2*P1dota)*P2
91
 
    Ap.compact()
92
 
    Am.compact()
93
 
    Ap_test.compact()
94
 
    Am_test.compact()
95
 
    assert Ap == Ap_test
96
 
    assert Am == Am_test
97
 
    Ap2 = Ap*Ap
98
 
    Am2 = Am*Am
99
 
    Ap2.compact()
100
 
    Am2.compact()
101
 
    assert Ap2 == ZERO
102
 
    assert Am2 == ZERO
103
 
 
104
 
def test_metrics():
105
 
    """
106
 
    Test specific metrics (diagpq, arbitrary_metric, arbitrary_metric_conformal)
107
 
    """
108
 
    from sympy.galgebra.GA import diagpq, arbitrary_metric, arbitrary_metric_conformal
109
 
    metric = diagpq(3)
110
 
    p1, p2, p3 = MV.setup('p1 p2 p3', metric, debug=0)
111
 
    MV.set_str_format(1)
112
 
    x1, y1, z1 = sympy.symbols('x1 y1 z1')
113
 
    x2, y2, z2 = sympy.symbols('x2 y2 z2')
114
 
    v1 = x1*p1 + y1*p2 + z1*p3
115
 
    v2 = x2*p1 + y2*p2 + z2*p3
116
 
    prod1 = v1*v2
117
 
    prod2 = (v1|v2) + (v1^v2)
118
 
    diff = prod1 - prod2
119
 
    diff.compact()
120
 
    assert diff == ZERO
121
 
    metric = arbitrary_metric(3)
122
 
    p1, p2, p3 = MV.setup('p1 p2 p3', metric, debug=0)
123
 
    v1 = x1*p1 + y1*p2 + z1*p3
124
 
    v2 = x2*p1 + y2*p2 + z2*p3
125
 
    prod1 = v1*v2
126
 
    prod2 = (v1|v2) + (v1^v2)
127
 
    diff = prod1 - prod2
128
 
    diff.compact()
129
 
    assert diff == ZERO
130
 
    metric = arbitrary_metric_conformal(3)
131
 
    p1, p2, p3 = MV.setup('p1 p2 p3', metric, debug=0)
132
 
    v1 = x1*p1 + y1*p2 + z1*p3
133
 
    v2 = x2*p1 + y2*p2 + z2*p3
134
 
    prod1 = v1*v2
135
 
    prod2 = (v1|v2) + (v1^v2)
136
 
    diff = prod1 - prod2
137
 
    diff.compact()
138
 
    assert diff == ZERO
139
 
 
140
 
def test_geometry():
141
 
    """
142
 
    Test conformal geometric description of circles, lines, spheres, and planes.
143
 
    """
144
 
    metric = '1 0 0 0 0,' + \
145
 
             '0 1 0 0 0,' + \
146
 
             '0 0 1 0 0,' + \
147
 
             '0 0 0 0 2,' + \
148
 
             '0 0 0 2 0'
149
 
 
150
 
    e0, e1, e2, n, nbar = MV.setup('e0 e1 e2 n nbar', metric, debug=0)
151
 
    e = n + nbar
152
 
    #conformal representation of points
153
 
 
154
 
    A = F(e0, n, nbar)    # point a = (1,0,0)  A = F(a)
155
 
    B = F(e1, n, nbar)    # point b = (0,1,0)  B = F(b)
156
 
    C = F(-1*e0, n, nbar)  # point c = (-1,0,0) C = F(c)
157
 
    D = F(e2, n, nbar)    # point d = (0,0,1)  D = F(d)
158
 
    x0, x1, x2 = sympy.symbols('x0 x1 x2')
159
 
    X = F(MV([x0, x1, x2], 'vector'), n, nbar)
160
 
 
161
 
    Circle = A ^ B ^ C ^ X
162
 
    Line = A ^ B ^ n ^ X
163
 
    Sphere = A ^ B ^ C ^ D ^ X
164
 
    Plane = A ^ B ^ n ^ D ^ X
165
 
 
166
 
    #Circle through a, b, and c
167
 
    Circle_test = -x2*(e0 ^ e1 ^ e2 ^ n) + x2*(
168
 
        e0 ^ e1 ^ e2 ^ nbar) + HALF*(-1 + x0**2 + x1**2 + x2**2)*(e0 ^ e1 ^ n ^ nbar)
169
 
    diff = Circle - Circle_test
170
 
    diff.compact()
171
 
    assert diff == ZERO
172
 
 
173
 
    #Line through a and b
174
 
    Line_test = -x2*(e0 ^ e1 ^ e2 ^ n) + \
175
 
        HALF*(-1 + x0 + x1)*(e0 ^ e1 ^ n ^ nbar) + \
176
 
        (HALF*x2)*(e0 ^ e2 ^ n ^ nbar) + \
177
 
        (-HALF*x2)*(e1 ^ e2 ^ n ^ nbar)
178
 
    diff = Line - Line_test
179
 
    diff.compact()
180
 
    assert diff == ZERO
181
 
 
182
 
    #Sphere through a, b, c, and d
183
 
    Sphere_test = HALF*(1 - x0**2 - x1**2 - x2**2)*(e0 ^ e1 ^ e2 ^ n ^ nbar)
184
 
    diff = Sphere - Sphere_test
185
 
    diff.compact()
186
 
    assert diff == ZERO
187
 
 
188
 
    #Plane through a, b, and d
189
 
    Plane_test = HALF*(1 - x0 - x1 - x2)*(e0 ^ e1 ^ e2 ^ n ^ nbar)
190
 
    diff = Plane - Plane_test
191
 
    diff.compact()
192
 
    assert diff == ZERO
193
 
 
194
 
 
195
 
def test_extract_plane_and_line():
196
 
    """
197
 
    Show that conformal trivector encodes planes and lines. See D&L section
198
 
    10.4.2
199
 
    """
200
 
    metric = '# # # 0 0,' + \
201
 
             '# # # 0 0,' + \
202
 
             '# # # 0 0,' + \
203
 
             '0 0 0 0 2,' + \
204
 
             '0 0 0 2 0'
205
 
 
206
 
    p1, p2, p3, n, nbar = MV.setup('p1 p2 p3 n nbar', metric, debug=0)
207
 
    MV.set_str_format(1)
208
 
 
209
 
    P1 = F(p1, n, nbar)
210
 
    P2 = F(p2, n, nbar)
211
 
    P3 = F(p3, n, nbar)
212
 
 
213
 
    #Line through p1 and p2
214
 
    L = P1 ^ P2 ^ n
215
 
    delta = (L | n) | nbar
216
 
    delta_test = 2*p1 - 2*p2
217
 
    diff = delta - delta_test
218
 
    diff.compact()
219
 
    assert diff == ZERO
220
 
 
221
 
    #Plane through p1, p2, and p3
222
 
    C = P1 ^ P2 ^ P3
223
 
    delta = ((C ^ n) | n) | nbar
224
 
    delta_test = 2*(p1 ^ p2) - 2*(p1 ^ p3) + 2*(p2 ^ p3)
225
 
    diff = delta - delta_test
226
 
    diff.compact()
227
 
    assert diff == ZERO
228
 
 
229
 
 
230
 
def test_reciprocal_frame():
231
 
    """
232
 
    Test of formula for general reciprocal frame of three vectors.
233
 
    Let three independent vectors be e1, e2, and e3. The reciprocal
234
 
    vectors E1, E2, and E3 obey the relations:
235
 
 
236
 
    e_i.E_j = delta_ij*(e1^e2^e3)**2
237
 
    """
238
 
    metric = '1 # #,' + \
239
 
             '# 1 #,' + \
240
 
             '# # 1,'
241
 
 
242
 
    e1, e2, e3 = MV.setup('e1 e2 e3', metric)
243
 
    E = e1 ^ e2 ^ e3
244
 
    Esq = (E*E)()
245
 
    Esq_inv = 1/Esq
246
 
    E1 = (e2 ^ e3)*E
247
 
    E2 = (-1)*(e1 ^ e3)*E
248
 
    E3 = (e1 ^ e2)*E
249
 
    w = (E1 | e2)
250
 
    w.collect(MV.g)
251
 
    w = w().expand()
252
 
    w = (E1 | e3)
253
 
    w.collect(MV.g)
254
 
    w = w().expand()
255
 
    assert w == 0
256
 
    w = (E2 | e1)
257
 
    w.collect(MV.g)
258
 
    w = w().expand()
259
 
    assert w == 0
260
 
    w = (E2 | e3)
261
 
    w.collect(MV.g)
262
 
    w = w().expand()
263
 
    assert w == 0
264
 
    w = (E3 | e1)
265
 
    w.collect(MV.g)
266
 
    w = w().expand()
267
 
    assert w == 0
268
 
    w = (E3 | e2)
269
 
    w.collect(MV.g)
270
 
    w = w().expand()
271
 
    assert w == 0
272
 
    w = (E1 | e1)
273
 
    w = w().expand()
274
 
    Esq = Esq.expand()
275
 
    assert w/Esq == 1
276
 
    w = (E2 | e2)
277
 
    w = w().expand()
278
 
    assert w/Esq == 1
279
 
    w = (E3 | e3)
280
 
    w = w().expand()
281
 
    assert w/Esq == 1
282
 
 
283
 
 
284
 
def test_derivative():
285
 
    coords = x, y, z = sympy.symbols('x y z')
286
 
    e_x, e_y, e_z = MV.setup('e', '1 0 0, 0 1 0, 0 0 1', coords=coords)
287
 
    X = x*e_x + y*e_y + z*e_z
288
 
    a = MV('a', 'vector')
289
 
 
290
 
    assert ((X | a).grad()) == a
291
 
    assert ((X*X).grad()) == 2*X
292
 
    assert (X*X*X).grad() == 5*X*X
293
 
    assert X.grad_int() == 3
294
 
 
295
 
 
296
 
def test_str():
297
 
    e_1, e_2, e_3 = MV.setup('e_1 e_2 e_3', '1 0 0, 0 1 0, 0 0 1')
298
 
 
299
 
    X = MV('x')
300
 
    assert str(X) == 'x+x__0*e_1+x__1*e_2+x__2*e_3+x__01*e_1e_2+x__02*e_1e_3+x__12*e_2e_3+x__012*e_1e_2e_3'
301
 
    Y = MV('y', 'spinor')
302
 
    assert str(Y) == 'y+y__01*e_1e_2+y__02*e_1e_3+y__12*e_2e_3'
303
 
    Z = X + Y
304
 
    assert str(Z) == 'x+y+x__0*e_1+x__1*e_2+x__2*e_3+(x__01+y__01)*e_1e_2+(x__02+y__02)*e_1e_3+(x__12+y__12)*e_2e_3+x__012*e_1e_2e_3'
305
 
    assert str(e_1 | e_1) == '1'
306
 
 
307
 
 
308
 
def test_metric():
309
 
    MV.setup('e_1 e_2 e_3', '[1,1,1]')
310
 
    assert str(MV.metric) == '[[1 0 0]\n [0 1 0]\n [0 0 1]]'
311
 
 
312
 
 
313
 
def test_constructor():
314
 
    """
315
 
    Test various multivector constructors
316
 
    """
317
 
    e_1, e_2, e_3 = MV.setup('e_1 e_2 e_3', '[1,1,1]')
318
 
    x = sympy.symbols('x')
319
 
    assert str(S(1)) == '1'
320
 
    assert str(S(x)) == 'x'
321
 
    assert str(MV('a', 'scalar')) == 'a'
322
 
    assert str(MV('a', 'vector')) == 'a__0*e_1+a__1*e_2+a__2*e_3'
323
 
    assert str(MV('a', 'pseudo')) == 'a*e_1e_2e_3'
324
 
    assert str(MV('a', 'spinor')) == 'a+a__01*e_1e_2+a__02*e_1e_3+a__12*e_2e_3'
325
 
    assert str(MV('a')) == 'a+a__0*e_1+a__1*e_2+a__2*e_3+a__01*e_1e_2+a__02*e_1e_3+a__12*e_2e_3+a__012*e_1e_2e_3'
326
 
    assert str(
327
 
        MV([2, 'a'], 'grade')) == 'a__01*e_1e_2+a__02*e_1e_3+a__12*e_2e_3'
328
 
    assert str(MV('a', 'grade2')) == 'a__01*e_1e_2+a__02*e_1e_3+a__12*e_2e_3'
329
 
 
330
 
 
331
 
def test__print_Mul_Add():
332
 
    from sympy.galgebra.latex_ex import LatexPrinter
333
 
    from sympy import symbols
334
 
    n, m = symbols('n,m', negative=True)
335
 
    z = symbols('z')
336
 
    l = LatexPrinter()
337
 
    assert l._print_Mul(n*m) == 'm n'
338
 
    assert l._print_Mul(-2*m) == '- 2 m'
339
 
    assert l._print_Mul(2*m) == '2 m'
340
 
    assert l._print_Add(-5 + 4*z) == '-5 + 4 z'
341
 
    assert l._print_Add(-5 - 4*z) == '-5 - 4 z'
342
 
    assert l._print_Add(n - 2) == '-2 + n'