Meilleur moyen de résoudre une équation non linéaire paramétrique stupidement compliquée avec réduire / résoudre, etc.

Aug 16 2020

J'ai une fonction non linéaire paramétrique qui est littéralement un cauchemar. Je sais que les racines existent, elles sont réelles et les deux paramètres p,esont tous deux positifs. Ce que j'attendais de Mathematica était d'obtenir une solution (sous forme de racine, pas de forme fermée) mais même de laisser le programme fonctionner toute la nuit, j'ai abandonné. Je ne peux pas comprendre si c'est moi qui pose le problème de manière non efficace, c'est mon PC qui a besoin de mises à niveau sérieuses ou le problème qui est tout simplement trop difficile pour des méthodes comme Réduire ou Résoudre. Si le cas est le dernier, je suppose que je suis condamné ... Une allusion aux deux autres? Merci pour l'aide.

Mes tentatives et l'équation:

f[x_]:=(1/(8 (p^2+x^2)^3))p^2 (-2 p^6+p^5 (4-8 x)+2 p^3 (3-8 x) x^2-6 p x^4+p^4 (80000+2 x-9 x^2)+2 p^2 x (40000+60000 x+x^2-5 x^3)-3 x^3 (-80000+40000 x+x^3)+(4 Sqrt[10] e x (p^2+x^2)^2 (2 p x^3+p^4 (-2+3 x)+2 p^3 x (-3+4 x)+x^2 (-80000+40000 x+x^3)+2 p^2 (40000-60000 x-x^2+2 x^3)))/Sqrt[-e p^2 (-1+x) x^2 (p^2+x^2)^2 (-40000+p^2+2 p x+x^2)])
Reduce[f[x]==0 && x>=0 &&p>=0 && e>=0,x,Reals] (*stuck running*)
Solve[f[x]==0 && x>=0 &&p>=0 && e>=0,x,Reals] (*stuck running*)

Réponses

2 CarlWoll Aug 18 2020 at 04:31

Si vous êtes satisfait d'une réponse approximative, vous pouvez essayer d'utiliser NDSolveValue. Votre fonction:

f[x_] := (1/(8 (p^2+x^2)^3))p^2 (-2 p^6+p^5 (4-8 x)+2 p^3 (3-8 x) x^2-6 p x^4+p^4 (80000+2 x-9 x^2)+2 p^2 x (40000+60000 x+x^2-5 x^3)-3 x^3 (-80000+40000 x+x^3)+(4 Sqrt[10] e x (p^2+x^2)^2 (2 p x^3+p^4 (-2+3 x)+2 p^3 x (-3+4 x)+x^2 (-80000+40000 x+x^3)+2 p^2 (40000-60000 x-x^2+2 x^3)))/Sqrt[-e p^2 (-1+x) x^2 (p^2+x^2)^2 (-40000+p^2+2 p x+x^2)])

Pour l'utiliser NDSolveValue, nous devons connaître une condition aux limites. Par exemple, voici la valeur de xquand pest 1:

x1 = x /. Block[{p=1}, First @ Solve[f[x] == 0, x]]

Root [256006399839996 + 1023948800640 e + (255942399200020 - 3071999998080 e) # 1 + (511955203840004 + 2304217598880 e) # 1 ^ 2 + (1279846402079976 - 2048025605760 e) # 1 ^ 3 + (-960279984200093) (192129617159615 + 4095897593600 e) # 1 ^ 5 + (-384151995680463-512486397600 e) # 1 ^ 6 + (-3455678391520375 + 2047846414080 e) # 1 ^ 7 + (2880427177039798 - 332788480 + e) ​​1 5 7 899480 e) 1024102385280 e) # 1 ^ 9 + (43199520578 + 256064008800 e) # 1 ^ 10 + (-14402879738 - 25598080 e) # 1 ^ 11 + (-359829 + 12802240 e) # 1 ^ 12 + 360087 # 1 ^ 13 + ( 9 + 160 e) # 1 ^ 14 + 9 # 1 ^ 15 &, 1]

Maintenant, nous pouvons utiliser NDSolveValue:

sol = NDSolveValue[
    {
    D[f[x[p,e]]==0, p], x[1, e] == x1},
    x,
    {p,.1,100},
    {e,.1,10000},
    MaxStepFraction->.0005,
    PrecisionGoal->10
]; //AbsoluteTiming

{19.2292, Null}

Vérifiez quelques échantillons aléatoires:

Block[{p = 50, e = 200}, f[sol[p, e]]]
Block[{p = 10, e = 2000}, f[sol[p, e]]]

6,42413 * 10 ^ -9

8,0893 * 10 ^ -9

Visualisation:

Plot3D[sol[p,e], {p,.1,100}, {e,.1,10000}]

2 MichaelE2 Aug 17 2020 at 04:15

FWIW, voici une demi-solution: prenez le numérateur, rationalisez-le pour qu'il soit un polynôme et trouvez les racines. Il ne reste plus qu'à sélectionner ceux qui sont positifs quand pet qui esont positifs. Cette étape prend du temps (si cela peut être fait), sauf lorsque des valeurs numériques spécifiques sont données à pet e.

num = Simplify[ff, x >= 0 && p >= 0 && e >= 0] // Together // 
      Numerator // Simplify[#, x >= 0 && p >= 0 && e >= 0] & // 
    FactorList // #[[-1, 1]] & // Simplify
(*
p Sqrt[-e (-1 + x) (-40000 + p^2 + 2 p x + x^2)] (2 p^6 + 6 p x^4 + 
    p^5 (-4 + 8 x) + 2 p^3 x^2 (-3 + 8 x) + 
    p^4 (-80000 - 2 x + 9 x^2) + 3 x^3 (-80000 + 40000 x + x^3) + 
    2 p^2 x (-40000 - 60000 x - x^2 + 5 x^3)) - 
 4 Sqrt[10]
   e (2 p x^5 + 4 p^3 x^3 (-1 + 2 x) + p^6 (-2 + 3 x) + 
    2 p^5 x (-3 + 4 x) + p^2 x^3 (-80000 - 2 x + 5 x^2) + 
    x^4 (-80000 + 40000 x + x^3) + 
    p^4 (80000 - 120000 x - 4 x^2 + 7 x^3))
*)

Vérifiez qu'il existe deux termes (le premier contenant évidemment le radical):

Length@num
(*  2  *)
rat = num*MapAt[-# &, num, 1] // Expand // Simplify
(*
e (p^2 (-1 + x) (-40000 + p^2 + 2 p x + x^2) (2 p^6 + 6 p x^4 + 
      p^5 (-4 + 8 x) + 2 p^3 x^2 (-3 + 8 x) + 
      p^4 (-80000 - 2 x + 9 x^2) + 3 x^3 (-80000 + 40000 x + x^3) + 
      2 p^2 x (-40000 - 60000 x - x^2 + 5 x^3))^2 + 
   160 e (2 p x^5 + 4 p^3 x^3 (-1 + 2 x) + p^6 (-2 + 3 x) + 
      2 p^5 x (-3 + 4 x) + p^2 x^3 (-80000 - 2 x + 5 x^2) + 
      x^4 (-80000 + 40000 x + x^3) + 
      p^4 (80000 - 120000 x - 4 x^2 + 7 x^3))^2)
*)
Solve[rat == 0, x]
(* <15 Root objects> *)

Ces racines contiennent des solutions étrangères, et on semble devoir substituer des valeurs numériques pour que les paramètres fonctionnent avec elles. Si une telle approche est utile, il serait peut-être préférable de substituer les paramètres f[x]et de traiter l'équation résultante f[x] == 0.