0

I want to solve this non linear equation: f100 = omega_nf_eq

where: f100 : numerical costant, defined as a variable for now.

omega_nf_eq: equation.

Firstly I've tried to solve it sybolically and my code was:

import sympy as sym

K_u, K_m = sym.symbols('K_u, K_m', real = True)
J_p1, J_p2, J_g1, J_g2, J_r, J_u, J_m, J_p12, J_g12, J_gb, J_2, J_1, J_p = sym.symbols('J_p1, J_p2, J_g1, J_g2, J_r, J_u, J_m, J_p12, J_g12, J_gb, J_2, J_1, J_p', real = True)
tau_1, tau_2 = sym.symbols('tau_1, tau_2', real = True)
omega_nf, f100 = sym.symbols('omega_nf, f100', real = True)


omega_nf_eq = sym.Eq(omega_nf, sym.sqrt(2)*sym.sqrt(K_m/(J_g2*tau_2**2 + J_p1 + J_p2) + K_u/(J_g2*tau_2**2 + J_p1 + J_p2) + K_u/(tau_2**2*(J_g1 + J_u)) + K_m/J_m - sym.sqrt(J_m**2*K_m**2*tau_2**4*(J_g1 + J_u)**2 + 2*J_m**2*K_m*K_u*tau_2**4*(J_g1 + J_u)**2 - 2*J_m**2*K_m*K_u*tau_2**2*(J_g1 + J_u)*(J_g2*tau_2**2 + J_p1 + J_p2) + J_m**2*K_u**2*tau_2**4*(J_g1 + J_u)**2 + 2*J_m**2*K_u**2*tau_2**2*(J_g1 + J_u)*(J_g2*tau_2**2 + J_p1 + J_p2) + J_m**2*K_u**2*(J_g2*tau_2**2 + J_p1 + J_p2)**2 + 2*J_m*K_m**2*tau_2**4*(J_g1 + J_u)**2*(J_g2*tau_2**2 + J_p1 + J_p2) - 2*J_m*K_m*K_u*tau_2**4*(J_g1 + J_u)**2*(J_g2*tau_2**2 + J_p1 + J_p2) - 2*J_m*K_m*K_u*tau_2**2*(J_g1 + J_u)*(J_g2*tau_2**2 + J_p1 + J_p2)**2 + K_m**2*tau_2**4*(J_g1 + J_u)**2*(J_g2*tau_2**2 + J_p1 + J_p2)**2)/(J_m*tau_2**2*(J_g1 + J_u)*(J_g2*tau_2**2 + J_p1 + J_p2)))/2)



solution = sym.solve(f100 - omega_nf_eq.args[1], J_u, dict = True) 

But this gave me just this result: [ ].

I've also tried to substitute all variable value except for J_u, which is the one i want. So now the omega_nf equation is:

omega_nf_eq = sym.Eq(omega_nf, sym.sqrt(2)*sym.sqrt(76019006.3529542 - 84187769.0684942*sym.sqrt(0.813040126459949*J_u**2 - 4.69199504596906e-5*J_u + 1.03236146920168e-9)/J_u + 2704.98520837442/J_u)/2)

So to solve now i've tried:

solution = sym.solve( 942.5 - omega_nf_eq.args[1], J_u,, dict = True, force=True, manual=True, set=True)

It works now, but it requires a couple of minutes.

So I've tried to solve it numerically, to speed up the process, with sympy.nsolve(); this is the code:

omega_nf_eq = sym.Eq(omega_nf, sym.sqrt(2)*sym.sqrt(76019006.3529542 - 84187769.0684942*sym.sqrt(0.813040126459949*J_u**2 - 4.69199504596906e-5*J_u + 1.03236146920168e-9)/J_u + 2704.98520837442/J_u)/2)

eq_solution = sym.nsolve(942.5 - omega_nf_eq, J_u, 0.0071, verify=False)

But i do not obtain the right result, which is: J_u = 0.00717865789803973.

What I'm doing wrong? There's a smarter way to use sympy?

2 Answers2

0

There is not J_u in your first symbolic equation so that's why you got [] for the solution. When you attempted a numerical solution you used omega_nf_eq (which is an Equality); I think you meant 'nsolve(942.5 - omega_nf_eq.rhs, J_u, .0071)'. But even still, that won't find a solution for you since the equation, as written, is ill-behaved with J_u in the denominator. If you use sympy.solvers.solvers.unrad to give you to radical-free expression whose roots will contain, as a subset, those you are interested in you will find that you need only solve a quadratic in J_u...and that will be fast.

>>> unrad(942.5 - omega_nf_eq.rhs)
(1.0022170762796e+15*J_u**2 - 2936792314038.5*J_u + 2.04890966415405e-7, [])
>>> solve(_[0])
[6.97669240810738e-20, 0.00293029562511584]

I would recommend that you revist your first symbolic expression and unrad that -- or even just try solve that -- after identifying which variable corresponds to J_u.

smichr
  • 16,948
  • 2
  • 27
  • 34
  • Thanks a lot. I made a mistake, I corrected it, but now, even after half an hour, sympy has not found a solution. This is a problem. I have to use it also to solve other equations on this project and I can't get the solutions, damn it. – cicciobombacannoniere Jun 06 '19 at 21:29
  • Please update the expression in your opening paragraph and someone might be able to help. – smichr Jun 06 '19 at 21:36
  • Your symbolic equation does not match the form of your equation once you have substituted in the numbers. Hopefully you are using a substitution dictionary like `reps = {tau_2: -0.28, K_m: 0.23, J_p1: 0.94, J_g2: -0.0062, J_m: 0.32, K_u: -0.79, J_p2 : 0.29, J_g1: 0.22}` -- those are just random values -- and then doing `omega_nf_eq.subs(reps)` to get the substitution done for you. After doing such a thing, you will end up (after feeding this expression to `unrad`) with a quartic equation in `J_u` but that is very easy for `real_roots` to manage. – smichr Jun 07 '19 at 02:39
  • I've solved using sympy.solveset(942.5 - omega_nf_eq.rhs, J_u). Faster than light. – cicciobombacannoniere Jun 12 '19 at 09:54
0

I have solved using :

sympy.solveset(942.5 - omega_nf_eq.rhs, J_u)

I link the sympy.solveset() docs : Now is pretty fast.