A estranha atração do mapa logístico

21

O objetivo do desafio é plotar aproximadamente o atrator do mapa logístico em função de seu parâmetro r (também chamado de diagrama de bifurcação ) ou uma sub-região dele. A aparência do gráfico pode ser vista na seguinte imagem da Wikipedia:

insira a descrição da imagem aqui

fundo

O mapa logístico é uma função matemática que pega uma entrada x k e mapeia para uma saída x k + 1 definida como

             x k + 1 = r x k (1− x k )

onde r é o parâmetro do mapa, assumido como estando no intervalo [0, 4].

Dado r em [0,4], e um valor inicial x 0 , no intervalo [0,1], é interessante para aplicar repetidamente a função para um grande número N de iterações, produzindo um valor final x N . Note que x N também estará necessariamente em [0,1].

Como exemplo, considere r = 3,2, N = 1000. O valor inicial x 0 = 0,01 fornece x 1000 = 0,5130. Para x 0 = 0,02, o resultado é x 0 = 0,7995. Para quaisquer outros valores iniciais x 0, os valores finais x 1000 são extremamente próximos de 0,5130 ou 0,7995. Isso é visto no gráfico como a altura das duas linhas na posição horizontal r = 3,2.

Isso não significa que para r = 3,2 cada sequência converja para um desses dois valores. De fato, para os dois valores iniciais considerados acima, as seqüências são (observe o comportamento oscilante):

             x 0 = 0,01, ..., x 1000 = 0,5130, x 1001 = 0,7995, x 1002 = 0,5130, ...
             x 0 = 0,02, ..., x 1000 = 0,7995, x 1001 = 0,5130, x 1002 = 0,7995 , ...

O que é verdade é que, para N suficientemente grande e para quase todos os valores iniciais x 0 , o termo x N estará próximo de um dos elementos do conjunto {0,5130, 0,7995}. Esse conjunto é chamado de atrator para esse r específico .

Para outros valores do parâmetro r, o tamanho do conjunto de atratores, ou seus elementos, será alterado. O gráfico plota os elementos no atrator para cada r .

O atrator para um r específico pode ser estimado por

  1. testando uma ampla gama de valores iniciais x 0 ;
  2. deixando o sistema evoluir para um grande número de N de iterações; e
  3. tomando nota dos valores finais x N que são obtidos.

O desafio

Entradas

  • N : número de iterações.

  • r 1 , r 2 e s . Eles definem o conjunto R de valores de r , ou seja, R = { r 1 , r 1 + s , r 1 + 2 s , ..., r 2 }.

Procedimento

O conjunto X dos valores iniciais x 0 é fixo: X = {0,01, 0,02, ..., 0,99}. Opcionalmente, 0 e 1, podem também ser incluídos em X .

Para cada R , em R e cada x 0 em X , iterar a logística mapa N vezes para produzir x N . Registre as tuplas obtidas ( r , x N ).

Saída

Traçar cada tuplo ( r , x N ) como um ponto no plano com R como eixo horizontal e x N como eixo vertical. A saída deve ser gráfica (não arte ASCII).

Regras adicionais

  • O procedimento indicado define o resultado necessário, mas não é imposto. Qualquer outro procedimento que processe o mesmo conjunto de tuplas ( r , x N ) pode ser usado.
  • A entrada é flexível, como de costume.
  • Erros de ponto flutuante não serão mantidos contra o atendedor.
  • A saída gráfica é necessária, em qualquer um dos formatos aceitos . Em particular, a saída pode ser exibida na tela ou um arquivo gráfico pode ser produzido ou uma matriz de valores RGB pode ser impressa. Se estiver produzindo um arquivo ou uma matriz, publique um exemplo da aparência quando exibida.
  • Os gráficos podem ser vetoriais ou rasterizados. Para gráficos raster, o tamanho da imagem deve ser de pelo menos 400 × 400 pixels.
  • Cada ponto deve ser mostrado como um único pixel ou como uma marca com tamanho da ordem de um pixel (caso contrário, o gráfico fica rapidamente confuso).
  • A faixa do eixo deve ser [0,4] para r (eixo horizontal) e [0,1] para x N (eixo vertical); ou pode ser menor, desde que inclua todos os pontos obtidos.
  • As escalas de eixo são arbitrárias. Em particular, a escala não precisa ser a mesma para os dois eixos.
  • Linhas de grade, rótulos de eixos, cores e elementos semelhantes são aceitáveis, mas não são necessários.
  • O menor código em bytes vence.

Casos de teste

Clique em cada imagem para obter uma versão de alta resolução.

N = 1000; r1 = 2.4; r2 = 4; s = 0.001;

insira a descrição da imagem aqui

N = 2000; r1 = 3.4; r2 = 3.8; s = 0.0002;

insira a descrição da imagem aqui

N = 10000; r1 = 3.56; r2 = 3.59; s = 0.00002;

insira a descrição da imagem aqui

Reconhecimento

Obrigado a @FryAmTheEggman e @AndrasDeak por seus comentários úteis enquanto o desafio estava na caixa de areia.

Luis Mendo
fonte
O que não há solução python ?!
@Lembik Eu tenho uma implementação de referência em Python (e em Matlab), mas eu não quero me responder
Luis Mendo
Você pode responder suas próprias perguntas sobre PPCG (talvez surpreendentemente).
@Lembik Eu sei, mas eu prefiro ter respostas dos outros
Luis Mendo

Respostas:

13

MATL, 32 30 28 27 bytes

4 bytes salvos graças a @Luis

3$:0:.01:1!i:"tU-y*]'.'3$XG

O formato de entrada é r1, s, r2, eN

Experimente no MATL Online

insira a descrição da imagem aqui

Explicação

        % Implicitly grab the first three inputs
3$:     % Take these three inputs and create the array [r1, r1+s, ...]
0:.01:1 % [0, 0.01, 0.02, ... 1]
!       % Transpose this array
i       % Implicitly grab the input, N
:"      % For each iteration
  tU    % Duplicate and square the X matrix
  -     % Subtract from the X matrix (X - X^2) == x * (1 - x)
  y     % Make a copy of R array
  *     % Multiply the R array by the (X - X^2) matrix to yield the new X matrix
]       % End of for loop
'.'    % Push the string literal '.' to the stack (specifies that we want
        % dots as markers)
3$XG    % Call the 3-input version of PLOT to create the dot plot
Suever
fonte
8

Mathematica, 65 bytes

Graphics@Table[Point@{r,Nest[r#(1-#)&,x,#]},{x,0,1,.01},{r,##2}]&

Função pura tomando os argumentos N, r1, r2, s nessa ordem. Nest[r#(1-#)&,x,N]itera a função logística r#(1-#)&um total de Nvezes, começando em x; aqui o primeiro argumento para a função ( #) é o Nem questão; Point@{r,...}produz um Pointque Graphicsterá prazer em traçar. Table[...,{x,0,1,.01},{r,##2}]cria um monte desses pontos, com o xvalor variando de 0para 1em incrementos de .01; o ##2in {r,##2}denota todos os argumentos da função original a partir do segundo e, portanto, se {r,##2}expande para o {r,r1,r2,s}qual define corretamente o intervalo e o incremento para r.

Saída de amostra, no segundo caso de teste: a entrada

Graphics@Table[Point@{r,Nest[r#(1-#)&,x,#]},{x,0,1,.01},{r,##2}]&[2000,3.4,3.8,0.0002]

produz os gráficos abaixo.

insira a descrição da imagem aqui

Greg Martin
fonte
1
59 bytes ListPlot @ Table [{r, Nest [r # (1 - #) &, x, #]}, {x, 0,1, .01}, {r, ## 2}] &
J42161217
Esclareci no desafio que o procedimento indicado se destina a definir o resultado necessário, mas o procedimento em si não é imposto. Você pode usar qualquer outro procedimento que dê o mesmo resultado. Desculpe se isso não ficou claro no início
Luis Mendo
Não é um problema, temos várias boas respostas!
Greg Martin
1
Você não vai usar esses -6 bytes. Você acha que algo está errado com esta solução?
J42161217
Oh, eu pensei que sua resposta foi a publicação de (a versão) o código do seu comentário ....
Greg Martin
5

Mathematica, 65 bytes

Usei alguns truques de Greg Martin e esta é a minha versão sem usar gráficos

ListPlot@Table[{r,NestList[#(1-#)r&,.5,#][[-i]]},{i,99},{r,##2}]&

entrada

[1000, 2,4, 4, 0,001]

saída

insira a descrição da imagem aqui

entrada

[2000, 3,4, 3,8, 0,0002]

saída

insira a descrição da imagem aqui

J42161217
fonte
1
Primeiro resposta que escolhe para evitar valores iniciais 0 ou 1 e o ( x linha = 0 geram) :-)
Luis Mendo
Você deve adicionar uma explicação do que seu código faz, pois ele não segue o procedimento especificado. O OP pode decidir se o resultado de aparência precisa justifica o método alternativo.
Greg Martin
O procedimento especificado não é imposto. Tudo o que der o mesmo resultado, por qualquer outro meio, é permitido (vou esclarecer isso). Independentemente disso, estou curioso para ver a explicação.
Luis Mendo
Os pontos que você precisa plotar para cada r já existem em cada "Ninho". esse é um código original e foi minha primeira abordagem (há um tempo atrás) na plotagem desse diagrama.
J42161217
@Luis Mendo Eu tenho uma versão ainda mais curta (que faz um registro para o mathematica) .58 bytes, mas você deve inserir apenas 3 entradas [N, r1, r2]. Leva tempo, mas funciona. Plot [Table [NestList [# ( 1 - #) r &, 5, #] [[- i]], {i, 99}], {r, ## 2}] &
J42161217
2

TI-Básico, 85 bytes

Prompt P,Q,S,N
P→Xmin:Q→Xmax
0→Ymin:1→Ymax
For(W,.01,1,.01
For(R,P,Q,S
W→X
For(U,1,N
R*X*(1-X→X
End
Pt-On(R,X
End
End

Um programa TI-Basic completo que recebe as entradas na ordem r1,r2,s,Ne mostra a saída em tempo real na tela do gráfico. Observe que isso tende a ser incrivelmente lento .

Aqui está uma amostra de saída incompleta gerada após cerca de 2,5 horas para a entrada 3,4,0.01,100:

insira a descrição da imagem aqui

R. Kap
fonte
Você não precisa dos *sinais.
precisa saber é o seguinte
1

ProcessandoJS, 125 123 120 bytes

Obrigado ao Kritixi Lithos por salvar 3 bytes.

var f(n,q,r,s){size(4e3,1e3);for(i=0;i<1;i+=.01)for(p=q;p<=r;p+=s){x=i;for(j=0;j<n;j++)x*=p-p*x;point(p*1e3,1e3-x*1e3)}}

Experimente online! Ligue usandof(N, r_1, r_2, s);

Somente ASCII
fonte
Eu acho que você pode substituir voidcom varporque de Processamento JS
Kritixi Lithos
E x*=p*(1-x)pode se tornarx*=p-p*x
Kritixi Lithos
Reorganizando o loop for, fico var f(n,q,r,s){size(4e3,1e3);for(i=0;i<1;i+=.01)for(p=q;x=i,p<=r;point(p*1e3,1e3-x*1e3),p+=s)for(j=0;j<n;j++)x*=p-p*x;}em 119 bytes
Kritixi Lithos
1

GEL , 158 bytes

`(N,r,t,s)=(LinePlotWindow=[r,t,0,1];for i=r to t by s do(p=.;for w=0to 1by 0.01do(x=w;for a=0to N do(x=i*x*(1-x););p=[p;q=[i,x]];);LinePlotDrawPoints(p);););

Pode não ser o mais curto, mas atrai em tempo real, embora possa ser incrivelmente lento com entradas enormes. De qualquer forma, esta é uma função anônima que recebe entrada no formato (N,r1,r2,s)e exibe o gráfico em uma nova janela. Observe que isso deve ser executado com a versão GNOME do Genius.

Saída de amostra

R. Kap
fonte
1

R, 159 147 bytes

pryr::f({plot(NA,xlim=c(a,b),ylim=0:1);q=function(r,n,x=1:99/100){for(i in 1:n)x=r*x*(1-x);x};for(i in seq(a,b,s))points(rep(i,99),q(i,n),cex=.1)})

Qual produz a função

function (a, b, n, s) 
{
    plot(NA, xlim = c(a, b), ylim = 0:1)
    q = function(r, n, x = 1:99/100) {
        for (i in 1:n) x = r * x * (1 - x)
        x
    }
    for (i in seq(a, b, s)) points(rep(i, 99), q(i, n), cex = 0.1)
}

plot(NA,...)cria uma tela vazia com as dimensões corretas. qé a função que faz a iteração. Leva um valor de re, em seguida, faz niterações para todos os pontos de partida entre 0.01e 0.99. Em seguida, retorna o vetor resultante.

A for-loop aplica a função qpara a sequência apara bo passo s. Em vez de retornar os valores, ele os adiciona como pontos ao gráfico. Se o ponto de atração é um valor, todos os pontos se sobrepõem e são exibidos como um ponto. cex=.1é uma adição necessária para tornar os pontos o menor possível.

insira a descrição da imagem aqui

JAD
fonte