# Solving a system of linear equations with complex numbers yields false solution

I have equations for the analysis of an electric circuit. They contain at several places the term w*I (I is the imaginary unit). I get a solution for the currents I1, I2, I3, I4. Then I substitute the solution into the equations, but I see that two of the four equations are not fulfilled:

```
var('L1 L2 L3 L4 C1 C2 C3 C4 M23 I1 I2 I3 I4 U w p RL')
lsg = solve([I1*(w*I*L1 + 1/(w*I*C1)) + I2*(-1/(w*I*C1) ) == U, \
I2*(w*I*L2 + 1/(w*I*C1) + 1/(w*I*C2)) + I3*(w*I*M23 - 1/(w*I*C2)) + I1*(- 1/(w*I*C1)) == 0, \
I3*(1/(w*I*C3) + w*I*L3 + 1/(w*I*C2)) + I2*(w*I*M23 - 1/(w*I*C2)) - I4/(w*I*C3) == 0, \
I4/(w*I*C3) + I4*RL - I3/(w*I*C3) == 0], \
[I1, I2, I3, I4])
param = [w==1, U==1, C1==1, C2==1, C3==1, C4==1, L1==1, L2==1, L3==1, L4==1, RL==1, M23==0]
I1 = I1.subs(lsg).subs(param)
I2 = I2.subs(lsg).subs(param)
I3 = I3.subs(lsg).subs(param)
I4 = I4.subs(lsg).subs(param)
eqn = [I1*(w*I*L1 + 1/(w*I*C1)) + I2*(-1/(w*I*C1) ) == U, \
I2*(w*I*L2 + 1/(w*I*C1) + 1/(w*I*C2)) + I3*(w*I*M23 - 1/(w*I*C2)) + I1*(- 1/(w*I*C1)) == 0, \
I3*(1/(w*I*C3) + w*I*L3 + 1/(w*I*C2)) + I2*(w*I*M23 - 1/(w*I*C2)) - I4/(w*I*C3) == 0, \
I4/(w*I*C3) + I4*RL - I3/(w*I*C3) == 0]
print("I1=", I1)
print("I2=", I2)
print("I3=", I3)
print("I4=", I4)
[eq.subs(param) for eq in eqn]
```

Output (unexpected):

```
I1= 1
I2= -I
I3= -2
I4= I - 1
[1 == 1, (-I - 1) == 0, I == 0, 0 == 0]
```

I have copied and pasted the equations from within the solve() command into the eqn = ... statement. So they are guaranteed to be equal.

Because I have set M23=0 in the parameters, I could remove the terms with M23 in the equations. Then I get the correct solution. I don't understand why not when M23 is present.

Another way to get the correct solution is replacing w*I with p in the equations and setting p=I in the parameters:

```
var('L1 L2 L3 L4 C1 C2 C3 C4 M23 I1 I2 I3 I4 U w p RL')
lsg = solve([I1*(p*L1 + 1/(p*C1)) + I2*(-1/(p*C1) ) == U, \
I2*(p*L2 + 1/(p*C1) + 1/(p*C2)) + I3*(p*M23 - 1/(p*C2)) + I1*(- 1/(p*C1)) == 0, \
I3*(1/(p*C3) + p*L3 + 1/(p*C2)) + I2*(p*M23 - 1/(p*C2)) - I4/(p*C3) == 0, \
I4/(p*C3) + I4*RL - I3/(p*C3) == 0], \
[I1, I2, I3, I4])
param = [p==I, U==1, C1==1, C2==1, C3==1, C4==1, L1==1, L2==1, L3==1, L4==1, RL==1, M23==0]
I1 = I1.subs(lsg).subs(param)
I2 = I2.subs(lsg).subs(param)
I3 = I3.subs(lsg).subs(param)
I4 = I4.subs(lsg).subs(param)
eqn = [I1*(p*L1 + 1/(p*C1)) + I2*(-1/(p*C1) ) == U, \
I2*(p*L2 + 1/(p*C1) + 1/(p*C2)) + I3*(p*M23 - 1/(p*C2)) + I1*(- 1/(p*C1)) == 0, \
I3*(1/(p*C3) + p*L3 + 1/(p*C2)) + I2*(p*M23 - 1/(p*C2)) - I4/(p*C3) == 0, \
I4/(p*C3) + I4*RL - I3/(p*C3) == 0]
print("I1=", I1)
print("I2=", I2)
print("I3=", I3)
print("I4=", I4)
[eq.subs(param) for eq in eqn]
```

Output as expected:

```
I1= 1
I2= -I
I3= -I - 1
I4= -1
[1 == 1, 0 == 0, 0 == 0, 0 == 0]
```

Does anyone know why it yields false results in the first case (taking w*I instead of p) ?

If you deal with linear equations, it's worth to employ linear algebra machinery (matrices/vectors/etc) rather than symbolic

`solve`

function. It'll give you more control over what's going on behind the scene.Maybe you are right. But solving a system of linear equations with constants (real or complex) and variables is a standard task. And since sagemath allows me to do it also without matrices as I did, I would expect to get correct results. Maybe there is some documentation about sagemath that I haven't read and that could explain the behaviour (and in fact I have read only parts so far).