1

(0 odpowiedzi, napisanych Mathematica: numeryka i grafika)

Podczas składania funkcji
f[x_]:=x^2
dla liczb o ustawionej przez użytkownika precyzji (x0=N[1,prec]), precyzja wyniku spada bardzo gwałtownie. Jak obejść ten błąd? Przy precyzji maszyny ten błąd się nie kumuluje. Przykład poniżej:

s1 = N[1];
s2 = N[1, 10];
TableForm[
 N /@ Precision /@ NestList[#^2 &, #, 20] & /@ {s1, s2},
 TableDirections -> Row,
 TableHeadings -> {{"Precision[x0]=MachinePrecision", 
    "Precision[x0]=10"}}]
ListPlot[%]

2

(1 odpowiedzi, napisanych Mathematica: numeryka i grafika)

Rozwiązaniem problemu jest użycie RuleDelayed, czyli

t2 = Select[Table[i, {i, 1, 40}], OddQ];
ReplacePart[t2, i_ :> (1 + t2[[i]])/2]

Różnica polega na kolejności obliczania prawej strony zastąpienia.

Rule (i_->t2[[ i ]]) zwraca błąd ponieważ prawa strona jest obliczana dla i=i.
RuleDelayed (i_:>t2[[ i ]]) nie oblicza prawej strony. Różnica jest analogiczna jak dla Set (=) i SetDelayed (:=), w manualu Mathematici jest to szczegółowo opisane.

Eksperymentowałem z podstawianiem explicite wzorów na pierwiastki za funkcję Root (używając ToRadicals), ale z marnym skutkiem. Zupełnie nie rozumiem dlaczego Mathematica nie radzi sobie z uproszczeniem tego wyrażenia.

Jest to prosta wersja skryptu, i jako taka wymaga podania funkcji wszystkich argumentów tj. "D[f[x,y],x]". Po takiej modyfikacji zamiana zmiennych działa prawidłowo.

bieg = {#1 Cos[#2] &, #1  Sin[#2] &};
rules = ChangeVariable[f, {x, y}, bieg, {r, th}, g];
D[f[x,y], x] /. rules // Simplify

Wynikiem jest wtedy

-((Sin[th]*Derivative[0, 1][g][r, th])/r) + 
 Cos[th]*Derivative[1, 0][g][r, th]

Z pewnością wystarczy drobna modyfikacja by program działał dla funkcji także o dowolnych argumentach (np z parametrem), pomyślę jeszcze o tej poprawce.

Z drugiej strony, co jeśli użytkownik poda argumenty w sposób niespójny? np najpierw f[x,y,t] a w innym miejscu f[t,x,y]? Program może wykryć położenia argumentów i inteligentnie pozamieniać zmienne, ale co jeśli użytkownik niedoczyta dokumentacji a ma na myśli równanie funkcyjne (np f[x,y] w jednym miejscu a f[y,x] w drugim)? Skrypt bez zająknięcia poda wynik (fałszywy), z którego nieświadomy użytkownik będzie korzystał. Chyba lepiej trzymać się jednej konwencji i pisać wszystkie argumenty i na właściwych miejscach, a podany przez Pana przypadek dodać do kontroli błędów.

Podniże znajduje się zbiór funkcji pozwalający na zamianę zmiennych w prostych (niefunkcyjnych) wyrażeniach zawierających pochodne cząstkowe dowolnego rzędu i w dowolnym wymiarze. Jest to wersja pozbawiona kontroli błędów, więc proszę używać na własne ryzyko. Funkcja nie wymaga podawania funkcji odwrotnej do funkcji definiującej zamianę zmiennych. Sposób użycia to np:

spherical = {# Cos[#3] Sin[#2] &, # Sin[#3] Sin[#2] &, # Cos[#2] &};
rules=ChangeVariable[f, {x, y, z}, spherical, {r, th, ph}, g];
f[x, y, z] + D[f[x, y, z], x, x] + D[f[x, y, z], y, y] + 
    D[f[x, y, z], z, z] /.rules // Simplify

"f" to nagłówek funkcji zastępowanej, {x,y,z} to nazwy dotychczasowych zmiennych, w miejscu "spherical" należy wstawić listę funkcji definiujących stare współrzędne ({x(r,th,ph),y(r,th,ph),...}), dalej nazwy nowych współrzędnych (np. {r,th,phi}) a na końcu nagłówek funkcji nowych zmiennych (może to być właściwie dowolny niezdefiniowany symbol).

Clear[Grad, Jac, InvJacGrad, ChangeVariable];
Grad[f_, args_] := Module[{IdM},
  IdM = IdentityMatrix[Length[args]];
  Evaluate[Derivative[Sequence @@ #][f]] @@ args & /@ IdM
  ]
Jac[yOf_, args_] := Module[{IdM},
  IdM = IdentityMatrix[Length[args]];
  Transpose[
   Outer[Derivative[Sequence @@ #][#2] @@ args &, IdM, yOf, 1]]
  ]
InvJacGrad[f_, x_, args_, l_] := Module[{ytemp},
  If[Length[l] == 0, f,
   ytemp = Unique[args];
   Function[
    Evaluate[ytemp],
    Evaluate[(
       Transpose[Inverse[Jac[x, ytemp]]].Grad[
         InvJacGrad[f, x, ytemp, Rest[l]],
         ytemp](*Inverse Jacobian used!*)
       )[[First[l]]]
     ]]
   ]
  ]
ChangeVariable[f_, args_, xOf_, newArgs_, g_] := Module[{},
  {
   Derivative[l__][f] @@ args :> Module[{k = List[l]},
     InvJacGrad[g, xOf, newArgs,
       Flatten[ConstantArray[#, k[[#]]] &
         /@ Range[Length[k]]]] @@ newArgs
     ],
   f @@ args -> g @@ newArgs,
   Sequence @@ Thread[args -> Through[xOf @@ newArgs]]
   }
  ]

Problem polega na tym, że nie istnieją wzory na ogólny pierwiastek wielomianu trzeciego stopnia, których Mathematica potrzebuje by uprościć wynik (chodzi o funkcje Root[], Wygląda na to, że odpowiednie wybranie stałych całkowania (C[1] i C[2]) mogłoby umożliwić przedstawienie rozwiązania w bardziej zwartej formie (C[1] występuje jako współczynnik wielomianu wewnątrz funkcji Root[], więc można by ją tak dobrać, by obniżyć stopień lub strywializować jeden z pierwiastków, a wtedy pozostanie do rozwiązania jedynie wielomian drugiego rzędu).

Poza tym według Wikipedii to równanie ma jawne rozwiązanie dla n=5 - http://en.wikipedia.org/wiki/Lane%E2%80 … n_equation .