Acessibilidade / Reportar erro

Caos determinístico: uma abordagem com o software Insight Maker focada em discussões sobre determinismo e previsibilidade

Deterministic Chaos: an approach with Insight Maker software focused on discussions about determinism and predictability

Resumos

A dinâmica não-linear, demonstrando a existência de eventos hipersensíveis às condições iniciais, sacudiu as bases da filosofia da Ciência, principalmente no que remete ao poder preditivo do conhecimento científico. No contexto educacional, no entanto, o caos determinístico continua sendo pouco explorado, mesmo em cursos de graduação, apesar do seu potencial para favorecer a aprendizagem tanto de conceitos científicos quanto de elementos centrais da filosofia da Ciência. Propondo uma alternativa para a abordagem do tema em cursos de formação de professores, apresentamos neste artigo uma abordagem didática para o ensino de tópicos de dinâmica não-linear utilizando um aplicativo online e open source chamado Insight Maker, que permite a integração numérica e geração de gráficos de forma automatizada, requerendo poucos conhecimentos de programação. Debatemos as noções de determinismo e previsibilidade analisando dois modelos teóricos de pêndulos: o pêndulo simples e o amortecido forçado. Na primeira situação, a integração numérica das equações de movimento pelo Insight Maker mostra que o fenômeno é em princípio previsível, já que suas soluções têm pouca sensibilidade às pequenas variações nas condições iniciais impostas. O mesmo não ocorreu com o pêndulo amortecido forçado. Utilizando a mesma ferramenta, mostramos que esse fenômeno pode apresentar regime caótico em função de variações de um parâmetro do sistema. Isso foi feito tanto pelo teste de sensibilidade às condições iniciais quanto pela obtenção indireta do expoente de Lyapunov. Evidenciamos as chamadas “ilhas de ordem” e demos particular enfoque à rota para o caos via cascata de duplicação do período, expondo conceitualmente o papel do número de Feigenbaum em sistemas não-lineares.

Palavras-chave:
Caos; Dinâmica não-linear; Previsibilidade; Simulação computacional; Insight Maker; Pêndulo amortecido forçado


Demonstrating the existence of events that are hypersensitive to initial conditions, non-linear dynamics shook the bases of philosophy of Science, mainly in what refers to the predictive power of scientific knowledge. In educational context, however, deterministic chaos continues to be little explored, even in undergraduate courses, despite its potential to favor the learning of both scientific knowledge and central elements of philosophy of science. Proposing an alternative approach to the subject in teacher training courses, this article presents a didactic approach for teaching nonlinear dynamics using an online and open source application called Insight Maker, which allows the solving of differential equations by numerical integration and automated generation of graphs, requiring little programming knowledge. We discuss the notions of causality, determinism and predictability by analyzing two theoretical models of pendulums: the simple pendulum and the forced damped one. In the first situation, numerical integration of the equations of motion in Insight Maker shows that the phenomenon is, in principle, predictable, since its solutions have little sensitivity to small variations in the imposed initial conditions. The same did not occur with the forced damped pendulum. Using the same tool, we show that this phenomenon can present chaotic behavior due to variations in a system parameter. This was done both by applying a sensitivity test considering variations in the initial conditions and the indirect calculation of the Lyapunov exponent. We showed the so-called “islands of order” and highlighted the route to chaos via the period doubling cascade, conceptually exposing the role of the Feigenbaum number in nonlinear systems.

Keywords:
Chaos; Nonlinear dynamics; Predictability; Computer simulation; Insight Maker; Forced damped pendulum


1. Introdução

Uma das principais características do conhecimento científico é a sua dinamicidade. Como um organismo vivo, a Ciência está em permanente transformação, amadurecendo a partir de intensos debates entre cientistas, fundamentados tanto em resultados experimentais e teorizações quanto em perspectivas filosóficas e tensões sociais. A despeito da vivacidade da Ciência, o ensino de Mecânica Clássica costuma ocupar longos períodos da educação formal nos quais é abordado como um conhecimento estático, assumindo-se, implícita ou explicitamente, que tal campo científico foi completamente estabilizado com os trabalhos de Galileu, Newton e Euler ainda no século XVIII, nada mais de interesse restando para analisar.

Os desdobramentos da Mecânica Clássica nos séculos XIX e XX, no entanto, foram fundamentais na evolução da Física e da Filosofia. O caos determinístico, que será explorado neste artigo, sacudiu as bases da Ciência, promovendo profundas reflexões sobre o poder preditivo dos modelos científicos. Conceitos relacionados com a Natureza da Ciência, como os de causalidade, determinismo e previsibilidade, estão no seio das discussões promovidas a partir dos avanços decorrentes desse campo de investigação. Cientes das potencialidades desse tema para instigar o imaginário popular, diferentes campos artísticos têm desenvolvido produções culturais pautadas por conceitos relacionados com o estudo do caos, como filmes (por exemplo, o “Efeito Borboleta”, premiada produção de Eric Bress e Mackye Gruber) e livros (como o “Problema de três corpos”, best-seller deCixin Liu).

O caos evidencia a Ciência como uma entidade viva, em permanente mudança, onde disciplinas se interconectam em suas mútuas influências. Foi investigando a estabilidade do sistema solar, no final do século XIX, que Henri Poincaré (1854–1912) reconheceu o caos determinístico. No entanto, foi com o desenvolvimento de outro campo, a informática, e a popularização dos computadores que a teoria do caos pôde alcançar avanços rápidos. Desde então, a dinâmica não-linear passou por importantes progressos, sendo popularizada pelos trabalhos de Edward Lorenz sobre a dinâmica atmosférica, cujos efeitos foram difundidos sob a alcunha de “efeito borboleta”. Recentemente, a medalha Fields concedida ao matemático brasileiro Artur Avila em 2014 pelo seu trabalho com sistemas dinâmicos atesta que a dinâmica não-linear continua sendo uma área aberta, valorizada e efervescente. Consolidando-se como um campo interdisciplinar, a importância dessa área vem crescendo por suas implicações em muitas áreas da Ciência além da física e da matemática, como economia, medicina e biologia [1[1] J.E. Skinner, M. Molnar, T. Vybiral e M. Mitra, Integr. Physiol. Behav. Sci. 27, 1 (1992)., 2[2] S. Vlad, P. Pascu e N. Morariu, Journal of computing 2, 1 (2010).].

Apesar da importância histórica e atual da dinâmica não-linear, ela é pouco explorada no ensino de Física, sendo em muitos casos omitida. A constatação desse fato vem de longa data: Gleick [3[3] J. Gleick, Caos: a criação de uma nova ciência (Campus, Rio de Janeiro, 1989).] já criticava a omissão do tópico de dinâmica não-linear nos cursos de formação científica na década de 70. Mas por que isso ocorre até hoje? De fato, um programa de estudos em dinâmica não-linear costuma demandar dos discentes certo domínio tanto de matemática quanto de programação. Ferrari et al. [4[4] P.C. Ferrari, J.A.P. Angotti e M.H.R. Tragtenberg, Ciência & Educação 15, 1 (2009).] constataram que esses fatores implicam dificuldades significativas para futuros docentes em Física. É possível que esta exigência contribua de forma significativa para a omissão desse tema nos cursos de graduação. No entanto, softwares modernos focados em modelagem computacional com fins didáticos possibilitam que as dificuldades relacionadas com o domínio do formalismo matemático e com os procedimentos de programação vinculados ao ensino de caos sejam significativamente amenizadas.

No presente artigo, procurando explorar tanto as potencialidades do uso de ferramentas computacionais no ensino de Física como oferecer uma alternativa para promover um ensino instigante, focado em temas que estimulam o imaginário das pessoas e centrado em conceitos basilares da filosofia da Ciência, apresentamos uma proposta de abordagem que permite a introdução de tópicos sobre caos determinístico em contextos educacionais. Acreditamos que a alternativa apresentada pode ser empregada no contexto de formação de professores, tanto no ensino público quanto privado. Para que isso seja possível, é necessário proporcionar ao futuro docente os conhecimentos que o tornam capaz da sua implementação. Com este intuito, conduziremos aqui um debate sobre o tema para uso em cursos de formação de professores de Física ou Ciências. A abordagem é centrada no estudo de um modelo da Física que possibilita a exploração de conceitos importantes no campo da dinâmica não-linear: o pêndulo amortecido forçado. Esse modelo é adequado para introduzir futuros docentes nessa área porque possibilita, com moderada complexidade matemática, o aprofundamento da aprendizagem de conhecimentos caros ao campo da Mecânica Clássica (e.g., relações entre forças e movimentos), assim como da aprendizagem de conceitos importantes no estudo de caos determinístico (e.g., rota para o caos, expoente de Lyapunov, número de Feigenbaum). Além disso, a resolução numérica tem interpretação gráfica intuitiva, vinculada a um evento real. O estudo quantitativo é feito por meio de um software chamado Insight Maker1 1 Disponível em: https://insightmaker.com/. Acesso em: 29/11/2022. , a partir do qual é possível resolver equações diferenciais de forma online com poucos conhecimentos de programação. Os gráficos que evidenciam a sensibilidade às condições iniciais e a duplicação do período na rota para o caos são gerados automaticamente.

Na seção seguinte, expomos considerações úteis para a abordagem proposta versando sobre os conceitos de causalidade, determinismo e previsibilidade. Para familiarizar o(a) leitor(a) com os conceitos e com o recurso computacional utilizado, fazemos isso a partir da exploração de um modelo conhecido mesmo entre os iniciantes em Física: o modelo de pêndulo simples. Na seção 3 3. Desafiando a Previsibilidade da Física: O Pêndulo Amortecido Forçado O pêndulo amortecido e forçado se diferencia do pêndulo simples pela inclusão de duas forças: i. uma relacionada com o arrasto devido ao ar, que assumimos ser sempre contrária à direção do movimento e diretamente proporcional à velocidade (F→a⁢r=-b⁢v→=-b⁢L⁢ϕ.⁢v^, onde b é uma constante de proporcionalidade que depende das características do corpo suspenso e do fluido em que está imerso), e ii. uma força periódica externa F→⁢(t), na direção do movimento. Admitimos que essa força é senoidal, com frequência angular ω e amplitude de força motriz F0, isto é, F→⁢(t)=F0⁢cos⁡ω⁢t⁢v^9. O esquema representando as forças é sintetizado na Figura 3. Figura 3: Esquema representando o pêndulo amortecido forçado e as forças nele atuantes. Observe que a direção da força externa pode agir tanto contra como a favor do movimento, já que sua frequência angular ω pode não ser ressonante com a oscilação do pêndulo. Pela Figura 3, constata-se que o torque em relação ao eixo de rotação é dado por τ→=[LF ⁢(t)-LFa⁢r-LP ⁢sin⁡ϕ]⁢j^. Considerando que I⁢ϕ¨=τ, onde I=m⁢L2 é o momento de inércia do sistema em relação ao eixo de sustentação do pêndulo, e que o torque resultante tem componentes apenas ao longo do eixo j^, a equação de movimento pode ser escrita, na forma escalar, como: (2) ϕ ¨ = 1 m ⁢ L ⁢ [ F ⁢ ( t ) - F a ⁢ r - P ⁢ sin ⁡ ϕ ] . Também vamos introduzir um parâmetro adimensional γ, definido como γ=F0/mg, que representa o quão intensa é a força externa em comparação com o peso do corpo. Ou seja, quando γ≪1, podemos dizer que a força externa é desprezível em relação ao peso do corpo; quando γ≈1, essa força é comparável com o peso do corpo. Esse é o parâmetro que vamos mudar nas simulações da rota para o caos descrita no artigo. A seguir, exploramos o modelo do pêndulo amortecido forçado resolvendo a Eq. (2) no Insight Maker. 3.1. Implementando o pêndulo amortecido forçado no Insight Maker e a constatação da imprevisibilidade A construção do diagrama de blocos no Insight Maker para o pêndulo amortecido forçado se dá pela inclusão de mais alguns blocos no diagrama usado para o pêndulo simples. Esses blocos são as novas forças atuantes, bem como os parâmetros a elas associados10. Para realizar o teste de sensibilidade, vamos considerar apenas variações no parâmetro ϕ, que dessa vez vai iniciar na vertical: ⟨ϕ⁢(0)⟩=0∘ e σϕ⁢(0)=0,5∘. Como estamos com intenções didáticas, ou seja, não é nossa prioridade termos precisão nos resultados, focando sobretudo em aspectos conceituais, optamos por rodar a simulação com menos exigências computacionais, limitando-nos à análise das implicações das incertezas apenas no ângulo inicial, e não nas outras grandezas do modelo. Esse procedimento será suficiente para as discussões conduzidas no artigo. Novamente, executamos o teste de sensibilidade realizando N=500 simulações. Na primeira análise, vamos considerar γ=0,2, isto é, a força motriz é pequena comparada com o peso do corpo. Os outros parâmetros escolhidos são: ϕ.⁢(0)=0rad/s, ω=2⁢π, m=1kg, g=9,759m/s2, L=g/(9⁢π2) e b=3⁢π/2. Na Figura 4 estão exibidos os resultados da simulação. Figura 4: (a) Resultado da simulação com o pêndulo amortecido forçado para γ=0,2 e (b) o correspondente teste de sensibilidade. Na Figura 4 percebemos que o sistema inicia uma fase de transiente, que leva cerca de dois segundos para se ajustar em uma trajetória periódica, com período de T=1s, medido diretamente a partir do gráfico. A análise de sensibilidade mostra que, mesmo variando o ângulo inicial levemente em torno da posição ϕ=0∘, o que pode ser entendido como uma incerteza no valor da posição inicial do pêndulo, não há variação significativa na predição do modelo. De fato, a predição parece ser melhor até que a do pêndulo simples. Na análise seguinte, repetimos a simulação, porém agora com γ=1,5. Esse valor implica que o valor máximo da força motriz é 50% maior do que o peso do corpo preso à extremidade da haste. Todos os outros parâmetros foram mantidos como antes. O resultado da simulação está expresso na Figura5. Figura 5: (a) Resultados da simulação com o pêndulo amortecido forçado para γ=1,5 e (b) o correspondente teste de sensibilidade. Começamos por observar um ponto fundamental que distingue as Figuras 4 e 5: enquanto a primeira exibe comportamento periódico, a segunda não o faz. Isso pode ser verificado aumentando o tempo de simulação para os parâmetros selecionados: não importa quanto tempo se observe, não é possível identificar um padrão. Observamos também que o modelo forneceu previsões idênticas até aproximadamente sete segundos mesmo com uma pequena variabilidade na condição inicial do ângulo. A partir daí, no entanto, os resultados divergiram tanto quanto mais afastados da média estavam inicialmente. Aparentemente, o parâmetro γ controla em algum nível a previsibilidade do sistema. De fato, as Figuras 4 e 5 evidenciam um ponto de transição da ordem para o caos que passamos a estudar em maiores detalhes agora. 3.2. A rota para o caos: duplicação do período Como visto, a condição necessária para a existência de soluções caóticas para uma equação diferencial são duas [17, p. 20]: a equação precisa conter termos não-lineares e envolver ao menos três variáveis dinâmicas (três dimensões). Na Eq. (2), o termo não-linear é sin⁡ϕ e as variáveis dinâmicas são ϕ, ϕ. e t11. Conforme mencionado, podemos identificar a rota para o caos pela escolha conveniente dos parâmetros. Cabe destacar, no entanto, que não há apenas uma rota para o caos no pêndulo amortecido e forçado. Algumas rotas encontradas até então são: a cascata de duplicação de período, tanto para soluções restritas e não-restritas12; as transições aleatórias entre dois estados restritos; a transição intermitente entre estados caóticos não-restritos com rotação alternada em sentido horário e anti-horário [27]. Neste artigo, vamos explorar a rota para o caos através da cascata de duplicação do período. Em 1975, o físico Mitchell J. Feigenbaum descobriu de forma empírica uma relação matemática genérica que relaciona a duplicação do período T de sistemas caóticos com um parâmetro ζ que o descreve. Se {T,2⁢T,22⁢T,…,2n-1⁢T}n=1∞ é o conjunto dos valores de período associados aos parâmetros {ζn}n=1∞, então se obtém o chamado número de Feigenbaum δ pela seguinte expressão: (3) δ = lim n → ∞ ⁡ ζ n - ζ n - 1 ζ n + 1 - ζ n ≈ 4 , 6692016 ⁢ … . Até o presente momento não sabemos explicar sua causa [28]. Sabemos, no entanto, que a cascata de duplicação do período ocorre respeitando essa lei para uma grande variedade de sistemas não-lineares, atestando a universalidade desse princípio. Quando trazido à tona em contextos de ensino, tal aspecto evidencia a contemporaneidade da dinâmica não-linear, bem como a importância de estudá-la para buscarmos uma melhor compreensão do mundo. Mais adiante, quando evidenciarmos a duplicação no período no pêndulo amortecido e forçado pelo Insight Maker, vamos obter um valor aproximado para o número de Feigenbaum. Iniciamos a rota para o caos buscando valores de γ maiores que γ=0,2. Utilizando as mesmas condições iniciais do exemplo anterior, aumentamos esse parâmetro sucessivamente até encontrar um valor limite a partir do qual são produzidas soluções periódicas com período duplicado (T=2s). Encontramos esse resultado para γ=1,0650 com uma precisão de 0,0001, isto é, diminuir o valor de γ em 0,0001 faz com que o período volte a ser de um segundo. Destacamos que esse é um valor próximo ao encontrado na literatura, que é de γ=1,0663[29, p. 474]. O resultado da simulação está exposto na Figura 6. Figura 6: Evolução do pêndulo amortecido forçado para γ=1,0650. A simulação foi executada até 40 segundos, buscando simultaneamente eliminar resíduos transientes na evolução temporal, bem como amenizar o erro associado aos arredondamentos realizados na integração numérica. Como no resultado apresentado na Figura 4, essa simulação apresenta uma fase transiente, que nesse caso dura até logo após 20 segundos. Depois disso, o sistema apresenta regime periódico, porém com período duplicado (T=2s) em relação à simulação com o parâmetro γ=0,2, quando o período era T=1s nas mesmas condições iniciais. Neste ponto, destacamos que o período não varia apenas em função de γ, tampouco apenas de forma a se duplicar. De fato, se a simulação fosse realizada com γ=1,077, seria observada uma evolução periódica com T=3s (período triplicado); ainda, para esse valor de γ, porém com ϕ⁢(0)=-π/2, obteríamos um período de T=1s. Por conta disso, é importante manter as condições iniciais fixas em cada simulação. Aumentando mais o parâmetro γ, obtemos a próxima duplicação do período em γ=1,0818, com a mesma precisão do exemplo anterior. O resultado é apresentado na Figura 7. Figura 7: Evolução do pêndulo amortecido forçado para γ=1,0818. A Figura 7 evidencia movimento periódico com período T=4s. Continuando esse processo, encontramos que, para γ=1,0825, o período é novamente duplicado (T=8s). Para γ=1,0827 temos mais uma duplicação de período (T=16s). Os resultados obtidos até agora estão condensados na Tabela 1. Tabela 1: Evolução do período do movimento em função do parâmetro γ. γ T (s) Δ ⁢ γ 1,0650 2 – 1,0818 4 0,0168 1,0825 8 0,0007 1,0827 16 0,0002 Pelo exposto, vê-se que o período duplica para incrementos Δ⁢γ cada vez menores do parâmetro γ, indicando que essa sequência é convergente. Os resultados da Tabela 1 podem ser usados para fazer uma primeira aproximação do número de Feigenbaum pela Eq. (3), apesar de que, a rigor, a relação só é valida para {γn};n→∞. Fazendo a conta para os três últimos valores desse parâmetro, obtemos (4) δ = 1 , 0825 - 1 , 0818 1 , 0827 - 1 , 0825 = 0 , 007 0 , 0002 = 3 , 5 . Apesar da falta de precisão, essa estimativa concreta ajuda na compreensão do significado da constante de Feigenbaum. Como a sequência de valores de Δ⁢γ é convergente, podemos concluir que, com a variação de γ, estamos duplicando o período do pêndulo diversas vezes, até um ponto em que Δ⁢γ tende a zero, significando que fazemos infinitas duplicações do período em uma variação Δ⁢γ tendendo a zero. Infinitas duplicações do período implicam em um período infinito, ou seja, caos. A literatura aponta que o limite da sequência da Tabela 1 ocorre para γc=1,0829, onde γc é o valor crítico ou limiar para a ocorrência de caos no sistema do pêndulo amortecido forçado [29, p. 475]13. Nesse caso, o período tem um valor infinito T=∞, de forma que o movimento não se repete. 3.3. A sensibilidade às condições iniciais em um coeficiente: o expoente de Lyapunov Vimos até aqui que sistemas caóticos possuem hipersensibilidade a variações nas condições iniciais, mesmo que essas variações sejam muito pequenas. Se compararmos a trajetória de dois sistemas idênticos que apenas diferem muito pouco nas suas condições iniciais, podemos quantificar sua divergência pelo expoente de Lyapunov, designado pelo símboloλ. Para compreendermos a lógica subjacente a esse coeficiente, começamos expressando sua definição em sistemas discretos. Um sistema dinâmico discreto é representado por uma função matemática que mapeia seu estado em certo instante para o seu estado no instante seguinte, em um conjunto discreto de valores possíveis. Um exemplo famoso é o mapa logístico, muito utilizado na modelagem de sistemas que descrevem a dinâmica populacional. Nesses casos, vamos supor que a diferença d entre as trajetórias dos sistemas comparados, a cada iteração n, seja proporcional à diferença na iteração imediatamente anterior, ou seja, |dn+1|=a⁢|dn|, onde a é uma constante. Desse modo, podemos escrever, por exemplo, |d2|=a⁢|d1|, |d1|=a⁢|d0|, |d2|=a2⁢|d0|, ou seja, |dn|=an⁢|d0|. Usando propriedades de logaritmos, podemos escrever |dn|=eλ⁢n⁢|d0|, onde o expoente de Lyapunov é dado por λ=ln⁡a[30, p. 358]. Um valor de λ>0 indica que o sistema é sensível a condições iniciais, pois as pequenas diferenças entre dois sistemas crescem exponencialmente, podendo haver a ocorrência de caos; caso λ≤0, o sistema pode ser entendido como estacionário ou periodicamente estável. A determinação do expoente de Lyapunov para sistemas contínuos descritos por equações diferenciais é mais complexa, especialmente quando a solução analítica não é conhecida, como no nosso caso. Nessa situação, o expoente de Lyapunov mede o comportamento assintótico médio da taxa de expansão/contração local dos eixos de uma hiperesfera infinitesimal de estados centrada nas condições iniciais do sistema [31, p. 245]. Está além do escopo deste artigo propor métodos de integração compatíveis com essa definição. Entretanto, como temos fins didáticos, podemos fazer uma aproximação com o caso discreto mencionado no parágrafo anterior. Tal aproximação nos permitirá determinar o sinal do expoente de Lyapunov, não seu valor propriamente dito, porém essa informação é suficiente para evidenciar a possibilidade de caos no sistema. Em analogia com o caso discreto, nossa definição tem a seguinte forma (já empregando as variáveis do pêndulo amortecido forçado): (5) Δ ⁢ ϕ ⁢ ( t ) = Δ ⁢ ϕ ⁢ ( 0 ) ⁢ e λ ⁢ t , onde Δ⁢ϕ⁢(t)=ϕ1⁢(t)-ϕ2⁢(t) é a diferença na posição angular de dois sistemas que diferem inicialmente pela quantidade pequena Δ⁢ϕ⁢(0). Linearizando a Eq. (5), obtemos: (6) x ⁢ ( t ) = λ ⁢ t + ln ⁡ ( Δ ⁢ ϕ ⁢ ( 0 ) ) , Onde x⁢(t)=ln⁡(ϕ1⁢(t)-ϕ2⁢(t)). Observe que, no formato da Eq. (6), λ é o coeficiente angular e ln⁡(Δ⁢ϕ⁢(0)), o coeficiente linear de uma reta. Nesse gráfico, uma inclinação positiva/negativa indicará que o expoente de Lyapunov é positivo/negativo. Como os sistemas estudados são altamente sensíveis às condições iniciais, esse método só tem sentido para pequenos intervalos de tempo, já que o próprio erro numérico decorrente dos arredondamentos realizados pelo computador pode ser entendido como uma pequena variação nas condições iniciais. Desse modo, o método que usaremos para avaliar λ não é adequado para pesquisas que demandam precisão. Existem métodos avançados de integração numérica usados por cientistas para resolver esse problema [31, 32]. Nosso objetivo aqui é apenas proporcionar uma compreensão conceitual dos efeitos desse expoente em situações de ensino de Física. Podemos usar o Insight Maker para simular as trajetórias ϕ1⁢(t) e ϕ2⁢(t) de dois sistemas que diferem apenas pelas condições iniciais duplicando o diagrama de blocos usados para estudar o pêndulo amortecido forçado14. Nossa primeira simulação foi feita usando as condições iniciais dos exemplos anteriores, porém com o parâmetro γ na região de caos, isto é: γ=1,2>γc. Além disso, usamos uma pequena variação de Δ⁢ϕ⁢(0)=0,25∘. Os resultados são mostrados na Figura 8. Figura 8: (a) Comparação das trajetórias ϕ1⁢(t) e ϕ2⁢(t) com γ=1,2; (b)x⁢(t) com inclinação positiva, indicando que o expoente de Lyapunov tem esse mesmo sinal. A Figura 8(a) indica que uma pequena divergência na condição inicial dos sistemas provocou diferenças significativas na previsão dos modelos em pouco tempo de simulação. Esse resultado era esperado porque o parâmetro γ é maior que o limite para haver caos. Se não soubéssemos esse fato de antemão, precisaríamos esperar cerca de 10 segundos de simulação para perceber que as trajetórias dos sistemas divergem. Entretanto, o gráfico na Figura 8(b) nos mostra desde o início da simulação que esse resultado aconteceria. Isso fica claro porque x⁢(t) tem inclinação positiva até esse instante, indicando que o coeficiente de Lyapunov é positivo. Vimos que a partir de γc=1,0829 o pêndulo amortecido passa a apresentar caos. Entretanto, isso não ocorre para todos os valores γ>γc. De fato, sabe-se que sistemas caóticos alternam intervalos de caos com outros de movimento periódico. Podemos evidenciar isso com o Insight Maker testando para os valores γ=1,5 e γ=1,6 e usando o restantes dos parâmetros como na simulação anterior. O resultado está na Figura 9. Figura 9: Simulações com γ=1,5 (a) e γ=1,6(b). A Figura 9(a) evidencia que, aumentando o parâmetro até γ=1,5, o sistema permanece apresentando caos, com expoente de Lyapunov positivo, como na Figura 8. Entretanto, aumentando o parâmetro mais um pouco, obtemos soluções periódicas (T=3s) estacionárias, de forma que o parâmetro γ=1,6 e seu entorno se constitui em uma “ilha de ordem”, como evidenciado pela Figura 9(b). Prosseguindo aumentando o parâmetro γ, seguem-se mais regiões alternadas de caos com regiões periódicas, que se tornam cada vez mais raras à medida que o parâmetro aumenta. , ampliamos a discussão explorando o pêndulo amortecido forçado no Insight Maker, evidenciando os problemas que surgem relativos à previsibilidade do modelo e exibindo sua rota para o caos. Na última seção tecemos as considerações finais.

2. As Noções de Causalidade, Determinismo e Previsibilidade: Uma Análise a Partir do Pêndulo Simples

O pêndulo tem um papel muito importante na história da física, assim como da Ciência e da cultura da sociedade ocidental [5[5] B.S. Hall, Annals of Science 35, 1 (1978)., 6[6] M.R. Matthews, in: International Handbook of Research in History, Philosophy and Science Teaching, editado por M.R. Matthews (Springer, Dordrecht, 2014).]. O estudo do pêndulo deu contribuições importantes tanto para o conhecimento humano, como para o estabelecimento de formas precisas de medir intervalos de tempo, o que na época das grandes navegações foi fundamental para a determinação da longitude das embarcações em alto-mar; também para a determinação experimental da aceleração da gravidade, mostrando que seu valor muda de acordo com a latitude no globo terrestre, que passa a ter formato reconhecidamente esferóide oblato; ainda, para o estudo das colisões e o consequente estabelecimento de leis de conservação na Física, consolidando a mecânica celeste newtoniana e o consequente fim da dicotomia entre o mundo sublunar e o supralunar [7[7] D.B. Meli, Thinking with Objects (The Johns Hopkins University Press, Baltimore, 2006)., p. 206].

Diversos notáveis cientistas se debruçaram sobre o estudo de pêndulos. Galileu, por exemplo, fortemente embebido no paradigma platônico, argumentava que as trajetórias dos pêndulos formavam segmentos de círculos – figura geométrica associada aos movimentos perfeitos. Por isso, dizia que os pêndulos eram isócronos e, a despeito dos seus resultados empíricos, assumia, equivocadamente, que seus períodos não dependiam da amplitude do movimento [6[6] M.R. Matthews, in: International Handbook of Research in History, Philosophy and Science Teaching, editado por M.R. Matthews (Springer, Dordrecht, 2014)., p. 21]. Foi apenas com Huygens, no entanto, que o período de pequenas oscilações do que hoje chamamos de pêndulo simples foi caracterizado pela equação T=2πL/g, onde L é o comprimento da haste e g a intensidade da aceleração da gravidade. As deduções de Huygens são fundamentadas em argumentos geométricos que hoje consideramos muito complicados (consulte-os em [8[8] Some Mathematical Works of the 17th & 18th Centuries, including Newton’s Principia, Euler’s Mechanica, Introductio in Analysin, etc., translated mainly from Latin into English, disponível em: http://www.17centurymaths.com/.
http://www.17centurymaths.com/...
]). Atualmente, no entanto, a dedução das equações associadas ao modelo de pêndulo simples é tradicionalmente apresentada para estudantes dos primeiros semestres dos cursos de Física sem significativas dificuldades.

Fenomenologicamente, o movimento de um pêndulo formado por uma pequena bolinha de massa m oscilando presa em uma haste é um dos mais simples de ser explicado. Trata-se de um evento que se repete, marcado por oscilações de períodos aproximadamente constantes. O modelo de pêndulo simples, no entanto, representa esse evento apenas de forma aproximada. Na dedução das suas predições, que não está no escopo deste artigo, consideramos que a haste que o sustenta é rígida e inextensível; que o corpo na sua extremidade é pontual; que o campo gravitacional tem módulo constante; e que a força de arrasto com o ar é desprezível. Com essas hipóteses, concluímos que a intensidade do torque total sobre o pêndulo simples em relação ao seu eixo de sustentação é dado apenas em função da força peso, de intensidade P=mg, que age sobre o corpo, de forma que a equação de movimento resulta em [9[9] D.K. Randall, Física – Uma Abordagem Estratégica (Bookman, Porto Alegre, 2009), v. 2, 2 ed.]2 2 Por termos fins didáticos, ao invés de suprimir o peso P, que desaparece da equação após a simplificação com a massa m no denominador da Eq. (1), optamos por considerá-lo de forma explícita, o que também será feito com as forças no modelo do pêndulo amortecido forçado. Tal decisão é tomada em função das intenções didáticas do artigo, para que, no momento em que formos construir os modelos no software Insigth Maker, as forças sejam apresentadas explicitamente, favorecendo a análise dos eventos do ponto de vista físico (avaliação das forças atuantes no sistema), e não puramente matemático (exclusivamente como taxas de variação de grandezas). :

(1) ϕ ¨ = - P m L sin ϕ ,

onde ϕ é o ângulo (crescendo no sentido anti-horário) entre a haste do pêndulo e a direção vertical. De acordo com o chamado determinismo clássico, que é a postura usualmente propagada em instituições de ensino, uma vez conhecida a posição e velocidade (ϕ.) em certo instante, podemos conhecer as mesmas grandezas num instante de tempo posterior (ou anterior) utilizando a Eq. (1). Ou seja, podemos prever o movimento futuro da bolinha. Mas podemos mesmo? É comum que, apressadamente, as pessoas respondam a essa pergunta considerando que, se podemos determinar a evolução de um sistema, podemos prevê-lo. Vamos aprofundar as diferenças entre esses aspectos neste artigo analisando as diferenças entre os conceitos de causalidade, determinismo e previsibilidade.

Historicamente, a noção de que sistemas físicos são previsíveis se desenvolveu a partir das ideias de causalidade e determinismo [10[10] M. Paty, Scientiae Studia 2, 1 (2004)., 11[11] M. Paty, Scientiae Studia 2, 4 (2004).]. A noção de causalidade é um dos princípios racionais do pensamento segundo o qual todo fenômeno tem uma causa, razão ou motivo. Essa concepção dá força à ideia de necessidade nas leis da natureza, as quais o pensamento racional do ser humano seria capaz de apreender. Nesse caso, o estado do pêndulo em um instante é entendido como a causa das suas variações de aceleração, ou seja, é o que causa suas oscilações, estando essas informações engendradas na própria forma da equação diferencial (Eq. (1)). Essa concepção de causalidade foi difundida no século XVIII com D’Alembert e Lagrange no contexto do desenvolvimento da mecânica analítica, tornando-se a concepção hegemônica na Física [12[12] J.R. Alembert, Alembert, J.R. d’ Alembert, em: Encyclopédie ou Dictionnaire raisonné des sciences, des arts et des métiers, editado por D. Diderot (Primeur – Libraire, Chez Pellet, 1752).].

O determinismo pode ser entendido como uma generalização da ideia de causalidade, referindo-se à ideia de fatalidade ou inexorabilidade dos eventos físicos. Também está vinculado com a crença na necessidade das leis naturais, implicando a capacidade de se determinar antecipadamente os eventos futuros. No caso do pêndulo, reforça-se a noção determinista de que seus movimentos futuros são plenamente estabelecidos pelo seu estado inicial e pelas leis de Newton. A noção atual de determinismo é decorrente da estabelecida por Laplace no século XIX. Segundo ele, uma inteligência que fosse capaz de conhecer todas as condições iniciais de um sistema poderia, em tese, deduzir seu passado e futuro, de forma que tudo é previsível e determinado. No entanto, Laplace não acreditava que o homem pudesse conhecer tais condições iniciais. De fato, toda medição científica é acompanhada de uma incerteza; por menor que seja a incerteza de uma medição, sempre caracterizamos o valor de uma grandeza como um intervalo, e não como um valor individual [13[13] J.P. Lima, M.T.X. Silva, F.L. Silveira e E.A. Veit, Laboratório de mecânica: subsídios para o ensino de Física experimental (Editora da UFRGS, Porto Alegre, 2013).]. Laplace sugere então, como paliativo a ignorância, o desenvolvimento da teoria das probabilidades [14[14] P.S. Laplace, Théorie analytique des probabilités (Imprimeur-Libraire, Paris, 1812).].

Hoje temos clareza de que a previsibilidade dos eventos físicos não é limitada apenas pela dificuldade em se conhecer as condições iniciais. Toda previsão da Física é obtida por meio de modelo(s) que descreve(m) a realidade de forma aproximada. De fato, para descrever a realidade precisamos necessariamente nos afastar dela [15[15] M. Bunge, Teoria e realidade (Perspectiva, São Paulo, 1974).]. Por exemplo, foi isso que fizemos quando admitimos que a haste que sustenta o pêndulo é rígida e inextensível e o corpo em sua extremidade é pontual. Em geral dizemos que um modelo é adequado/inadequado quando descreve certo fenômeno com margem de erro menor/maior que o grau de precisão admitido [16[16] L.A. Heidemann, I.S. Araujo e E.A. Veit, Investigações em Ensino de Ciências 23, 2, 2018.]. É ingênua, portanto, a noção de que os modelos da Física representam fenômenos com absoluta precisão para qualquer instante de tempo. Entretanto, mesmo se os modelos descrevessem a realidade de forma exata, teríamos outro fator que limita a previsibilidade.

Sistemas físicos descritos por equações diferenciais não-lineares podem apresentar comportamento caótico. Pequenas variações nas condições iniciais desses sistemas são suficientes para modificar totalmente seu comportamento futuro, o que compromete de forma irrevogável o resultado da previsão. Como nossa capacidade de determinar precisamente tais condições é limitada – como já havia percebido Laplace –, não podemos assumir que sistemas não-lineares sejam previsíveis (a priori). Existem duas condições necessárias para a existência de soluções caóticas de uma equação diferencial. A equação precisa conter termos não-lineares e envolver ao menos três variáveis dinâmicas (três dimensões) [17[17] M. Gitterman, The Chaotic Pendulum (World Scientific Publishing, New Jersey, 2010)., p. 20]. Nesse caso, a solução só pode ser obtida de forma numérica, e não há uma forma geral de determinar se um dado sistema apresentará ou não comportamento caótico.

Para compreender de forma intuitiva como um sistema determinista pode produzir resultados imprevisíveis, podemos fazer uma analogia com o conceito de equilíbrio mecânico: a imprevisibilidade de sistemas caóticos decorre do fato de que as soluções das equações diferenciais que os representam contém regiões de equilíbrio instável, tornando pequenas variações nas condições iniciais muito significativas [18[18] R. Duit e M. Komorek, Research in Science Education 27, 3 (1997).]. Essa analogia é útil para entender como a evolução de um sistema pode ser imprevisível; porém, o fenômeno do caos não se limita à mecânica, podendo ocorrer em diversos campos, como o da ótica. Por exemplo, nos lasers de diodo de Fabry-Perot (um tipo comum), a corrente elétrica deve ser modulada na frequência correta para garantir a emissão periódica de luz, isto é, fora da região caótica [19[19] S. Iezekiel, S. Bennett e C. M. Snowden, Workshop on High Performance Electron Devices for Microwave and Optoelectronic Applications 1, 1 (1997).].

Antes de explorarmos o pêndulo amortecido forçado, vamos analisar aqui as implicações do nosso desconhecimento sobre as condições iniciais em um sistema mais familiar aos estudantes: o pêndulo simples. Faremos isso com dois objetivos didáticos: i. familiarizar o(a) leitor(a) com o uso do software Insight Maker através da análise de um evento mais conhecido pelos(as) estudantes; e ii. introduzir os procedimentos relacionados com os conceitos de causalidade, determinismo e previsibilidade que abordaremos com esse software a partir de um sistema mais simples – que não apresenta comportamento caótico. Seremos dirigidos pela seguinte questão: Quais as implicações do desconhecimento das condições iniciais no caso do pêndulo simples? Guiados por questão semelhante, mostraremos, na seção 3 3. Desafiando a Previsibilidade da Física: O Pêndulo Amortecido Forçado O pêndulo amortecido e forçado se diferencia do pêndulo simples pela inclusão de duas forças: i. uma relacionada com o arrasto devido ao ar, que assumimos ser sempre contrária à direção do movimento e diretamente proporcional à velocidade (F→a⁢r=-b⁢v→=-b⁢L⁢ϕ.⁢v^, onde b é uma constante de proporcionalidade que depende das características do corpo suspenso e do fluido em que está imerso), e ii. uma força periódica externa F→⁢(t), na direção do movimento. Admitimos que essa força é senoidal, com frequência angular ω e amplitude de força motriz F0, isto é, F→⁢(t)=F0⁢cos⁡ω⁢t⁢v^9. O esquema representando as forças é sintetizado na Figura 3. Figura 3: Esquema representando o pêndulo amortecido forçado e as forças nele atuantes. Observe que a direção da força externa pode agir tanto contra como a favor do movimento, já que sua frequência angular ω pode não ser ressonante com a oscilação do pêndulo. Pela Figura 3, constata-se que o torque em relação ao eixo de rotação é dado por τ→=[LF ⁢(t)-LFa⁢r-LP ⁢sin⁡ϕ]⁢j^. Considerando que I⁢ϕ¨=τ, onde I=m⁢L2 é o momento de inércia do sistema em relação ao eixo de sustentação do pêndulo, e que o torque resultante tem componentes apenas ao longo do eixo j^, a equação de movimento pode ser escrita, na forma escalar, como: (2) ϕ ¨ = 1 m ⁢ L ⁢ [ F ⁢ ( t ) - F a ⁢ r - P ⁢ sin ⁡ ϕ ] . Também vamos introduzir um parâmetro adimensional γ, definido como γ=F0/mg, que representa o quão intensa é a força externa em comparação com o peso do corpo. Ou seja, quando γ≪1, podemos dizer que a força externa é desprezível em relação ao peso do corpo; quando γ≈1, essa força é comparável com o peso do corpo. Esse é o parâmetro que vamos mudar nas simulações da rota para o caos descrita no artigo. A seguir, exploramos o modelo do pêndulo amortecido forçado resolvendo a Eq. (2) no Insight Maker. 3.1. Implementando o pêndulo amortecido forçado no Insight Maker e a constatação da imprevisibilidade A construção do diagrama de blocos no Insight Maker para o pêndulo amortecido forçado se dá pela inclusão de mais alguns blocos no diagrama usado para o pêndulo simples. Esses blocos são as novas forças atuantes, bem como os parâmetros a elas associados10. Para realizar o teste de sensibilidade, vamos considerar apenas variações no parâmetro ϕ, que dessa vez vai iniciar na vertical: ⟨ϕ⁢(0)⟩=0∘ e σϕ⁢(0)=0,5∘. Como estamos com intenções didáticas, ou seja, não é nossa prioridade termos precisão nos resultados, focando sobretudo em aspectos conceituais, optamos por rodar a simulação com menos exigências computacionais, limitando-nos à análise das implicações das incertezas apenas no ângulo inicial, e não nas outras grandezas do modelo. Esse procedimento será suficiente para as discussões conduzidas no artigo. Novamente, executamos o teste de sensibilidade realizando N=500 simulações. Na primeira análise, vamos considerar γ=0,2, isto é, a força motriz é pequena comparada com o peso do corpo. Os outros parâmetros escolhidos são: ϕ.⁢(0)=0rad/s, ω=2⁢π, m=1kg, g=9,759m/s2, L=g/(9⁢π2) e b=3⁢π/2. Na Figura 4 estão exibidos os resultados da simulação. Figura 4: (a) Resultado da simulação com o pêndulo amortecido forçado para γ=0,2 e (b) o correspondente teste de sensibilidade. Na Figura 4 percebemos que o sistema inicia uma fase de transiente, que leva cerca de dois segundos para se ajustar em uma trajetória periódica, com período de T=1s, medido diretamente a partir do gráfico. A análise de sensibilidade mostra que, mesmo variando o ângulo inicial levemente em torno da posição ϕ=0∘, o que pode ser entendido como uma incerteza no valor da posição inicial do pêndulo, não há variação significativa na predição do modelo. De fato, a predição parece ser melhor até que a do pêndulo simples. Na análise seguinte, repetimos a simulação, porém agora com γ=1,5. Esse valor implica que o valor máximo da força motriz é 50% maior do que o peso do corpo preso à extremidade da haste. Todos os outros parâmetros foram mantidos como antes. O resultado da simulação está expresso na Figura5. Figura 5: (a) Resultados da simulação com o pêndulo amortecido forçado para γ=1,5 e (b) o correspondente teste de sensibilidade. Começamos por observar um ponto fundamental que distingue as Figuras 4 e 5: enquanto a primeira exibe comportamento periódico, a segunda não o faz. Isso pode ser verificado aumentando o tempo de simulação para os parâmetros selecionados: não importa quanto tempo se observe, não é possível identificar um padrão. Observamos também que o modelo forneceu previsões idênticas até aproximadamente sete segundos mesmo com uma pequena variabilidade na condição inicial do ângulo. A partir daí, no entanto, os resultados divergiram tanto quanto mais afastados da média estavam inicialmente. Aparentemente, o parâmetro γ controla em algum nível a previsibilidade do sistema. De fato, as Figuras 4 e 5 evidenciam um ponto de transição da ordem para o caos que passamos a estudar em maiores detalhes agora. 3.2. A rota para o caos: duplicação do período Como visto, a condição necessária para a existência de soluções caóticas para uma equação diferencial são duas [17, p. 20]: a equação precisa conter termos não-lineares e envolver ao menos três variáveis dinâmicas (três dimensões). Na Eq. (2), o termo não-linear é sin⁡ϕ e as variáveis dinâmicas são ϕ, ϕ. e t11. Conforme mencionado, podemos identificar a rota para o caos pela escolha conveniente dos parâmetros. Cabe destacar, no entanto, que não há apenas uma rota para o caos no pêndulo amortecido e forçado. Algumas rotas encontradas até então são: a cascata de duplicação de período, tanto para soluções restritas e não-restritas12; as transições aleatórias entre dois estados restritos; a transição intermitente entre estados caóticos não-restritos com rotação alternada em sentido horário e anti-horário [27]. Neste artigo, vamos explorar a rota para o caos através da cascata de duplicação do período. Em 1975, o físico Mitchell J. Feigenbaum descobriu de forma empírica uma relação matemática genérica que relaciona a duplicação do período T de sistemas caóticos com um parâmetro ζ que o descreve. Se {T,2⁢T,22⁢T,…,2n-1⁢T}n=1∞ é o conjunto dos valores de período associados aos parâmetros {ζn}n=1∞, então se obtém o chamado número de Feigenbaum δ pela seguinte expressão: (3) δ = lim n → ∞ ⁡ ζ n - ζ n - 1 ζ n + 1 - ζ n ≈ 4 , 6692016 ⁢ … . Até o presente momento não sabemos explicar sua causa [28]. Sabemos, no entanto, que a cascata de duplicação do período ocorre respeitando essa lei para uma grande variedade de sistemas não-lineares, atestando a universalidade desse princípio. Quando trazido à tona em contextos de ensino, tal aspecto evidencia a contemporaneidade da dinâmica não-linear, bem como a importância de estudá-la para buscarmos uma melhor compreensão do mundo. Mais adiante, quando evidenciarmos a duplicação no período no pêndulo amortecido e forçado pelo Insight Maker, vamos obter um valor aproximado para o número de Feigenbaum. Iniciamos a rota para o caos buscando valores de γ maiores que γ=0,2. Utilizando as mesmas condições iniciais do exemplo anterior, aumentamos esse parâmetro sucessivamente até encontrar um valor limite a partir do qual são produzidas soluções periódicas com período duplicado (T=2s). Encontramos esse resultado para γ=1,0650 com uma precisão de 0,0001, isto é, diminuir o valor de γ em 0,0001 faz com que o período volte a ser de um segundo. Destacamos que esse é um valor próximo ao encontrado na literatura, que é de γ=1,0663[29, p. 474]. O resultado da simulação está exposto na Figura 6. Figura 6: Evolução do pêndulo amortecido forçado para γ=1,0650. A simulação foi executada até 40 segundos, buscando simultaneamente eliminar resíduos transientes na evolução temporal, bem como amenizar o erro associado aos arredondamentos realizados na integração numérica. Como no resultado apresentado na Figura 4, essa simulação apresenta uma fase transiente, que nesse caso dura até logo após 20 segundos. Depois disso, o sistema apresenta regime periódico, porém com período duplicado (T=2s) em relação à simulação com o parâmetro γ=0,2, quando o período era T=1s nas mesmas condições iniciais. Neste ponto, destacamos que o período não varia apenas em função de γ, tampouco apenas de forma a se duplicar. De fato, se a simulação fosse realizada com γ=1,077, seria observada uma evolução periódica com T=3s (período triplicado); ainda, para esse valor de γ, porém com ϕ⁢(0)=-π/2, obteríamos um período de T=1s. Por conta disso, é importante manter as condições iniciais fixas em cada simulação. Aumentando mais o parâmetro γ, obtemos a próxima duplicação do período em γ=1,0818, com a mesma precisão do exemplo anterior. O resultado é apresentado na Figura 7. Figura 7: Evolução do pêndulo amortecido forçado para γ=1,0818. A Figura 7 evidencia movimento periódico com período T=4s. Continuando esse processo, encontramos que, para γ=1,0825, o período é novamente duplicado (T=8s). Para γ=1,0827 temos mais uma duplicação de período (T=16s). Os resultados obtidos até agora estão condensados na Tabela 1. Tabela 1: Evolução do período do movimento em função do parâmetro γ. γ T (s) Δ ⁢ γ 1,0650 2 – 1,0818 4 0,0168 1,0825 8 0,0007 1,0827 16 0,0002 Pelo exposto, vê-se que o período duplica para incrementos Δ⁢γ cada vez menores do parâmetro γ, indicando que essa sequência é convergente. Os resultados da Tabela 1 podem ser usados para fazer uma primeira aproximação do número de Feigenbaum pela Eq. (3), apesar de que, a rigor, a relação só é valida para {γn};n→∞. Fazendo a conta para os três últimos valores desse parâmetro, obtemos (4) δ = 1 , 0825 - 1 , 0818 1 , 0827 - 1 , 0825 = 0 , 007 0 , 0002 = 3 , 5 . Apesar da falta de precisão, essa estimativa concreta ajuda na compreensão do significado da constante de Feigenbaum. Como a sequência de valores de Δ⁢γ é convergente, podemos concluir que, com a variação de γ, estamos duplicando o período do pêndulo diversas vezes, até um ponto em que Δ⁢γ tende a zero, significando que fazemos infinitas duplicações do período em uma variação Δ⁢γ tendendo a zero. Infinitas duplicações do período implicam em um período infinito, ou seja, caos. A literatura aponta que o limite da sequência da Tabela 1 ocorre para γc=1,0829, onde γc é o valor crítico ou limiar para a ocorrência de caos no sistema do pêndulo amortecido forçado [29, p. 475]13. Nesse caso, o período tem um valor infinito T=∞, de forma que o movimento não se repete. 3.3. A sensibilidade às condições iniciais em um coeficiente: o expoente de Lyapunov Vimos até aqui que sistemas caóticos possuem hipersensibilidade a variações nas condições iniciais, mesmo que essas variações sejam muito pequenas. Se compararmos a trajetória de dois sistemas idênticos que apenas diferem muito pouco nas suas condições iniciais, podemos quantificar sua divergência pelo expoente de Lyapunov, designado pelo símboloλ. Para compreendermos a lógica subjacente a esse coeficiente, começamos expressando sua definição em sistemas discretos. Um sistema dinâmico discreto é representado por uma função matemática que mapeia seu estado em certo instante para o seu estado no instante seguinte, em um conjunto discreto de valores possíveis. Um exemplo famoso é o mapa logístico, muito utilizado na modelagem de sistemas que descrevem a dinâmica populacional. Nesses casos, vamos supor que a diferença d entre as trajetórias dos sistemas comparados, a cada iteração n, seja proporcional à diferença na iteração imediatamente anterior, ou seja, |dn+1|=a⁢|dn|, onde a é uma constante. Desse modo, podemos escrever, por exemplo, |d2|=a⁢|d1|, |d1|=a⁢|d0|, |d2|=a2⁢|d0|, ou seja, |dn|=an⁢|d0|. Usando propriedades de logaritmos, podemos escrever |dn|=eλ⁢n⁢|d0|, onde o expoente de Lyapunov é dado por λ=ln⁡a[30, p. 358]. Um valor de λ>0 indica que o sistema é sensível a condições iniciais, pois as pequenas diferenças entre dois sistemas crescem exponencialmente, podendo haver a ocorrência de caos; caso λ≤0, o sistema pode ser entendido como estacionário ou periodicamente estável. A determinação do expoente de Lyapunov para sistemas contínuos descritos por equações diferenciais é mais complexa, especialmente quando a solução analítica não é conhecida, como no nosso caso. Nessa situação, o expoente de Lyapunov mede o comportamento assintótico médio da taxa de expansão/contração local dos eixos de uma hiperesfera infinitesimal de estados centrada nas condições iniciais do sistema [31, p. 245]. Está além do escopo deste artigo propor métodos de integração compatíveis com essa definição. Entretanto, como temos fins didáticos, podemos fazer uma aproximação com o caso discreto mencionado no parágrafo anterior. Tal aproximação nos permitirá determinar o sinal do expoente de Lyapunov, não seu valor propriamente dito, porém essa informação é suficiente para evidenciar a possibilidade de caos no sistema. Em analogia com o caso discreto, nossa definição tem a seguinte forma (já empregando as variáveis do pêndulo amortecido forçado): (5) Δ ⁢ ϕ ⁢ ( t ) = Δ ⁢ ϕ ⁢ ( 0 ) ⁢ e λ ⁢ t , onde Δ⁢ϕ⁢(t)=ϕ1⁢(t)-ϕ2⁢(t) é a diferença na posição angular de dois sistemas que diferem inicialmente pela quantidade pequena Δ⁢ϕ⁢(0). Linearizando a Eq. (5), obtemos: (6) x ⁢ ( t ) = λ ⁢ t + ln ⁡ ( Δ ⁢ ϕ ⁢ ( 0 ) ) , Onde x⁢(t)=ln⁡(ϕ1⁢(t)-ϕ2⁢(t)). Observe que, no formato da Eq. (6), λ é o coeficiente angular e ln⁡(Δ⁢ϕ⁢(0)), o coeficiente linear de uma reta. Nesse gráfico, uma inclinação positiva/negativa indicará que o expoente de Lyapunov é positivo/negativo. Como os sistemas estudados são altamente sensíveis às condições iniciais, esse método só tem sentido para pequenos intervalos de tempo, já que o próprio erro numérico decorrente dos arredondamentos realizados pelo computador pode ser entendido como uma pequena variação nas condições iniciais. Desse modo, o método que usaremos para avaliar λ não é adequado para pesquisas que demandam precisão. Existem métodos avançados de integração numérica usados por cientistas para resolver esse problema [31, 32]. Nosso objetivo aqui é apenas proporcionar uma compreensão conceitual dos efeitos desse expoente em situações de ensino de Física. Podemos usar o Insight Maker para simular as trajetórias ϕ1⁢(t) e ϕ2⁢(t) de dois sistemas que diferem apenas pelas condições iniciais duplicando o diagrama de blocos usados para estudar o pêndulo amortecido forçado14. Nossa primeira simulação foi feita usando as condições iniciais dos exemplos anteriores, porém com o parâmetro γ na região de caos, isto é: γ=1,2>γc. Além disso, usamos uma pequena variação de Δ⁢ϕ⁢(0)=0,25∘. Os resultados são mostrados na Figura 8. Figura 8: (a) Comparação das trajetórias ϕ1⁢(t) e ϕ2⁢(t) com γ=1,2; (b)x⁢(t) com inclinação positiva, indicando que o expoente de Lyapunov tem esse mesmo sinal. A Figura 8(a) indica que uma pequena divergência na condição inicial dos sistemas provocou diferenças significativas na previsão dos modelos em pouco tempo de simulação. Esse resultado era esperado porque o parâmetro γ é maior que o limite para haver caos. Se não soubéssemos esse fato de antemão, precisaríamos esperar cerca de 10 segundos de simulação para perceber que as trajetórias dos sistemas divergem. Entretanto, o gráfico na Figura 8(b) nos mostra desde o início da simulação que esse resultado aconteceria. Isso fica claro porque x⁢(t) tem inclinação positiva até esse instante, indicando que o coeficiente de Lyapunov é positivo. Vimos que a partir de γc=1,0829 o pêndulo amortecido passa a apresentar caos. Entretanto, isso não ocorre para todos os valores γ>γc. De fato, sabe-se que sistemas caóticos alternam intervalos de caos com outros de movimento periódico. Podemos evidenciar isso com o Insight Maker testando para os valores γ=1,5 e γ=1,6 e usando o restantes dos parâmetros como na simulação anterior. O resultado está na Figura 9. Figura 9: Simulações com γ=1,5 (a) e γ=1,6(b). A Figura 9(a) evidencia que, aumentando o parâmetro até γ=1,5, o sistema permanece apresentando caos, com expoente de Lyapunov positivo, como na Figura 8. Entretanto, aumentando o parâmetro mais um pouco, obtemos soluções periódicas (T=3s) estacionárias, de forma que o parâmetro γ=1,6 e seu entorno se constitui em uma “ilha de ordem”, como evidenciado pela Figura 9(b). Prosseguindo aumentando o parâmetro γ, seguem-se mais regiões alternadas de caos com regiões periódicas, que se tornam cada vez mais raras à medida que o parâmetro aumenta. , que as implicações desse desconhecimento para o pêndulo amortecido forçado são muito distintas das para o pêndulo simples. De partida, podemos afirmar que o pêndulo simples não apresenta caos já que, apesar de a equação diferencial que o representa possuir termos não-lineares, contém apenas duas variáveis dinâmicas (ϕeϕ.). Já o pêndulo amortecido forçado – um sistema levemente modificado em relação ao pêndulo simples, como veremos – apresenta esse comportamento. Na próxima subseção, continuaremos estudando o pêndulo simples com o objetivo de consolidar nosso entendimento sobre o conceito de causalidade e determinismo, bem como introduzir o software Insight Maker, que será usado para explorar o caos determinista no pêndulo amortecido forçado.

2.1. Implementando o modelo de pêndulo simples no Insight Maker

O Insight Maker é um software online, gratuito e de código aberto utilizado para construir e compartilhar modelos que descrevem sistemas dinâmicos [20[20] S. Fortmann, Simulation Modelling Practice and Theory 47, 1 (2014).]. Os modelos são construídos graficamente, de forma que as relações causais entre os elementos do sistema são visualmente explicitadas. A apresentação do modelo é facilitada por um sistema de diagramação em blocos que separa os elementos do sistema dinâmico em variáveis, taxas de variação e constantes. A solução numérica da simulação computacional pode ser feita pelo método de Euler ou pelo método de Runge-Kutta de quarta ordem, com passo ajustável pelo usuário. Uma vantagem da utilização do Insight Maker é que os estudantes podem explorar modelos científicos mesmo dispondo de poucos conhecimentos de programação, demandando apenas domínio sobre o conceito de derivada, relacionando-o com taxas de variação.

Não está no escopo deste artigo expor um tutorial sobre o uso do Insight Maker, que pode ser consultado em [21[21] M. Lauschner, L.A. Heidemann e E.A. Veit, Tutorial para o uso do software Insight Maker, UFRGS, RS (2020), disponível em: https://cref-ufrgs.github.io/TutorialIM/, acessado em 13/01/2023.
https://cref-ufrgs.github.io/TutorialIM/...
]. Nos limitaremos a apresentar os resultados da simulação com a ferramenta. Os detalhes da construção de cada diagrama usado no artigo se encontram no Material Suplementar. Na Figura 1, explicitamos o resultado da simulação, que é uma integração numérica da Eq. (1) considerando ϕ(0)=30, ϕ.(0)=0rad/s, L=1m, g=9,81m/s2 e m=1kg.

Figura 1:
Resultados da simulação do modelo do pêndulo simples no software Insight Maker. Em (a), é exposta a evolução temporal da posição angular. Em (b), o diagrama de fase do sistema, isto é, ϕ.(ϕ).

2.2. Previsibilidade no modelo de pêndulo simples: análise da sensibilidade às condições iniciais

Uma forma de auferir a previsibilidade de um modelo é avaliar a sensibilidade dos resultados da simulação frente a pequenas variações nas condições iniciais. Em sistemas reais, tais variações são inerentes ao processo de medição e não podemos eliminá-las [23[23] L. Pigosso, Um estudo exploratório sobre atividades investigativas com enfoque no processo de medição no ensino fundamental. Dissertação de Mestrado, Universidade Federal do Rio Grande do Sul, Porto Alegre (2022).]. No caso de um pêndulo, sempre teremos incerteza na medição do seu comprimento L, da sua massa m, da aceleração da gravidade g e do seu ângulo inicial ϕ(0). Como essas incertezas (desconhecimentos das condições iniciais do evento) afetam a previsão do modelo?

Para avaliar a sensibilidade dos resultados da simulação frente as pequenas variações nas condições iniciais, vamos repetir N=500 vezes a simulação no Insight Maker3 3 Fazemos isso com a ferramenta Sensitivity Testing, explicada no Material Suplementar. , de forma que, a cada nova simulação, os valores das variáveis do sistema têm condições iniciais diferentes, seguindo uma distribuição normal4 4 No software, também é possível assumir outras distribuições como a Binomial e de Poisson. . Nesse processo, assumimos as incertezas experimentais típicas como um desvio padrão [13[13] J.P. Lima, M.T.X. Silva, F.L. Silveira e E.A. Veit, Laboratório de mecânica: subsídios para o ensino de Física experimental (Editora da UFRGS, Porto Alegre, 2013).].

Na simulação que segue, utilizamos ϕ(0)=30, σϕ(0)=0,5 e ϕ.(0)=0rad/s. Consideramos os seguintes valores médios e desvios padrão para os outros parâmetros do modelo de pêndulo simples: L=1m e σL=0,5mm5 5 Simulando uma medição feita com uma régua cuja menor escala é 1 mm. ; m=1kg e σm=1g6 6 Simulando uma medição feita com uma balança de precisão digital cujo menor valor não nulo indicado pelo mostrador é 1 g. Com isso temos uma medida conservadora da incerteza [13, p. 36]. ; g=9,759m/s2 e σg=0,003m/s2.7 7 Valor medido em Porto Alegre, segundo Silveira [24]. O resultado do teste de sensibilidade pode ser observado na Figura2.

Figura 2:
Resultado do teste de sensibilidade a variações nas condições iniciais para o pêndulo simples.

A Figura 2 se parece muito com a Figura 1(a), indicando que, mesmo com os quatro parâmetros variando a cada início de simulação, o resultado da predição não se altera significativamente. No destaque apresentado na figura, pode-se ver que, na cor amarela, estão representados aproximadamente 68% dos resultados da simulação; já na cor azul, 95% dos resultados8 8 Na distribuição normal, cerca de 68% dos valores distam até 1 desvio padrão em relação à média, enquanto cerca de 95% distam até dois desvios. .

Concluímos essa seção destacando que, por mais que os resultados obtidos difiram entre si por pequenas variações nas condições iniciais, eles reforçam a noção de que o evento descrito pelo modelo é previsível, já que, pelos resultados da simulação, a diferença no ângulo foi relativamente baixa. Tomando 95% dos dados resultantes das 500 simulações realizadas com variações em diversas das condições iniciais, identificamos incertezas máximas propagadas de 3º nas predições, um valor pequeno dependendo dos objetivos que estão sendo considerados. Assim, apesar das limitações na medida, podemos concluir que é possível prever a evolução de um pêndulo de forma aproximada usando o modelo de pêndulo simples. Na próxima seção faremos uma análise semelhante, porém com o modelo do pêndulo amortecido forçado, mostrando que tal previsibilidade é bastante distinta neste caso.

3. Desafiando a Previsibilidade da Física: O Pêndulo Amortecido Forçado

O pêndulo amortecido e forçado se diferencia do pêndulo simples pela inclusão de duas forças: i. uma relacionada com o arrasto devido ao ar, que assumimos ser sempre contrária à direção do movimento e diretamente proporcional à velocidade (Far=-bv=-bLϕ.v^, onde b é uma constante de proporcionalidade que depende das características do corpo suspenso e do fluido em que está imerso), e ii. uma força periódica externa F(t), na direção do movimento. Admitimos que essa força é senoidal, com frequência angular ω e amplitude de força motriz F0, isto é, F(t)=F0cosωtv^9 9 Essa forma para a força foi escolhida porque gera o efeito do caos que queremos estudar no artigo. Fisicamente, ela poderia representar, por exemplo, o “embalo” fornecido pela mãe que empurra seu filho no balanço ou a tração que controla o movimento de um “barco viking” em um parque de diversões. . O esquema representando as forças é sintetizado na Figura 3.

Figura 3:
Esquema representando o pêndulo amortecido forçado e as forças nele atuantes. Observe que a direção da força externa pode agir tanto contra como a favor do movimento, já que sua frequência angular ω pode não ser ressonante com a oscilação do pêndulo.

Pela Figura 3, constata-se que o torque em relação ao eixo de rotação é dado por τ=[LF (t)-LFar-LP sinϕ]j^. Considerando que Iϕ¨=τ, onde I=mL2 é o momento de inércia do sistema em relação ao eixo de sustentação do pêndulo, e que o torque resultante tem componentes apenas ao longo do eixo j^, a equação de movimento pode ser escrita, na forma escalar, como:

(2) ϕ ¨ = 1 m L [ F ( t ) - F a r - P sin ϕ ] .

Também vamos introduzir um parâmetro adimensional γ, definido como γ=F0/mg, que representa o quão intensa é a força externa em comparação com o peso do corpo. Ou seja, quando γ1, podemos dizer que a força externa é desprezível em relação ao peso do corpo; quando γ1, essa força é comparável com o peso do corpo. Esse é o parâmetro que vamos mudar nas simulações da rota para o caos descrita no artigo. A seguir, exploramos o modelo do pêndulo amortecido forçado resolvendo a Eq. (2) no Insight Maker.

3.1. Implementando o pêndulo amortecido forçado no Insight Maker e a constatação da imprevisibilidade

A construção do diagrama de blocos no Insight Maker para o pêndulo amortecido forçado se dá pela inclusão de mais alguns blocos no diagrama usado para o pêndulo simples. Esses blocos são as novas forças atuantes, bem como os parâmetros a elas associados10 10 A simulação está disponível pelo link em [25]. . Para realizar o teste de sensibilidade, vamos considerar apenas variações no parâmetro ϕ, que dessa vez vai iniciar na vertical: ϕ(0)=0 e σϕ(0)=0,5. Como estamos com intenções didáticas, ou seja, não é nossa prioridade termos precisão nos resultados, focando sobretudo em aspectos conceituais, optamos por rodar a simulação com menos exigências computacionais, limitando-nos à análise das implicações das incertezas apenas no ângulo inicial, e não nas outras grandezas do modelo. Esse procedimento será suficiente para as discussões conduzidas no artigo.

Novamente, executamos o teste de sensibilidade realizando N=500 simulações. Na primeira análise, vamos considerar γ=0,2, isto é, a força motriz é pequena comparada com o peso do corpo. Os outros parâmetros escolhidos são: ϕ.(0)=0rad/s, ω=2π, m=1kg, g=9,759m/s2, L=g/(9π2) e b=3π/2. Na Figura 4 estão exibidos os resultados da simulação.

Figura 4:
(a) Resultado da simulação com o pêndulo amortecido forçado para γ=0,2 e (b) o correspondente teste de sensibilidade.

Na Figura 4 percebemos que o sistema inicia uma fase de transiente, que leva cerca de dois segundos para se ajustar em uma trajetória periódica, com período de T=1s, medido diretamente a partir do gráfico. A análise de sensibilidade mostra que, mesmo variando o ângulo inicial levemente em torno da posição ϕ=0, o que pode ser entendido como uma incerteza no valor da posição inicial do pêndulo, não há variação significativa na predição do modelo. De fato, a predição parece ser melhor até que a do pêndulo simples.

Na análise seguinte, repetimos a simulação, porém agora com γ=1,5. Esse valor implica que o valor máximo da força motriz é 50% maior do que o peso do corpo preso à extremidade da haste. Todos os outros parâmetros foram mantidos como antes. O resultado da simulação está expresso na Figura5.

Figura 5:
(a) Resultados da simulação com o pêndulo amortecido forçado para γ=1,5 e (b) o correspondente teste de sensibilidade.

Começamos por observar um ponto fundamental que distingue as Figuras 4 e 5: enquanto a primeira exibe comportamento periódico, a segunda não o faz. Isso pode ser verificado aumentando o tempo de simulação para os parâmetros selecionados: não importa quanto tempo se observe, não é possível identificar um padrão. Observamos também que o modelo forneceu previsões idênticas até aproximadamente sete segundos mesmo com uma pequena variabilidade na condição inicial do ângulo. A partir daí, no entanto, os resultados divergiram tanto quanto mais afastados da média estavam inicialmente.

Aparentemente, o parâmetro γ controla em algum nível a previsibilidade do sistema. De fato, as Figuras 4 e 5 evidenciam um ponto de transição da ordem para o caos que passamos a estudar em maiores detalhes agora.

3.2. A rota para o caos: duplicação do período

Como visto, a condição necessária para a existência de soluções caóticas para uma equação diferencial são duas [17[17] M. Gitterman, The Chaotic Pendulum (World Scientific Publishing, New Jersey, 2010)., p. 20]: a equação precisa conter termos não-lineares e envolver ao menos três variáveis dinâmicas (três dimensões). Na Eq. (2), o termo não-linear é sinϕ e as variáveis dinâmicas são ϕ, ϕ. e t11 11 A dependência explícita do tempo implica um grau de liberdade extra no sistema, já que o tempo pode ser entendido como uma variável dinâmica que obedece a relação t.=1 (veja, p. ex., (26, p. 5)). . Conforme mencionado, podemos identificar a rota para o caos pela escolha conveniente dos parâmetros. Cabe destacar, no entanto, que não há apenas uma rota para o caos no pêndulo amortecido e forçado. Algumas rotas encontradas até então são: a cascata de duplicação de período, tanto para soluções restritas e não-restritas12 12 Soluções restritas são aquelas em que o ângulo ϕ fica confinado em um intervalo. Quando esse não é o caso, chamamos as soluções de não-restritas. ; as transições aleatórias entre dois estados restritos; a transição intermitente entre estados caóticos não-restritos com rotação alternada em sentido horário e anti-horário [27[27] D. D’Humueres, M.R. Beasley, B.A. Huberman e A. Libchaber, Phys. Rev. A 26, 3483 (1982).]. Neste artigo, vamos explorar a rota para o caos através da cascata de duplicação do período.

Em 1975, o físico Mitchell J. Feigenbaum descobriu de forma empírica uma relação matemática genérica que relaciona a duplicação do período T de sistemas caóticos com um parâmetro ζ que o descreve. Se {T,2T,22T,,2n-1T}n=1 é o conjunto dos valores de período associados aos parâmetros {ζn}n=1, então se obtém o chamado número de Feigenbaum δ pela seguinte expressão:

(3) δ = lim n ζ n - ζ n - 1 ζ n + 1 - ζ n 4 , 6692016 .

Até o presente momento não sabemos explicar sua causa [28[28] L. A. Debnath, J. Appl. Comput. Math. 1, 1 (2015).]. Sabemos, no entanto, que a cascata de duplicação do período ocorre respeitando essa lei para uma grande variedade de sistemas não-lineares, atestando a universalidade desse princípio. Quando trazido à tona em contextos de ensino, tal aspecto evidencia a contemporaneidade da dinâmica não-linear, bem como a importância de estudá-la para buscarmos uma melhor compreensão do mundo. Mais adiante, quando evidenciarmos a duplicação no período no pêndulo amortecido e forçado pelo Insight Maker, vamos obter um valor aproximado para o número de Feigenbaum.

Iniciamos a rota para o caos buscando valores de γ maiores que γ=0,2. Utilizando as mesmas condições iniciais do exemplo anterior, aumentamos esse parâmetro sucessivamente até encontrar um valor limite a partir do qual são produzidas soluções periódicas com período duplicado (T=2s). Encontramos esse resultado para γ=1,0650 com uma precisão de 0,0001, isto é, diminuir o valor de γ em 0,0001 faz com que o período volte a ser de um segundo. Destacamos que esse é um valor próximo ao encontrado na literatura, que é de γ=1,0663[29[29] J.R. Taylor, Classical mechanics (University Science Books, USA, 2005)., p. 474]. O resultado da simulação está exposto na Figura 6.

Figura 6:
Evolução do pêndulo amortecido forçado para γ=1,0650. A simulação foi executada até 40 segundos, buscando simultaneamente eliminar resíduos transientes na evolução temporal, bem como amenizar o erro associado aos arredondamentos realizados na integração numérica.

Como no resultado apresentado na Figura 4, essa simulação apresenta uma fase transiente, que nesse caso dura até logo após 20 segundos. Depois disso, o sistema apresenta regime periódico, porém com período duplicado (T=2s) em relação à simulação com o parâmetro γ=0,2, quando o período era T=1s nas mesmas condições iniciais. Neste ponto, destacamos que o período não varia apenas em função de γ, tampouco apenas de forma a se duplicar. De fato, se a simulação fosse realizada com γ=1,077, seria observada uma evolução periódica com T=3s (período triplicado); ainda, para esse valor de γ, porém com ϕ(0)=-π/2, obteríamos um período de T=1s. Por conta disso, é importante manter as condições iniciais fixas em cada simulação.

Aumentando mais o parâmetro γ, obtemos a próxima duplicação do período em γ=1,0818, com a mesma precisão do exemplo anterior. O resultado é apresentado na Figura 7.

Figura 7:
Evolução do pêndulo amortecido forçado para γ=1,0818.

A Figura 7 evidencia movimento periódico com período T=4s. Continuando esse processo, encontramos que, para γ=1,0825, o período é novamente duplicado (T=8s). Para γ=1,0827 temos mais uma duplicação de período (T=16s). Os resultados obtidos até agora estão condensados na Tabela 1.

Tabela 1:
Evolução do período do movimento em função do parâmetro γ.

Pelo exposto, vê-se que o período duplica para incrementos Δγ cada vez menores do parâmetro γ, indicando que essa sequência é convergente. Os resultados da Tabela 1 podem ser usados para fazer uma primeira aproximação do número de Feigenbaum pela Eq. (3), apesar de que, a rigor, a relação só é valida para {γn};n. Fazendo a conta para os três últimos valores desse parâmetro, obtemos

(4) δ = 1 , 0825 - 1 , 0818 1 , 0827 - 1 , 0825 = 0 , 007 0 , 0002 = 3 , 5 .

Apesar da falta de precisão, essa estimativa concreta ajuda na compreensão do significado da constante de Feigenbaum. Como a sequência de valores de Δγ é convergente, podemos concluir que, com a variação de γ, estamos duplicando o período do pêndulo diversas vezes, até um ponto em que Δγ tende a zero, significando que fazemos infinitas duplicações do período em uma variação Δγ tendendo a zero. Infinitas duplicações do período implicam em um período infinito, ou seja, caos.

A literatura aponta que o limite da sequência da Tabela 1 ocorre para γc=1,0829, onde γc é o valor crítico ou limiar para a ocorrência de caos no sistema do pêndulo amortecido forçado [29[29] J.R. Taylor, Classical mechanics (University Science Books, USA, 2005)., p. 475]13 13 Na referência indicada, os valores para γ1,…,γ4 são levemente diferentes do valor que chegamos, possivelmente em função do método de integração numérica. . Nesse caso, o período tem um valor infinito T=, de forma que o movimento não se repete.

3.3. A sensibilidade às condições iniciais em um coeficiente: o expoente de Lyapunov

Vimos até aqui que sistemas caóticos possuem hipersensibilidade a variações nas condições iniciais, mesmo que essas variações sejam muito pequenas. Se compararmos a trajetória de dois sistemas idênticos que apenas diferem muito pouco nas suas condições iniciais, podemos quantificar sua divergência pelo expoente de Lyapunov, designado pelo símboloλ. Para compreendermos a lógica subjacente a esse coeficiente, começamos expressando sua definição em sistemas discretos. Um sistema dinâmico discreto é representado por uma função matemática que mapeia seu estado em certo instante para o seu estado no instante seguinte, em um conjunto discreto de valores possíveis. Um exemplo famoso é o mapa logístico, muito utilizado na modelagem de sistemas que descrevem a dinâmica populacional. Nesses casos, vamos supor que a diferença d entre as trajetórias dos sistemas comparados, a cada iteração n, seja proporcional à diferença na iteração imediatamente anterior, ou seja, |dn+1|=a|dn|, onde a é uma constante. Desse modo, podemos escrever, por exemplo, |d2|=a|d1|, |d1|=a|d0|, |d2|=a2|d0|, ou seja, |dn|=an|d0|. Usando propriedades de logaritmos, podemos escrever |dn|=eλn|d0|, onde o expoente de Lyapunov é dado por λ=lna[30[30] L.H. Monteiro, Sistemas Dinâmicos (Livraria da Física, São Paulo, 2002)., p. 358]. Um valor de λ>0 indica que o sistema é sensível a condições iniciais, pois as pequenas diferenças entre dois sistemas crescem exponencialmente, podendo haver a ocorrência de caos; caso λ0, o sistema pode ser entendido como estacionário ou periodicamente estável.

A determinação do expoente de Lyapunov para sistemas contínuos descritos por equações diferenciais é mais complexa, especialmente quando a solução analítica não é conhecida, como no nosso caso. Nessa situação, o expoente de Lyapunov mede o comportamento assintótico médio da taxa de expansão/contração local dos eixos de uma hiperesfera infinitesimal de estados centrada nas condições iniciais do sistema [31[31] N.F. Ferrara e C.P.C. Prado, Caos – Uma introdução (Edgar Blücher, São Paulo, 2000)., p. 245]. Está além do escopo deste artigo propor métodos de integração compatíveis com essa definição. Entretanto, como temos fins didáticos, podemos fazer uma aproximação com o caso discreto mencionado no parágrafo anterior. Tal aproximação nos permitirá determinar o sinal do expoente de Lyapunov, não seu valor propriamente dito, porém essa informação é suficiente para evidenciar a possibilidade de caos no sistema. Em analogia com o caso discreto, nossa definição tem a seguinte forma (já empregando as variáveis do pêndulo amortecido forçado):

(5) Δ ϕ ( t ) = Δ ϕ ( 0 ) e λ t ,

onde Δϕ(t)=ϕ1(t)-ϕ2(t) é a diferença na posição angular de dois sistemas que diferem inicialmente pela quantidade pequena Δϕ(0). Linearizando a Eq. (5), obtemos:

(6) x ( t ) = λ t + ln ( Δ ϕ ( 0 ) ) ,

Onde x(t)=ln(ϕ1(t)-ϕ2(t)). Observe que, no formato da Eq. (6), λ é o coeficiente angular e ln(Δϕ(0)), o coeficiente linear de uma reta. Nesse gráfico, uma inclinação positiva/negativa indicará que o expoente de Lyapunov é positivo/negativo. Como os sistemas estudados são altamente sensíveis às condições iniciais, esse método só tem sentido para pequenos intervalos de tempo, já que o próprio erro numérico decorrente dos arredondamentos realizados pelo computador pode ser entendido como uma pequena variação nas condições iniciais. Desse modo, o método que usaremos para avaliar λ não é adequado para pesquisas que demandam precisão. Existem métodos avançados de integração numérica usados por cientistas para resolver esse problema [31[31] N.F. Ferrara e C.P.C. Prado, Caos – Uma introdução (Edgar Blücher, São Paulo, 2000)., 32[32] A.M. Calvão, Estudos de sistemas dinâmicos não lineares: pêndulo duplo, batimentos cardíacos e coletivos de animais. Tese de Doutorado, Universidade Federal Fluminense, Niterói (2014).]. Nosso objetivo aqui é apenas proporcionar uma compreensão conceitual dos efeitos desse expoente em situações de ensino de Física.

Podemos usar o Insight Maker para simular as trajetórias ϕ1(t) e ϕ2(t) de dois sistemas que diferem apenas pelas condições iniciais duplicando o diagrama de blocos usados para estudar o pêndulo amortecido forçado14 14 O diagrama está disponível em [33]. . Nossa primeira simulação foi feita usando as condições iniciais dos exemplos anteriores, porém com o parâmetro γ na região de caos, isto é: γ=1,2>γc. Além disso, usamos uma pequena variação de Δϕ(0)=0,25. Os resultados são mostrados na Figura 8.

Figura 8:
(a) Comparação das trajetórias ϕ1(t) e ϕ2(t) com γ=1,2; (b)x(t) com inclinação positiva, indicando que o expoente de Lyapunov tem esse mesmo sinal.

A Figura 8(a) indica que uma pequena divergência na condição inicial dos sistemas provocou diferenças significativas na previsão dos modelos em pouco tempo de simulação. Esse resultado era esperado porque o parâmetro γ é maior que o limite para haver caos. Se não soubéssemos esse fato de antemão, precisaríamos esperar cerca de 10 segundos de simulação para perceber que as trajetórias dos sistemas divergem. Entretanto, o gráfico na Figura 8(b) nos mostra desde o início da simulação que esse resultado aconteceria. Isso fica claro porque x(t) tem inclinação positiva até esse instante, indicando que o coeficiente de Lyapunov é positivo.

Vimos que a partir de γc=1,0829 o pêndulo amortecido passa a apresentar caos. Entretanto, isso não ocorre para todos os valores γ>γc. De fato, sabe-se que sistemas caóticos alternam intervalos de caos com outros de movimento periódico. Podemos evidenciar isso com o Insight Maker testando para os valores γ=1,5 e γ=1,6 e usando o restantes dos parâmetros como na simulação anterior. O resultado está na Figura 9.

Figura 9:
Simulações com γ=1,5 (a) e γ=1,6(b).

A Figura 9(a) evidencia que, aumentando o parâmetro até γ=1,5, o sistema permanece apresentando caos, com expoente de Lyapunov positivo, como na Figura 8. Entretanto, aumentando o parâmetro mais um pouco, obtemos soluções periódicas (T=3s) estacionárias, de forma que o parâmetro γ=1,6 e seu entorno se constitui em uma “ilha de ordem”, como evidenciado pela Figura 9(b). Prosseguindo aumentando o parâmetro γ, seguem-se mais regiões alternadas de caos com regiões periódicas, que se tornam cada vez mais raras à medida que o parâmetro aumenta.

4. Considerações Finais

Neste artigo, exploramos o software Insight Maker como uma alternativa para estudar sistemas não-lineares com fins didáticos. Como exemplo de aplicação, estudamos o modelo do pêndulo amortecido forçado e exploramos a rota para o caos via cascata de duplicação do período. Isso foi feito iniciando por considerações sobre a noção de causalidade e determinismo, que propagam a concepção equivocada de que eventos naturais são sempre previsíveis à luz das leis da Física. As implicações dessas ideias foram gradualmente testadas na análise do movimento pendular, iniciando pelo teste de sensibilidade às condições iniciais do seu movimento quando descrito pelo modelo de pêndulo simples. Os resultados da simulação atestaram a previsibilidade do evento quando descrito por esse modelo. Entretanto, para um sistema apenas um pouco mais complexo, como o pêndulo sujeito ao arrasto do ar e uma força externa periódica (descrito pelo modelo do pêndulo amortecido forçado), verificamos que o movimento passa a não ser previsível sob certas condições, apesar de ser determinístico e obedecer às leis de Newton. Observamos que o sistema evolui para o caos duplicando seu período, de acordo com uma variação de seus parâmetros que obedecem de forma aproximada a relação empírica de Feigenbaum. Por meio de uma analogia com o caso discreto, avaliamos o sinal do expoente de Lyapunov em situações de movimento caótico e periódico, obtendo resultados coerentes.

Em termos do seu potencial no ensino de Física, a atividade se mostrou frutífera para fomentar discussões em nível qualitativo e quantitativo sobre as implicações da quebra da previsibilidade, mesmo em fenômenos simples como o movimento de um pêndulo. Em uma única atividade, fornecemos subsídios para que os estudantes aprendam ciências (conceitos de mecânica e de dinâmica não-linear) e sobre ciências (elementos de filosofia da Ciência, como causalidade e determinismo) [34[34] D. Hodson, International Journal of Science Education 36, 15 (2014).].

Por fim, os resultados atestaram a viabilidade e potencialidades do estudo de sistemas não-lineares pelo Insight Maker no contexto de ensino. Em particular, destacamos a facilidade de implementação dos códigos via construção de diagrama de blocos, que requerem poucos conhecimentos de programação, bem como a geração automatizada de resultados pelas séries temporais.

Material suplementar

O seguinte material suplementar está disponível online: Utilizando o Software Insight Maker na Prática.

Referências

  • [1]
    J.E. Skinner, M. Molnar, T. Vybiral e M. Mitra, Integr. Physiol. Behav. Sci. 27, 1 (1992).
  • [2]
    S. Vlad, P. Pascu e N. Morariu, Journal of computing 2, 1 (2010).
  • [3]
    J. Gleick, Caos: a criação de uma nova ciência (Campus, Rio de Janeiro, 1989).
  • [4]
    P.C. Ferrari, J.A.P. Angotti e M.H.R. Tragtenberg, Ciência & Educação 15, 1 (2009).
  • [5]
    B.S. Hall, Annals of Science 35, 1 (1978).
  • [6]
    M.R. Matthews, in: International Handbook of Research in History, Philosophy and Science Teaching, editado por M.R. Matthews (Springer, Dordrecht, 2014).
  • [7]
    D.B. Meli, Thinking with Objects (The Johns Hopkins University Press, Baltimore, 2006).
  • [8]
    Some Mathematical Works of the 17th & 18th Centuries, including Newton’s Principia, Euler’s Mechanica, Introductio in Analysin, etc., translated mainly from Latin into English, disponível em: http://www.17centurymaths.com/
    » http://www.17centurymaths.com/
  • [9]
    D.K. Randall, Física – Uma Abordagem Estratégica (Bookman, Porto Alegre, 2009), v. 2, 2 ed.
  • [10]
    M. Paty, Scientiae Studia 2, 1 (2004).
  • [11]
    M. Paty, Scientiae Studia 2, 4 (2004).
  • [12]
    J.R. Alembert, Alembert, J.R. d’ Alembert, em: Encyclopédie ou Dictionnaire raisonné des sciences, des arts et des métiers, editado por D. Diderot (Primeur – Libraire, Chez Pellet, 1752).
  • [13]
    J.P. Lima, M.T.X. Silva, F.L. Silveira e E.A. Veit, Laboratório de mecânica: subsídios para o ensino de Física experimental (Editora da UFRGS, Porto Alegre, 2013).
  • [14]
    P.S. Laplace, Théorie analytique des probabilités (Imprimeur-Libraire, Paris, 1812).
  • [15]
    M. Bunge, Teoria e realidade (Perspectiva, São Paulo, 1974).
  • [16]
    L.A. Heidemann, I.S. Araujo e E.A. Veit, Investigações em Ensino de Ciências 23, 2, 2018.
  • [17]
    M. Gitterman, The Chaotic Pendulum (World Scientific Publishing, New Jersey, 2010).
  • [18]
    R. Duit e M. Komorek, Research in Science Education 27, 3 (1997).
  • [19]
    S. Iezekiel, S. Bennett e C. M. Snowden, Workshop on High Performance Electron Devices for Microwave and Optoelectronic Applications 1, 1 (1997).
  • [20]
    S. Fortmann, Simulation Modelling Practice and Theory 47, 1 (2014).
  • [21]
    M. Lauschner, L.A. Heidemann e E.A. Veit, Tutorial para o uso do software Insight Maker, UFRGS, RS (2020), disponível em: https://cref-ufrgs.github.io/TutorialIM/, acessado em 13/01/2023.
    » https://cref-ufrgs.github.io/TutorialIM/
  • [22]
    Simulação no Software Insight Maker, disponível em: https://insightmaker.com/insight/6QMTJprV9hOxImZ8wqpbzy/P-ndulo-simples-Artigo, acessado em 13/01/2023.
    » https://insightmaker.com/insight/6QMTJprV9hOxImZ8wqpbzy/P-ndulo-simples-Artigo
  • [23]
    L. Pigosso, Um estudo exploratório sobre atividades investigativas com enfoque no processo de medição no ensino fundamental Dissertação de Mestrado, Universidade Federal do Rio Grande do Sul, Porto Alegre (2022).
  • [24]
    F.L. Silveira, Revista de Ensenãnza de la Física 10, 2 (1995).
  • [25]
    Simulação no Software Insight Maker, disponível em: https://insightmaker.com/insight/3NhBG2TW0SAOmkbSEVTCmG/P-ndulo-amortecido-e-for-ado-Artigo, acessado em 13/01/2023.
    » https://insightmaker.com/insight/3NhBG2TW0SAOmkbSEVTCmG/P-ndulo-amortecido-e-for-ado-Artigo
  • [26]
    R.M. Angelo, E.I. Duzzioni e A.D. Ribeiro, J. Phys. A: Math. Theor 45, 055101 (2012).
  • [27]
    D. D’Humueres, M.R. Beasley, B.A. Huberman e A. Libchaber, Phys. Rev. A 26, 3483 (1982).
  • [28]
    L. A. Debnath, J. Appl. Comput. Math. 1, 1 (2015).
  • [29]
    J.R. Taylor, Classical mechanics (University Science Books, USA, 2005).
  • [30]
    L.H. Monteiro, Sistemas Dinâmicos (Livraria da Física, São Paulo, 2002).
  • [31]
    N.F. Ferrara e C.P.C. Prado, Caos – Uma introdução (Edgar Blücher, São Paulo, 2000).
  • [32]
    A.M. Calvão, Estudos de sistemas dinâmicos não lineares: pêndulo duplo, batimentos cardíacos e coletivos de animais Tese de Doutorado, Universidade Federal Fluminense, Niterói (2014).
  • [33]
    Simulação no Software Insight Maker, disponível em: https://insightmaker.com/insight/7Dpw1ajC2gkLZNEhn0r9zV/P-ndulo-amortecido-e-for-ado-Lyapunov-Artigo, acessado em 13/01/2023.
    » https://insightmaker.com/insight/7Dpw1ajC2gkLZNEhn0r9zV/P-ndulo-amortecido-e-for-ado-Lyapunov-Artigo
  • [34]
    D. Hodson, International Journal of Science Education 36, 15 (2014).
  • 1
    Disponível em: https://insightmaker.com/. Acesso em: 29/11/2022.
  • 2
    Por termos fins didáticos, ao invés de suprimir o peso P, que desaparece da equação após a simplificação com a massa m no denominador da Eq. (1), optamos por considerá-lo de forma explícita, o que também será feito com as forças no modelo do pêndulo amortecido forçado. Tal decisão é tomada em função das intenções didáticas do artigo, para que, no momento em que formos construir os modelos no software Insigth Maker, as forças sejam apresentadas explicitamente, favorecendo a análise dos eventos do ponto de vista físico (avaliação das forças atuantes no sistema), e não puramente matemático (exclusivamente como taxas de variação de grandezas).
  • 3
    Fazemos isso com a ferramenta Sensitivity Testing, explicada no Material Suplementar.
  • 4
    No software, também é possível assumir outras distribuições como a Binomial e de Poisson.
  • 5
    Simulando uma medição feita com uma régua cuja menor escala é 1 mm.
  • 6
    Simulando uma medição feita com uma balança de precisão digital cujo menor valor não nulo indicado pelo mostrador é 1 g. Com isso temos uma medida conservadora da incerteza [13[13] J.P. Lima, M.T.X. Silva, F.L. Silveira e E.A. Veit, Laboratório de mecânica: subsídios para o ensino de Física experimental (Editora da UFRGS, Porto Alegre, 2013)., p. 36].
  • 7
    Valor medido em Porto Alegre, segundo Silveira [24[24] F.L. Silveira, Revista de Ensenãnza de la Física 10, 2 (1995).].
  • 8
    Na distribuição normal, cerca de 68% dos valores distam até 1 desvio padrão em relação à média, enquanto cerca de 95% distam até dois desvios.
  • 9
    Essa forma para a força foi escolhida porque gera o efeito do caos que queremos estudar no artigo. Fisicamente, ela poderia representar, por exemplo, o “embalo” fornecido pela mãe que empurra seu filho no balanço ou a tração que controla o movimento de um “barco viking” em um parque de diversões.
  • 10
    A simulação está disponível pelo link em [25[25] Simulação no Software Insight Maker, disponível em: https://insightmaker.com/insight/3NhBG2TW0SAOmkbSEVTCmG/P-ndulo-amortecido-e-for-ado-Artigo, acessado em 13/01/2023.
    https://insightmaker.com/insight/3NhBG2T...
    ].
  • 11
    A dependência explícita do tempo implica um grau de liberdade extra no sistema, já que o tempo pode ser entendido como uma variável dinâmica que obedece a relação t.=1 (veja, p. ex., (26, p. 5)[26] R.M. Angelo, E.I. Duzzioni e A.D. Ribeiro, J. Phys. A: Math. Theor 45, 055101 (2012).).
  • 12
    Soluções restritas são aquelas em que o ângulo ϕ fica confinado em um intervalo. Quando esse não é o caso, chamamos as soluções de não-restritas.
  • 13
    Na referência indicada, os valores para γ1,,γ4 são levemente diferentes do valor que chegamos, possivelmente em função do método de integração numérica.
  • 14
    O diagrama está disponível em [33[33] Simulação no Software Insight Maker, disponível em: https://insightmaker.com/insight/7Dpw1ajC2gkLZNEhn0r9zV/P-ndulo-amortecido-e-for-ado-Lyapunov-Artigo, acessado em 13/01/2023.
    https://insightmaker.com/insight/7Dpw1aj...
    ].

Datas de Publicação

  • Publicação nesta coleção
    17 Abr 2023
  • Data do Fascículo
    2023

Histórico

  • Recebido
    03 Fev 2023
  • Revisado
    28 Mar 2023
  • Aceito
    07 Mar 2023
Sociedade Brasileira de Física Caixa Postal 66328, 05389-970 São Paulo SP - Brazil - São Paulo - SP - Brazil
E-mail: marcio@sbfisica.org.br