momente şi schiţe de informatică şi matematică
To attain knowledge, write. To attain wisdom, rewrite.

Fractalii Newton ai rădăcinilor unităţii - partea a treia

limbajul R
2017 jun

6. Determinarea ciclurilor de lungime 2

Pentru unele puncte, orbita constituită prin aplicarea iterativă a transformării $N()$ - vezi [2] - este un şir periodic. De exemplu, dacă $N(z_1)=z_2$ şi $z_2 \ne z_1$, iar $N(z_2)=N(N(z_1))=z_1$, atunci orbita lui $z_1$ ciclează la nesfârşit între $z_1$ şi $z_2$.

Ciclurile de lungime 2 ale lui $N()$ sunt date de acele soluţii ale ecuaţiei $N(N(z))=z$ pentru care $N(z) \ne z$; avem $N(z)=z \iff \dfrac{(n-1)z^n+1}{nz^{n-1}}=z \iff z^n=1$ - astfel că, notând $t=z^n$, ajungem la ecuaţia acestor cicluri astfel:

$$N(N(z))=z \iff \dfrac{(n-1)N^n(z)+1}{nN^{n-1}(z)}=z \iff$$ $$(n-1)((n-1)t+1)^n-n^2t((n-1)t+1)^{n-1}+n^nt^{n-1}=0,\,\text{cu}\,\,t=z^n\ne 1 \,\,\,\,(*)$$

Excluzând soluţia evidentă $t=1$, rămân $n-1$ rădăcini $t$; radicalul de ordin $n$ al fiecăreia dintre acestea conduce la câte $n$ valori $z\in \mathbb{C}$ - prin urmare, mulţimea care este invariantă faţă de transformarea $N(N(z))$ conţine exact $n(n-1)$ puncte (repartizate câte $n$ pe $n-1$ cercuri cu centrul în origine), iar ciclurile de lungime 2 au capetele în această mulţime.

Ecuaţia obţinută poate fi tratată analitic doar pentru $n \le 4$ (şi ar trebui considerat fiecare caz în parte); preferăm desigur să trasăm numeric coeficienţii din (*), urmând apoi să folosim polyroot() pentru a obţine rădăcinile ecuaţiei.

Modelăm întâi înmulţirea a două polinoame (profitând aici de outer() şi tapply()):

poly_mul <- function(p1, p2) {  # vectorii coeficienţilor polinoamelor
  m2 <- outer(p1, p2)  # cu parametrul implicit FUN="*" (produsele de coeficienţi)
  return(as.vector(tapply(m2, row(m2) + col(m2), sum)))  # coeficienţii lui p1*p2
}

outer(p1, p2) produce matricea m2 în care liniile sunt asociate coeficienţilor p1, iar coloanele - coeficienţilor p2; valorile acestei matrici sunt produsele coeficienţilor respectivi, m2[i][j] = p1[i]*p2[j]. Adunând valorile m2[i][j] pentru i + j = d =constant, obţinem coeficientul monomului de gradul d din produsul celor două polinoame - operaţie montată prin tapply(): clasează valorile după suma indicilor de linie şi coloană şi apoi aplică sum() pe fiecare categorie de valori.

Folosind poly_mul(), putem obţine uşor polinomul din (*):

poly_2cycles <- function(n) {
  P <- p1 <- c(1, n-1)  # 1 + (n-1)t
  for(i in 1:(n-2))  # construieşte (1 + (n-1)t)^(n-1)
    P <- poly_mul(P, p1)
  P <- poly_mul(P, c(n-1, 1-2*n))
  P[n] <- P[n] + n^n  # ajustează penultimul coeficient
  for(i in n:2)  # elimină rădăcina 1 (schema lui Horner)
    P[i] <- P[i] + P[i+1]
  return(-P[2:(n+1)])  # polinomul (*), divizat prin (t-1)*(-1)
}

În final, am eliminat rădăcina $t=1$ şi am schimbat semnele coeficienţilor (având în vedere că suma coeficienţilor din (*) este $0$, iar coeficientul dominant ${\small (n-1)(n-1)^n-n^2(n-1)^{n-1}}$ este negativ şi este cel mai mare în valoare absolută).

Pentru $n=3$ obţinem:

> poly_2cycles(3)
[1]  2  5  20

adică $20t^2+5t+2$, cu zerorile $t_{1,2}=\frac{-5}{40} \pm \frac{\sqrt{135}}{40}i=\small -0.125 \pm 0.075\sqrt{5}i$.
Rezolvând ecuaţiile $z^3=t_{1,2}$, vom găsi cele 6 puncte ale mulţimii care este invariantă faţă de transformarea $N(z)$ şi care conţine cele 3 cicluri de lungime 2:

Nw <- function(n) {  # transformarea lui Newton pentru z^n-1
  function(z) ((n-1)*z + z^(1-n)) / n
}  # v. [1]
> p3 <- poly_2cycles(3)
> t12 <- polyroot(p3)
> z_2c <- as.vector(sapply(t12, function(t) polyroot(c(-t, 0, 0, 1))))
> print(z_2c)  # rădăcinile ecuaţiilor z^3 = t12
[1]  0.5386087+0.4172045i -0.6306140+0.2578466i  0.0920053-0.6750510i
[4]  0.0920053+0.6750510i -0.6306140-0.2578466i  0.5386087-0.4172045i
> sapply(z_2c, Nw(3))  # verificăm că mulţimea acestora este invariantă faţă de N()
[1]  0.5386087-0.4172045i  0.0920053+0.6750510i -0.6306140-0.2578466i
[4] -0.6306140+0.2578466i  0.0920053-0.6750510i  0.5386087+0.4172045i

Pentru $\small n=4$ găsim polinomul $\small 189t^3+41t^2+23t+3$ (şi desigur, putem adapta imediat secvenţa de mai sus pentru a găsi mulţimea care conţine ciclurile de lungime 2); pentru $\small n=5$ găsim $\small 2304t^4+459t^3+299t^2+59t+4$, etc. Putem sintetiza secvenţa de mai sus într-o funcţie care primind $n$, returnează mulţimea corespunzătoare a ciclurilor de lungime 2:

get_invariant <- function(n) {
  p3 <- poly_2cycles(n)  # polinomul din (*), divizat prin -(t-1)
  t12 <- polyroot(p3)  # soluţiile ecuaţiei (*), exceptând t=1
  as.vector(sapply(t12, function(t) polyroot(c(-t, rep(0, n-1), 1))))  # z^n-t = 0
}  # returnează punctele "2-periodice" faţă de N()

S-ar pune problema dacă există şi puncte 3-periodice, adică puncte distincte $z_1, z_2, z_3$ pentru care $z_2=N(z_1),\, z_3=N(z_2)$, iar $N(z_3)=z_1$. Pe de altă parte, avem de văzut care este locul punctelor 2-periodice faţă de transformarea $N()$: în "mulţimea Fatou", sau în "mulţimea Julia"?

vezi Cărţile mele (de programare)

docerpro | Prev | Next