¿Cómo puedo obtener la solución de una ecuación complicada?

Introducción

Actualización 6/2014

Originalmente, encontré tres soluciones; ahora son seis. Pensé que sabía que había tres porque la racionalización de la ecuación da como resultado un polinomio de grado 24 y NSolve encuentra 24 raíces y luego elimina las raíces que no eran ceros de la ecuación original. Resulta que la ecuación era problemática desde el punto de vista de los números y me perdí tres.

Otras modificaciones:

• Parece más natural usar Rationalize para convertir los coeficientes a números exactos en este caso que SetPrecision.

• Asimismo, sustituí First @ eq por Subtract @@ eq.

Racionalizar la ecuación

Se debe tener en cuenta que de acuerdo con la documentación sobre NSolve

NSolve se ocupa principalmente de ecuaciones lineales y polinómicas.

La ecuación de OP es una ecuación algebraica que se puede convertir en una ecuación polinómica, que NSolve puede resolver fácilmente. Dado que NSolve no hace esto por sí solo, tenemos que racionalizar ” a mano.”La racionalización tiende a introducir soluciones extrañas, por lo que tenemos que comprobar las respuestas devueltas por NSolve. Los coeficientes algo desagradables de eq lo hacen más difícil, ya que al racionalizar eq, el error de redondeo conduce a que algunas de las raíces cuadradas no se cancelen como deberían. Podemos arreglarlo ajustando la precisión de eq a Infinity.

Podemos encontrar las raíces cuadradas en una expresión expr con

DeleteDuplicates @ Cases, Infinity]

Luego para cada raíz cuadrada, podemos multiplicar expr por la expresión con la raíz cuadrada reemplazada por su negativo y simplificar.

A continuación se muestra el resultado de aplicar el método a la función obtenida reemplazando Equal en eq por Subtract por Apply (@@). Podemos resolver la ecuación exactamente.

eqExact = Rationalize;rationalized = Simplify @ Fold &, eqExact, DeleteDuplicates@Cases, Infinity]];rootsRatExact = Solve;rootsRat = NSolve(* 24 roots {{ζ -> -24559.6 + 24559.6 I}, <<22>>, {ζ -> 24559.6 - 24559.7 I}}*)

Hay algunas diferencias significativas entre las dos soluciones:

(ζ /. rootsRat) - (ζ /. rootsRatExact) // Abs // Max(* 3.15226*10^-10*)

Selección de las raíces de la ecuación dada

Actualización: Aquí es donde me perdí algunos ceros.

Utilicé Pick para seleccionar raíces en las que la ecuación original eq tiene un valor pequeño, menor que algún umbral pequeño, digamos 10^-1. Podemos comprobar más tarde que cada una es una raíz, pero no me echa para otras raíces. Las raíces originales eran:

rootsEq = Pick < 10^-1 /. rootsRat](* {{ζ -> -3.78042 - 5.50655 I}, {ζ -> -3.20562 + 5.39914 I}, {ζ -> 6.98478 + 0.493405 I}}*)

Estos corresponden a raíces 7, 9 y 20 en rootsRat:

Position < 10^-1 &)](* {{7}, {9}, {20}}*)

Si compruebo la ecuación exacta de ceros con PossibleZeroQ, obtengo más raíces:

Position(* {{1}, {7}, {9}, {14}, {20}, {21}}*)

lleva una eternidad en estas raíces potenciales, o al menos más de lo que estaba dispuesto a esperar.]

Curiosamente, estas raíces adicionales están conectadas a las diferencias entre los resultados devueltos por NSolve y Solve:

Position, _?Positive](* {{1}, {2}, {3}, {14}, {21}, {22}, {23}, {24}}*)

Dejemos que nuestras nuevas raíces sean las siguientes y las podemos comprobar a continuación:

rootsEqNew = Pick];rootsEqNew // N(* {{ζ -> -24559.6 + 24559.6 I}, {ζ -> -3.78042 - 5.50655 I}, {ζ -> -3.20562 + 5.39914 I}, {ζ -> -0.0832786 - 722.827 I}, {ζ -> 6.98478 + 0.493405 I}, {ζ -> 722.642 - 0.100823 I}}*)

Comprueba

Al principio, at no se ve demasiado bien para las nuevas incorporaciones:

eqExact /. rootsEqNew // N(* { 0. + 9.7614*10^12 I, (* {1} *) -1.49012*10^-8 + 3.91155*10^-8 I, (* {7} *) -1.39698*10^-8 - 3.63216*10^-8 I, (* {9} *) -2. + 1536. I, (* {14} *) 1.49012*10^-8 + 1.11759*10^-8 I, (* {20} *) 0. + 1024. I} (* {21} *)*)

El problema está en evaluar la ecuación con la precisión de la máquina.

N(* N::meprec warning about $MaxExtraPrecision being reached *)(* {0``69.62973568978663 + 0``69.69302870899077 I, 0``90.5174054423328 + 0``90.55817837498498 I, 0``90.50250822096415 + 0``90.54414468499085 I, 0``80.1824915073549 + 0``79.76650578675965 I, 0``90.47483650216002 + 0``90.49782363232914 I, 0``80.17292602755023 + 0``79.76710897249409 I}*)

La advertencia no parece ser significativo. Subir $MaxExtraPrecision a 2000 sigue produciendo la advertencia y los ceros con precisiones de entre 2020 y 2040. Si son ceros, N probablemente siempre producirá la advertencia, ya que no puede tener un Precision (que no sea MachinePrecision).