2

I have a function named F which looks like below: a two-well line and I want to solve equations to get two common tangent points (x1,F(x1)),(x2,F(x2)) of the function in a given range. the equations are:

F'(x1)=F'(x2)                     (1)

F'(x1)*(x2-x1)+F(x1)=F(x2)        (2)

and to solve these two equations I use sympy.solve. the code is:

x1 = Symbol('x1')
F1 = F(x1,AA, BB, CC, DD, EE, FF)
x2 = Symbol('x2')
F2 = F(x2,AA, BB, CC, DD, EE, FF)
dF1=diff(F1,x1)
dF2=diff(F2,x2)
print(solve([dF1-dF2, F2-F1-dF1*(x2-x1)], [x1, x2]))

but the output is [(x2, x2)]. it is quite confusing me

  1. how to let x1 != x2?
  2. how to get numerical output of the equations?
  3. how to let x1 and x2 in a certain range?

I relly do not know what is wrong here, the following solve command has correct output:

print(solve([y1+y2,y1-y2], [x1, x2]))
print(solve(dy1, x1))

but not work for binary differential equations(1)and (2)

thankyou!

@Lutz Lehmann:Print out the functions and equations to see if they actually are what you think they are:

-------------------------------------------------------------------
In [1]: x1 = Symbol('x1')
F1 = F(x1,AA, BB, CC, DD, EE, FF)
x2 = Symbol('x2')
F2 = F(x2,AA, BB, CC, DD, EE, FF)
dF1=diff(F1,x1)
dF2=diff(F2,x2)
-------------------------------------------------------------------
In [2]: F1
Out[2]: 1⋅(1−1)(−0.05017656665838731+0.244249682150038(21−1)5−0.535363152864316(21−1)4+0.197479965493092(21−1)3−0.366799918223659(21−1)2−0.203076238788926)
(Out [2] is F)
-------------------------------------------------------------------
In [3]:dF1
Out[3]:1⋅(1−1)(−2.934399345789271+2.44249682150038(21−1)4−4.28290522291453(21−1)3+1.18487979295855(21−1)2+1.41702310623625)−1(−0.05017656665838731+0.244249682150038(21−1)5−0.535363152864316(21−1)4+0.197479965493092(21−1)3−0.366799918223659(21−1)2−0.203076238788926)+(1−1)(−0.05017656665838731+0.244249682150038(21−1)5−0.535363152864316(21−1)4+0.197479965493092(21−1)3−0.366799918223659(21−1)2−0.203076238788926)
(Out [3] is correct)
-------------------------------------------------------------------
In [4]:print(solve([F1+F2,F1-F2], [x1, x2]))  
[(0.0, 0.0), (0.0, 1.00000000000000), (1.00000000000000, 0.0), (1.00000000000000, 1.00000000000000)]
(In [4] is correct)
-------------------------------------------------------------------
In [5]:print(solve(dF1, x1))
[0.143449671600321, 0.462174698289538, 0.744534434388284, 1.44258261390104, 0.573315370351492 - 0.261496971434863*I, 0.573315370351492 + 0.261496971434863*I]
(In [5] is correct)
-------------------------------------------------------------------
In [6]:print(solve([dF1-dF2, F2-F1-dF1*(x2-x1)], [x1, x2]))
[(x2, x2)]
(In [6] is what I confused??)
-------------------------------------------------------------------
Er g Ch
  • 21
  • 7
  • `dsolve` is for solving differential equations, what you have however is a system of equations constructed from a function and its derivative. – Lutz Lehmann Jun 16 '22 at 04:06
  • thank you, but `solve` does not work either – Er g Ch Jun 16 '22 at 04:35
  • Print out the functions and equations to see if they actually are what you think they are. – Lutz Lehmann Jun 16 '22 at 05:09
  • Thank you, I add them in the main text :) – Er g Ch Jun 16 '22 at 05:37
  • You might want to add one artificial equation with one auxiliary variable to avoid the diagonal `Y*(x2-x1)==1`. This trick is even named, but I do not quite remember the name (Rabinovich?). – Lutz Lehmann Jun 16 '22 at 05:49
  • Looks like `F` is a degree 5 polynomial. So both equations are degree 5 polynomials in `x1` and `x2`. Have you tried replacing the floats with their symbolic form in the hope that the x^5 coefficients cancel? Eg. instead of `np.sqrt(2)`, you use `sym.sqrt(S(2))` where applicable. Or worst case, `F(x1, Symbol("AA"), Symbol("BB"), ...)`. – Chris du Plessis Jun 16 '22 at 18:24

0 Answers0