Карта возврата Пуанкаре для логистической карты f_c (x) = cx (1-x)

Question

Предположим:

f[c_,x_]:= c x (1-x) 

Фиксировать интервал (a, b) в [0,1].Карта возврата Пуанкаре (к интервалу (a, b)) функции f является функцией

R [c, x] = f ^ k [c, x],

где k - первая итерация, такая, что a<f^k[c,x]<b (здесь f^k[c,x] означает k кратную композицию f с самим собой, то есть f[c,f[c,...f[c,x]...]])

Итак, я хотел бынаписать функцию (или модуль)

 R[f_,a_,b_,n_,x_]

, которая рассматривает первые n итераций f и возвращает значение итерации f[c,x], которое сначала попадает в интервал [a,b].

Вот что я пытался:

R[f_[x___ ], a_,b_, n_, x0_] := 
Module[{i, y=x0},
    Catch[
        For[i = 0, i <= n, i++,
            If[a < f[{x}[[1]], y] < b,
                Throw[f[{x}[[1]], y]], y = f[{x}[[1]], y]
              ]
           ]
         ]
      ]

Код не работает, потому что там, где написано f[{x}[[1]], y], f понимается как умножение на число {x}[[1]], а не как логистическая функция, определенная выше.

Обратите внимание, что я ищу простой фрагмент кода, и, если возможно, не меняйте количество входов для функции R в вашем ответе.

РЕДАКТИРОВАТЬ: Я хотел бынравится называть R следующим образом.

R[f[3.5, t], 0.4, 0.7, 100, 0.2]

Это должно вернуть значение итерации x0=0.2 при первом попадании в интервал (0.4,0.7) после применения функции f[3.5,x]=3.5x(1-x).n - это максимальные итерации, которые мы пробуем перед тем, как отказаться.

Answer 1

Проблема в том, что вы ничего не возвращаете с Module. Чтобы исправить ваш код, как написано, я бы использовал

R[f_[x___], a_, b_, n_, x0_] :=
 Module[{i, y = x0},
  Catch[
   For[i = 0, i <= n, i++, Print["i= ", i];
    If[
     a < f[{x}[[1]], y] < b,
     y= f[{x}[[1]], y]; Throw[y],
     y = f[{x}[[1]], y]
     ](*If*)
    ](*For*)
   ];(*Catch*)
   y  (* returns y *)
  ](*Module*)

Однако это можно переписать более кратко, как

g[c_][x_] := c x(1 - x)
(* I used Q to differentiate it from R, above. *)
Q[f_, a_, b_, n_, x0_] :=
 Module[{i},
  (* -- FIXED THIS, See below. -- *)
  NestWhile[f, f[x0], Not[a < # < b]&, 1, n]
 ](*Module*)

Обратите внимание, я изменил метод вызова f с f[c,x] на g[c][x]. Основным преимуществом является возможность передачи g[c] в Q вместо f[c,t], где t - фиктивная переменная. Тогда он вызывается как

Q[g[3.5], 0.4, 0.7, 100, 0.2]

, который работает так же, как R.

Редактировать : Я смотрел на расширение вышеприведенного кода и заметил недостаток. У меня было условие a < f[#] < b&, которое говорит, что цикл будет продолжаться только в том случае, если значение значения следующей итерации находится в границах. Вместо этого мы хотим продолжить, только если текущая итерация выходит за пределы диапазона, поэтому я изменил ее на Not[ a < # < b ]&.

---

Что касается изменений, которые я рассматривал, то иногда с таким расчетом приятно иметь возможность просматривать полный список итераций. Для этого нам нужно внести несколько небольших изменений в приведенный выше код.

Clear[Q]
Options[Q] = {FullList -> False};
Q[f_, a_, b_, n_, x0_, opts : OptionsPattern[]] :=
  Module[{i, nst},
   nst = If[OptionValue[FullList], 
              NestWhileList, NestWhile];
   nst[f, f[x0], Not[a < # < b] &, 1, n]
  ](*Module*)

, который вводит Option FullList, который при значении True, Q будет использовать NestWhileList вместо NestWhile.