symbolic.py

   1 #!/usr/bin/python
   2
   3 from sympy import *
   4
   5 import sys
   6
   7 import sys, codecs
   8 if sys.stdout.encoding is None:
   9   sys.stdout = codecs.open("/dev/stdout", "w", 'utf-8')
  10
  11 init_printing(use_unicode=True)
  12
  13 r, theta, s, la, mu = symbols('r theta s lambda mu')
  14
  15 # start      original formulation
  16 # rightvars  replaces
  17
  18 p_start = Matrix([
  19   r * (1 - cos(theta)),
  20   r * sin(theta),
  21   mu * s,
  22 ])
  23
  24 p_rightvars = p_start.subs( theta, s/r ).subs( r, 1/la )
  25
  26 def dbg(*args):
  27   for vn in args:
  28     print('\n    ' + vn + '\n')
  29     pprint(eval(vn))
  30
  31 dbg('p_rightvars')
  32
  33 p_dirn_rightvars = diff(p_rightvars, s)
  34
  35 dbg('p_dirn_rightvars')
  36
  37 zeta = Wild('zeta')
  38
  39 p_nosing = (p_rightvars
  40             .replace( 1-cos(zeta)  ,   2*sin(zeta/2)**2          )
  41             .replace( sin(zeta)**2 ,   zeta*sinc(zeta)*sin(zeta) )
  42             )
  43 p_nosing[1] = (p_nosing[1]
  44             .replace( sin(zeta) , zeta * sinc(zeta)                )
  45                )
  46
  47 dbg('p_nosing')
  48
  49 t = symbols('t')
  50
  51 q_owncoords = p_nosing.replace(s,t).replace(la,-la)
  52 q_dirn_owncoords = p_dirn_rightvars.replace(s,t).replace(la,-la)
  53
  54 dbg('q_owncoords','q_dirn_owncoords')
  55 dbg('q_owncoords.replace(t,0)','q_dirn_owncoords.replace(t,0)')
  56
  57 p2q_translate = p_nosing
  58 #p2q_rotate_2d = Matrix([ p_dirn_rightvars[0:2],
  59
  60 #p2q_rotate = eye(3)
  61 #p2q_rotate[0:2, 0] = Matrix([ p_dirn_rightvars[1], -p_dirn_rightvars[0] ])
  62 #p2q_rotate[0:2, 1] = p_dirn_rightvars[0:2]
  63
  64 p2q_rotate = Matrix([[  cos(theta), sin(theta), 0 ],
  65                      [ -sin(theta), cos(theta), 0 ],
  66                      [  0         , 0,          1 ]]).subs(theta,la*s)
  67 #p2q_rotate.add_col([0,0])
  68 #p2q_rotate.add_row([0,0,1])
  69
  70 dbg('p2q_rotate')
  71
  72 q_dirn_maincoords = p2q_rotate * q_dirn_owncoords;
  73 q_maincoords = p2q_rotate * q_owncoords + p2q_translate
  74
  75 dbg('diff(p_dirn_rightvars,s)')
  76 dbg('diff(q_dirn_maincoords,t)')
  77 dbg('diff(q_dirn_maincoords,t).replace(t,0)')
  78
  79 assert(Eq(p2q_rotate * Matrix([0,1,mu]), p_dirn_rightvars))
  80
  81 #for v in 's','t','la','mu':
  82 #  dbg('diff(q_maincoords,%s)' % v)
  83
  84 #print('\n eye3 subs etc.\n')
  85 #dbg('''Eq(eye(3) * Matrix([1,0,mu]),
  86 #     p_dirn_rightvars .cross(Matrix([0,0,1]) .subs(s,0)))''')
  87
  88 #dbg('''Eq(p2q_rotate * Matrix([1,0,mu]),
  89 #          p_dirn_rightvars .cross(Matrix([0,0,1])))''')
  90
  91 #eq = Eq(qmat * q_dirn_owncoords_0, p_dirn_rightvars)
  92 #print
  93 #pprint(eq)
  94 #solve(eq, Q)
  95
  96 dbg('q_maincoords.replace(t,0)','q_dirn_maincoords.replace(t,0)')
  97
  98 dbg('q_maincoords','q_dirn_maincoords')