I have a system of longish equations and when calling .solve() on it, I get
this traceback:
Traceback (most recent call last):
File "C:/Users/Andreas
Schuldei/PycharmProjects/lissajous-achse/hgü-kabel-detektion-symbolic.py",
line 109, in <module>
result = solve(equations, (CM0, theta_i, theta_j, theta_k, i, a,
B_earth))
File "C:\Users\Andreas
Schuldei\PycharmProjects\lissajous-achse\venv\lib\site-packages\sympy\solvers\solvers.py",
line 1096, in solve
solution = _solve_system(f, symbols, **flags)
File "C:\Users\Andreas
Schuldei\PycharmProjects\lissajous-achse\venv\lib\site-packages\sympy\solvers\solvers.py",
line 1730, in _solve_system
i, d = _invert(g, *symbols)
File "C:\Users\Andreas
Schuldei\PycharmProjects\lissajous-achse\venv\lib\site-packages\sympy\solvers\solvers.py",
line 3118, in _invert
rhs -= indep
File "C:\Users\Andreas
Schuldei\PycharmProjects\lissajous-achse\venv\lib\site-packages\sympy\core\numbers.py",
line 2194, in __sub__
return Rational.__sub__(self, other)
File "C:\Users\Andreas
Schuldei\PycharmProjects\lissajous-achse\venv\lib\site-packages\sympy\core\decorators.py",
line 89, in __sympifyit_wrapper
return func(a, b)
File "C:\Users\Andreas
Schuldei\PycharmProjects\lissajous-achse\venv\lib\site-packages\sympy\core\numbers.py",
line 1725, in __sub__
return Number.__sub__(self, other)
File "C:\Users\Andreas
Schuldei\PycharmProjects\lissajous-achse\venv\lib\site-packages\sympy\core\decorators.py",
line 89, in __sympifyit_wrapper
return func(a, b)
File "C:\Users\Andreas
Schuldei\PycharmProjects\lissajous-achse\venv\lib\site-packages\sympy\core\numbers.py",
line 733, in __sub__
return AtomicExpr.__sub__(self, other)
File "C:\Users\Andreas
Schuldei\PycharmProjects\lissajous-achse\venv\lib\site-packages\sympy\core\decorators.py",
line 251, in _func
return func(self, other)
File "C:\Users\Andreas
Schuldei\PycharmProjects\lissajous-achse\venv\lib\site-packages\sympy\core\decorators.py",
line 126, in binary_op_wrapper
return f(self)
File "C:\Users\Andreas
Schuldei\PycharmProjects\lissajous-achse\venv\lib\site-packages\sympy\core\decorators.py",
line 127, in binary_op_wrapper
return func(self, other)
File "C:\Users\Andreas
Schuldei\PycharmProjects\lissajous-achse\venv\lib\site-packages\sympy\vector\basisdependent.py",
line 352, in __rsub__
raise TypeError("Invalid argument types for subtraction")
TypeError: Invalid argument types for subtraction
Process finished with exit code 1
I stared at my code long and hard, and fixed all instances of wrong types
that i could think of. Also, in order to catch the error earlier, I tried
insterting .simplify() calls in places that might benefit from them - but
the process ran for hours, doing the simplyfy() without even reaching the
solve() call.
So I am asking for tricks to investigate the types and them fitting
together, earlier. You can find my code attached. (please critique it, I
want to learn!).
--
You received this message because you are subscribed to the Google Groups
"sympy" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to [email protected].
To view this discussion on the web visit
https://groups.google.com/d/msgid/sympy/68a215e3-db62-4c4a-b3e7-d028f6f54d7bn%40googlegroups.com.
import sympy
from sympy.vector import CoordSys3D, express
from sympy import symbols, Eq, solve
from sympy.vector import AxisOrienter
#import numpy as np
pi = sympy.pi
mu_0 = symbols("mu_0")
Sys_sensors = CoordSys3D('Sys_sensors')
O = Sys_sensors.origin
# Earths magnetic field along the sensors vector components
B_x, B_y, B_z = symbols("B_x, B_y, B_z")
B_earth = B_x * Sys_sensors.i + B_y * Sys_sensors.j + B_z * Sys_sensors.k
#sensor_coordinates = np.array([(.265, 0, .382), (0, .712, .764), (0, .712, 0),
(.752, .712, .382)])
equations = []
# Points of measurement 1-4 (Sensors) (all sensors are aligned and in the same
Frame of reference)
m01, m02, m03, m11, m12, m13, m21, m22, m23, m31, m32, m33 = symbols('m:4(1:4)')
M0 = O.locate_new('M0', m01 * Sys_sensors.i + m02 * Sys_sensors.j + m03 *
Sys_sensors.k)
M1 = O.locate_new('M1', m11 * Sys_sensors.i + m12 * Sys_sensors.j + m13 *
Sys_sensors.k)
M2 = O.locate_new('M2', m21 * Sys_sensors.i + m22 * Sys_sensors.j + m23 *
Sys_sensors.k)
M3 = O.locate_new('M3', m31 * Sys_sensors.i + m32 * Sys_sensors.j + m33 *
Sys_sensors.k)
# vector from the Origin to each sensor
OM0 = O.position_wrt(M0)
OM1 = O.position_wrt(M1)
OM2 = O.position_wrt(M2)
OM3 = O.position_wrt(M3)
# each of the above points, projected on to the cable
c01, c02, c03, c11, c12, c13, c21, c22, c23, c31, c32, c33 = symbols('c:4(1:4)')
C0 = O.locate_new('C0', c01 * Sys_sensors.i + c02 * Sys_sensors.j + c03 *
Sys_sensors.k)
C1 = O.locate_new('C1', c11 * Sys_sensors.i + c12 * Sys_sensors.j + c13 *
Sys_sensors.k)
C2 = O.locate_new('C2', c21 * Sys_sensors.i + c22 * Sys_sensors.j + c23 *
Sys_sensors.k)
C3 = O.locate_new('C3', c31 * Sys_sensors.i + c32 * Sys_sensors.j + c33 *
Sys_sensors.k)
# shortest vector from each Sensor to the cable
CM0 = C0.position_wrt(M0)
CM1 = C1.position_wrt(M1)
CM2 = C2.position_wrt(M2)
CM3 = C3.position_wrt(M3)
# vectors along the cable
C01 = C0.position_wrt(C1)
C02 = C0.position_wrt(C2)
C03 = C0.position_wrt(C3)
# Magnetic field vector measured in the sensors 0...3
B01, B02, B03, B11, B12, B13, B21, B22, B23, B31, B32, B33 = symbols('B:4(1:4)')
B0 = B01*Sys_sensors.i + B02*Sys_sensors.j + B03*Sys_sensors.k
B1 = B11*Sys_sensors.i + B12*Sys_sensors.j + B13*Sys_sensors.k
B2 = B21*Sys_sensors.i + B22*Sys_sensors.j + B23*Sys_sensors.k
B3 = B31*Sys_sensors.i + B32*Sys_sensors.j + B33*Sys_sensors.k
# this notation implies that the dot product is zero, the vectors stand
perpendicular on each other
equations.append(Eq(C01.dot(CM0), 0))
for (cable_vect, sensor_vect) in zip([CM1, CM2, CM3], [C01, C02, C03]):
equations.append(Eq(cable_vect.dot(sensor_vect), 0).simplify())
# angels that the cables's reference frame is tilted by against the measurement
coordinate system
theta_i, theta_j, theta_k = symbols("theta_i, theta_j, theta_k")
axis_orienter_i = AxisOrienter(theta_i, Sys_sensors.i)
axis_orienter_j = AxisOrienter(theta_j, Sys_sensors.j)
axis_orienter_k = AxisOrienter(theta_k, Sys_sensors.k)
# the cable is rotated by the angels theta_i, _j, _k and translated by the
vector CM0+OM0, from the origin
a, i = symbols("a, i") # a: diameter of the cable, as well as (fixed)
orientation of the dipol in the x-y plane,
# note: the rotation of the dipol field in relation to the sensor frame of
reference is encoded in the thetas above. The rotation around the cable's
z-axis is equivalent to the rotation of the a vector in the cable's xy-plane.
# i: current in the cable
# the a vector is in the x-y plane. for my purpose its ok if its fixed - say
along x-axis. a is in the cable's frame of reference, C
# a = C.x + 0*C.y + 0*C.z #.... i dont know how to do those vectors
# this is in the cable's frame of reference, C
# this function is utterly broken. regard as pseudo-code.
def B_el(r):
a1 = a
K = mu_0 * i / pi
r1 = r.coeff(C.i)
r2 = r.coeff(C.j)
r3 = r.coeff(C.k) # = 0, because
# - we chose our vectors that way
# - the cable's field is invariant in z-direction
# u = ((2 * K * r2 * (a1 * r1 + a2 * r2)) / (r1 * r1 + r2 * r2) ** 2 - (K *
a2) / (r1 * r1 + r2 * r2)) * C.i
# v = ((K * a1) / (r1 * r1 + r2 * r2) - (2 * K * r1 * (a1 * r1 + a2 * r2)) /
(r1 * r1 + r2 * r2) ** 2) * C.j
u = ((2 * K * r2 * (a1 * r1 )) / (r1 * r1 + r2 * r2) ** 2 ) * C.i
v = ((K * a1) / (r1 * r1 + r2 * r2) - (2 * K * r1 * (a1 * r1)) / (r1 * r1 +
r2 * r2) ** 2) * C.j
w = 0 * C.k
return u + v + w # =B
for (cable_vect, sensor_vect, magnetic_vect) in zip([CM0, CM1, CM2, CM3], [OM0,
OM1, OM2, OM3], [B0, B1, B2, B3]):
location_vect = cable_vect + sensor_vect
# the cable is rotated by the angels theta_i, _j, _k and translated by the
vector location_vect, from the origin
C = Sys_sensors.orient_new('C', [axis_orienter_i, axis_orienter_j,
axis_orienter_k], location=location_vect)
radius_vect = express(cable_vect, C)
# calculate magnetic field, transform back to S coordinate system and add
up all magnetic fields
B_cable_sys = B_el(radius_vect)
B_sensor_sys = express(B_cable_sys, Sys_sensors)
equations.append(Eq(B_sensor_sys + B_earth - magnetic_vect,
0*Sys_sensors.i+0*Sys_sensors.j+0*Sys_sensors.k))
print(equations)
result = solve(equations, (CM0, theta_i, theta_j, theta_k, i, a, B_earth))
print(result)
