Solve system of equations

[Sancho_Pansa]

I've tried with fsolve.
Here is my code:

import math import numpy as np from scipy.optimize import fsolve def f(LCR): #variables L = LCR[0] Cs = LCR[1] Rs = LCR[2] # global constants Pi = math.pi Rin = 10.0 Cl = 3.0e-12 # measurement-related constants f0 = 130e6 f1 = 140e6 f2 = 150e6 omg = [None]*3 out = [None]*3 omg[0] = 2*Pi*f0 omg[1] = 2*Pi*f1 omg[2] = 2*Pi*f2 out[0] = 6.658482971121064 out[1] = 25.89226207742419 out[2] = 7.013950018010094 # A and B arrays Ar = [None]*3 Aj = [None]*3 Br = [None]*3 Bj = [None]*3 equ = [None]*3 # equation building for i in range(3): Ar[i] = 1 - omg[i]**2*Cl*L Aj[i] = omg[i]*Cs*Rs Br[i] = 1 - omg[i]**2*(L*(Cs + Cl) + Cl*Cs*Rin*Rs) Bj[i] = omg[i]*(Cs*Rs + Cl*Rin + Cl*Rs - omg[i]**2*Cl*Cs*Rin*L) equ[i] = out[i] - math.sqrt(Ar[i]*Ar[i] + Aj[i]*Aj[i])/math.sqrt(Br[i]*Br[i] + Bj[i]*Bj[i]) return equ LCR0 = np.array([0.3e-6, 1e-12, 1]) LCR = fsolve(f, LCR0) LCR

Here is output:

Output:
Warning (from warnings module): File "/home/pavel47/.local/lib/python3.8/site-packages/scipy/optimize/minpack.py", line 175 warnings.warn(msg, RuntimeWarning) RuntimeWarning: The iteration is not making good progress, as measured by the improvement from the last five Jacobian evaluations. >>> LCR array([ 3.10756682e-07, 1.02921301e-12, -4.22520679e+00]) >>> 

While the first two values are correct, the third is completely wrong.
Any suggestions ?
Thanks.

[Sancho_Pansa]

Here is updated version of code:

import math import numpy as np from scipy.optimize import fsolve def f(LCR, *args): omg = np.array([args[0], args[1], args[2]]) out = np.array([args[3], args[4], args[5]]) # variables L = LCR[0] Cs = LCR[1] Rs = LCR[2] # global constants Rin = 10.0 Cl = 3.0e-12 # declaration A, B and equ arrays Ar, Aj, Br, Bj, equ = [np.zeros(3) for _ in range(5)] # equation building for i in range(3): Ar[i] = 1 - omg[i]**2*Cl*L Aj[i] = omg[i]*Cs*Rs Br[i] = 1 - omg[i]**2*(L*(Cs + Cl) + Cl*Cs*Rin*Rs) Bj[i] = omg[i]*(Cs*Rs + Cl*Rin + Cl*Rs - omg[i]**2*Cl*Cs*Rin*L) #equ[i] = out[i] - math.sqrt(Ar[i]*Ar[i] + Aj[i]*Aj[i])/math.sqrt(Br[i]*Br[i] + Bj[i]*Bj[i]) equ[i] = out[i]*math.sqrt(Br[i]*Br[i] + Bj[i]*Bj[i]) - math.sqrt(Ar[i]*Ar[i] + Aj[i]*Aj[i]) return equ # Define frequencies in MHz f1 = 130 f2 = 140 f3 = 150 omega = np.array([f1, f2, f3])*1.0e6*2*math.pi # Define outputs: out1 = 6.658482971121064 out2 = 25.89226207742419 out3 = 7.013950018010094 out = np.array([out1, out2, out3]) # Define guesses for L [uH], Cs [pF], Rs [Ohm] L0 = 0.3e-6 Cs0 = 1e-12 Rs0 = 5 LCR0 = np.array([L0, Cs0, Rs0]) LCR = fsolve(f, LCR0, args=tuple(np.concatenate((omega, out))), maxfev=10000) print(LCR)

With initial conditions closer to the solution, the result is much worse.

[Gribouillis]

Unfortunately, I don't have much time to help you now but I think you reduced too fast the complex equations to real numbers. Go back to the complex equations and try to interprete them in geometrical terms. It seems to me that with some work you could express Cs in terms of L and Rs, then design an optimisation form of the problem. When I have the time to do so, I may answer in more details but I'm afraid it won't be before a week because I'm going to be on a journey for a few days.

[Sancho_Pansa]

It seems to me that with some work you could express Cs in terms of L and Rs

Frankly speaking I didn't understand your idea. Well ... after hard work I could do that. And after that ... what?
I think a solution can come from some parameters passing in fsolve.

[Askic]

Sancho_Pansa,

can you please elaborate more about the problem itself. You said you have 3 equations with 3 unknowns. The unknown variables that need to be found are: Cs, L, Rs
So the first equation is the equation for Br, second is for Bj and the third equation is for TF.
When you say that frequency omega is different, you mean in the first equation is omega1, in second is omega2 and in the third one is omega3 or maybe something else?

[Sancho_Pansa]

Sancho_Pansa,

can you please elaborate more about the problem itself. You said you have 3 equations with 3 unknowns. The unknown variables that need to be found are: Cs, L, Rs
So the first equation is the equation for Br, second is for Bj and the third equation is for TF.
When you say that frequency omega is different, you mean in the first equation is omega1, in second is omega2 and in the third one is omega3 or maybe something else?

Well, at the beginning I badly formulated the problem.
Here is correct setup:

• there are 3 variables: Cs, L, Rs
  
• there are 3 equations: How these 3 equations differ - at each frequency there is a specific output
  

[Askic]

I see,
this looks like as if the constrained optimization would probably do better (since L,C and R cannot be negative)

import math from scipy.optimize import least_squares def f(LCR): # variables L = LCR[0] Cs = LCR[1] Rs = LCR[2] # global constants Pi = math.pi Rin = 10.0 Cl = 3.0e-12 # measurement-related constants f0 = 130e6 f1 = 140e6 f2 = 150e6 omg = [0, 0, 0] out = [0, 0, 0] omg[0] = 2*Pi*f0 omg[1] = 2*Pi*f1 omg[2] = 2*Pi*f2 out[0] = 6.658482971121064 out[1] = 25.89226207742419 out[2] = 7.013950018010094 # A and B arrays Ar = [0, 0, 0] Aj = [0, 0, 0] Br = [0, 0, 0] Bj = [0, 0, 0] equ = [0, 0, 0] # equation building for i in range(3): Ar[i] = 1 - omg[i]**2*Cl*L Aj[i] = omg[i]*Cs*Rs Br[i] = 1 - omg[i]**2*(L*(Cs + Cl) + Cl*Cs*Rin*Rs) Bj[i] = omg[i]*(Cs*Rs + Cl*Rin + Cl*Rs - omg[i]**2*Cl*Cs*Rin*L) equ[i] = out[i] - math.sqrt(Ar[i]*Ar[i] + Aj[i] * Aj[i])/math.sqrt(Br[i]*Br[i] + Bj[i]*Bj[i]) return equ LCR0 = [0.3e-6, 1e-12, 1] res = least_squares(f, (LCR0), bounds=((0, 0, 0), (1.0e-6, 1.0e-10, 2))) print(res)

Let me know if this is closer to the expected solution.

[Sancho_Pansa]

Thanks.
1st value (i.e. L) is Ok.
2nd value (i.e. C) is lower than expected of 30%
The last one (i.e. R) is 40% of expected.
Here is code, that generates expected values:

import math Pi = math.pi for f_raw in (130, 140, 150): Cs_raw = 1.0 Cl_raw = 3.0 L_raw = 0.32 Cs = Cs_raw*1e-12 Cl = Cl_raw*1e-12 L = L_raw*1e-6 Rs = 5.0 Rin = 10.0 f = f_raw*1e6 omega = 2*Pi*f Ar = 1 - omega**2*Cs*Cl Aj = omega*Cs*Rs A = math.sqrt(Ar**2 + Aj**2) Br = 1 - omega**2*Cs*L - omega**2*Cl*L - omega**2*Cl*Cs*Rin*Rs Bj = omega*Cs*Rs + omega*Cl*Rin - omega**3*Cl*Cs*Rin*L + omega*Cl*Rs B = math.sqrt(Br**2 + Bj**2) print(f_raw, A/B)

[Askic]

Hello Sancho,

please double check your equations. I'm looking your last code. Is it Ar = 1 - omega**2*Cs*Cl or
Ar = 1 - omega**2*Cl*L???

[Sancho_Pansa]

Hello Sancho,

please double check your equations. I'm looking your last code. Is it Ar = 1 - omega**2*Cs*Cl or
Ar = 1 - omega**2*Cl*L???

Hello Askic.
Here is the solution, which works fine.
Modifications that was made:

• I replaced the capacitance Cl by the resistor Rl which drastically simplifies the transfer function, which in turn simplifies the work of the fsolve function.
  
• The equations are now generated by transformations performed on the circuit, created by the Circuit function, which minimizes typos.
  

from lcapy import Circuit, j, omega import numpy as np from math import pi from scipy.optimize import fsolve def f(LCR, *args): omg = np.array([args[0], args[1], args[2]]) out = np.array([args[3], args[4], args[5]]) cct = Circuit(''' ... Rin 1 2 ... Cs 2 4 ... La 2 3 ... Rs 3 4 ... Rl 4 0''') cct1 = cct.subs({'Rin': 50, 'Cs':LCR[1], 'La': LCR[0], 'Rs': LCR[2], 'Rl': 100}) H = cct1.transfer(1,0,4,0) A = H(j*omega).simplify() Am = A.magnitude equ = np.zeros(3) # equation building for i in range(3): equ[i] = out[i] - Am.evaluate(omg[i]) return equ # Define frequencies in MHz f1 = 130 f2 = 140 f3 = 150 omg = np.array([f1, f2, f3])*1.0e6*2*pi # Define outputs: out = [0.27175539, 0.24606606, 0.22193977] # Define guesses for L [uH], Cs [pF], Rs [Ohm] L0 = 0.3e-6 Cs0 = 1e-12 Rs0 = 1 LCR0 = np.array([L0, Cs0, Rs0]) LCR = fsolve(f, LCR0, args=tuple(np.concatenate((omg, out))), maxfev=10000) print(LCR)

Here is output:

Output:
[3.20000543e-07 9.99996130e-13 4.99951198e+00]

Here are real values:
0.32e-6, 1.0e-12, 5.0

Sincerely,

Sancho.