Resumos
Apresentamos neste trabalho uma breve revisão sobre o cálculo fracionário e também o método de Adam-Bashforth fracionário para solução de equações diferenciais fracionárias. O método de Adam-Bashforth fracionário é devido a Diethelm et al. [1[1] K. Diethelm, N.J. Ford e A.D. Freed, Nonlinear Dynamics 29, 3 (2002).] e caracteriza-se pela sua eficiência, já que considera os efeitos de memória das derivadas fracionárias. A apresentação do tema é realizada de forma pedagógica, proporcionando a compreensão daqueles que são incipientes no assunto. Nesse contexto, exemplos são apresentados e algumas aplicações são discutidas.
Palavras-chaves:
Cálculo Fracionário; Derivada de Caputo; Método de Adam-Bashforth Fracionário
In this paper we present a brief review about fractional calculus and also the fractional Adam-Bashforth method for solving fractional differential equations. The fractional Adam-Bashforth method is due to Diethelm et al. [1[1] K. Diethelm, N.J. Ford e A.D. Freed, Nonlinear Dynamics 29, 3 (2002).] and is characterized by its efficiency, as it considers the memory effects of fractional differentiation. The presentation is performed in a pedagogical perspective, delivered to beginner in the subject. In this context, some examples are presented and applications are discussed.
Keywords
Fractional Calculus; Caputo’s Derivative; Fractional Adam-Bashforth Method
1. Introdução
Embora pouco popular, o Cálculo Fracionário é mais antigo do que imaginamos. A sua origem é contemporânea do cálculo diferencial de ordem inteira, tradicionalmente estudado nos cursos universitários. O indício do surgimento dessa generalização do cálculo de ordem inteira é datado no século dezessete numa carta que L’Hôspital enviou a Leibniz indagando sobre a possibilidade de uma derivada de ordem [2[2] B. Ross, Historia Mathematica 4, 75 (1977).]. A resposta dada por Leibniz culminou na introdução das primeiras definições de derivada e integral de ordem não-inteira. Em 1812, Laplace estabeleceu a primeira definição de derivada fracionária. A partir desse fato, eminentes matemáticos passaram a investigar as consequências da definição de derivadas fracionárias e aplicá-las na abordagem de problemas físicos. Nesse caminho, a primeira vez que a denominação “derivada fracionária” apareceu na literatura foi em 1819 nos trabalhos de Lacroix. No bojo das aplicações, um exemplo foi o matemático Abel, que em 1823 usou a derivada fracionária para estudar o conhecido problema da tautócrona. O estudo de Abel é considerado como a primeira aplicação prática do cálculo fracionário. Já em 1892, surge a definição de derivada e integral fracionária proposta por Riemann-Liouville. Tal definição é muito estudada até hoje, e tem sido aplicada na modelagem de diversos problemas. Contudo, a derivada fracionária de Riemann-Liouville mostrou-se um pouco obsoleta por alguns motivos, sendo um deles o fato da possibilidade da derivada fracionária de uma constante ser diferente de zero [3[3] R.F. Camargo e E.C. Oliveira, Cálculo fracionário (Livraria da Física, São Paulo, 2015).]. A partir desse problema, em 1967 emerge a derivada fracionária proposta por Michele Caputo [4[4] M. Caputo, Geophysical Journal of the Royal Astronomical Society 13, 529 (1967).], a qual soluciona parte dos desconfortos oriundos da versão de Riemann-Liouville, e principalmente por esse motivo, tem sido aplicada de forma abundante em problemas da física e da engenharia [3[3] R.F. Camargo e E.C. Oliveira, Cálculo fracionário (Livraria da Física, São Paulo, 2015)., 5[5] K.B. Oldhan e J. Spanier, The Fractional Calculus: Theory and Aplications, Differentiation and Integration to Arbitrary Order (Academic, New York, 1974)., 6[6] N. Varalta, A.V. Gomes e R.F. Camargo, TEMA 15, 211 (2014).]. A ideia inicial de Caputo era desenvolver uma abordagem da viscoelasticidade usando cálculo fracionário, tendo muito sucesso em sua proposta. Uma característica do cálculo fracionário é que a derivada fracionária pode ser definida de muitas formas distintas, e geralmente, o resultado de um definição geralmente não coincide com o de outras. As definições apresentadas por Riemann-Liouville e Caputo são as mais conhecidas e serão abordadas neste trabalho, com ênfase dada na formulação de Caputo devido a proposta deste trabalho. O cálculo fracionário tem se desenvolvido bastante a partir de 1967 e tem sido utilizado para modelar diversas situações nas mais variadas áreas, quais sejam: física [7[7] R. Herrmann, Fractional Calculus: An Introduction to Physicists (World Scientific, Singapura (2011).]; ciências mecânicas [8[8] P.J. Torvik e R.L. Bagley, Journal of Applied Mechanics 51, 94 (1984)., 9[9] W. Labecca, O. Guimarães e J.R.C. Piqueira, Mathematical Problems in Engineering, 591715 (2015)., 10[10] J. Vaz Jr e E.C. Oliveira, Mathematics in Engineering 4, 1 (2022).]; modelos matemáticos para ciências biológicas [11[11] A.D. Freed e K. Diethelm, Fractional Calculus and Applied Analysis 10, 219 (2007)., 12[12] T.W.J. Almeida, Sistema para análise viscoelástica de tecidos moles por ondas de cisalhamento usando excitação magnética e medida ultrassônica. Tese de Doutorado, Universidade de São Paulo, Ribeirão Preto (2015).]; dentre muitas outras. As principais características da derivada fracionária que ampliam o seu horizonte de aplicações são a não-localidade e o efeito de memória [3[3] R.F. Camargo e E.C. Oliveira, Cálculo fracionário (Livraria da Física, São Paulo, 2015).]. Nesse sentido, um dos campos de estudo teórico do cálculo fracionário são as equações diferenciais fracionárias. Relevantes avanços já foram alcançados com a demonstração de teoremas e o desenvolvimento de técnicas de solução [5[5] K.B. Oldhan e J. Spanier, The Fractional Calculus: Theory and Aplications, Differentiation and Integration to Arbitrary Order (Academic, New York, 1974).]; porém, há muito ainda por desenvolver sobretudo por haver equações diferenciais fracionárias relativamente simples para as quais não são conhecidas soluções analíticas. É justamente nesse ponto que surge a relevância deste artigo. Apresentaremos de forma pedagógica alguns elementos do cálculo fracionário e como principal constructo, abordaremos alguns métodos numéricos de solução de equações diferenciais fracionárias. Assim sendo, na seção 2 apresentamos os elementos matemáticos fundamentais para a compreensão do trabalho; além da definição de derivada fracionária segundo Caputo. Na seção 3 introduziremos as equações diferenciais fracionárias. Na seção 4 apresentaremos um método numérico para equações fracionárias e análogo ao método de Adams-Bashforth para equações diferenciais de ordem inteira. Três aplicações de áreas distintas serão apresentadas na seção 5. Enquanto na seção 6 elencamos as considerações finais e perspectivas.
2. Preliminares
Para que a derivada fracionária possa ser definida e aplicada será necessário definir e explorar algumas funções especiais que ajudarão a construir o conhecimento necessário para o desenvolvimento dessa. A seguir revisaremos a Funções Gama, a Função Beta e a Função de Mittag-Leffler.
2.1. Função Gama
A Função Gama, ou integral de Euler de segunda espécie, é uma das funções especiais mais importantes na matemática, tendo sido apresentada pela primeira vez pelo matemático Euler em 1730, como resultado de uma pesquisa sobre uma forma de interpolação do fatorial de um número [2[2] B. Ross, Historia Mathematica 4, 75 (1977).]. Foi estudada por diversos outros matemáticos, mas, Adrian Marie Legendre, no ano de 1809, quem a definiu na forma que trabalhamos nos dias atuais. Para , com , a Função Gama é definida por [13[13] J.P.S. Ramírez, Função Gama, disponível em: https://www.ime.unicamp.br/~ftorres/ENSINO/MONOGRAFIAS/JP_VC2_2015.pdf (2015).
https://www.ime.unicamp.br/~ftorres/ENSI...
]
Um dos fatos notáveis da Função Gama é que ela generaliza o fatorial de um número natural n. É imediato verificar a partir de uma integração por partes que
para qualquer z. Consequentemente, por indução, para temos que
2.2. Função Beta
A Função Beta foi introduzida pela primeira vez por Euler e também é conhecida por integral de Euler de primeira espécie, definida como [13[13] J.P.S. Ramírez, Função Gama, disponível em: https://www.ime.unicamp.br/~ftorres/ENSINO/MONOGRAFIAS/JP_VC2_2015.pdf (2015).
https://www.ime.unicamp.br/~ftorres/ENSI...
]
em que . A Função Beta apresenta diversas propriedades interessantes as quais podemos destacar: (i) a comutatividade
cuja demonstração é imediata, e (ii) a relação entre a Função Beta e a Função Gama
Para demonstrarmos essa identidade consideramos o produto
Introduzindo as mudanças de variáveis e , na Eq. (6), pode-se escrever,
Usando as coordenadas polares no plano e , em que o Jacobiano da transformação é r, então e
Aqui vamos notar que
Pois, fazendo a mudança de na Eq. (4), segue que,
em que e . Assim, temos que
Portanto,
Finalmente, fazendo a mudança de variáveis e 2rdr = dy em Eq. (7) temos
Logo, chegamos à identidade desejada
2.3. Função de Mittag-Leffler
A Função de Mittag-Leffler é conhecida como uma generalização da função exponencial. Essa função, conforme introduzida por Leffler, é uma função de variável complexa que depende de um parâmetro complexo , onde , e é definida na forma de uma série de potências [14[14] D.S. Oliveira, Derivada Fracionária e as Funções de Mittag-Leffler. Dissertação de Mestrado, Universidade Estadual de Campinas, Campinas (2014).]
Note que no caso em que , tem-se:
o que permite-nos observar que a Função de Mittag-Leffler admite, como caso particular, a função exponencial.
2.3.1. Função de Mittag-Leffler de dois parâmetros
A Função de Mittag-Leffler de dois parâmetros, , foi introduzida por Wiman [14[14] D.S. Oliveira, Derivada Fracionária e as Funções de Mittag-Leffler. Dissertação de Mestrado, Universidade Estadual de Campinas, Campinas (2014).] e é uma função de variável complexa que depende de dois parâmetros complexos, e , onde , dada pela equação a seguir,
A condição é exigida para garantir a convergência da série dada na Eq. (9).
É imediato observar que quando a Eq. (9) se reduz a Função de Mittag-Leffler de um parâmetro, ou seja,
2.3.2. Transformada de Laplace da Função deMittag-Leffler
Nas seções posteriores necessitaremos da transformada de Laplace da Função de Mittag-Leffler. Por essa razão, apresentamos neste momento os passos de sua obtenção. A transformada de Laplace de uma função , definida para , é dada pela integral
desde que a integral convirja. Com isso, podemos calcular a transformada de Laplace da seguinte generalização da Função de Mittag-Leffler com dois parâmetros
em que é um parâmetro qualquer.
Como a Transformada de Laplace de é , obtemos
Logo,
Se recordarmos da série geométrica que é dada por
quando , temos que
contanto que
Finalmente, reescrevendo a Eq. (10), temos que a transformada de Laplace da Função de Mittag-Leffler é dada por,
Este resultado será relevante na sequência do trabalho.
2.4. A derivada fracionária de Caputo
Nesta seção apresentaremos a definição de derivada fracionária segundo Caputo. Dedicaremos a nossa atenção a este tipo de derivada fracionária, em detrimento das outras possibilidades, devido a sua importância na análise de sistemas físicos [3[3] R.F. Camargo e E.C. Oliveira, Cálculo fracionário (Livraria da Física, São Paulo, 2015).].
Definição: A derivada fracionária de Caputo de ordem é definida por [3[3] R.F. Camargo e E.C. Oliveira, Cálculo fracionário (Livraria da Física, São Paulo, 2015).]
em que , e é a n-ésima derivada de ordem inteira. As derivadas de Caputo são válidas para , desde que , no entanto, aqui será considerado apenas . Daqui em diante, por não haver confusão com outra definição de derivada fracionária, escreveremos o operador derivada de Caputo na forma
Para fins didáticos, vamos calcular a derivada fracionária segundo Caputo da função potência , com . Por definição ela é dada por
em que Então,
Realizando a mudança de variáveis u = ty, tem-se que,
Assim,
Temos da Eq. (5) que
obtêm-se a derivada fracionária segundo Caputo,
É imediato notar que a derivada fracionária de uma constante, segundo a definição de Caputo, é sempre igual a zero, o que não ocorre no caso da derivada fracionária segundo Riemann-Liouville [15[15] A.A. Kilbas, H.M. Srivastava e J.H. Trujillo, Theory and Application of Fractional Differential Equations (Elsevier, Reino Unido, 2006)., 16[16] K.M. Owolabi, A. Atangana, Numerical Methods for Fractional Differentiation (Springer, Singapura, 2019).], por exemplo. Esse fato é uma das justificativas por se optar pela proposta de Caputo.
Na sequência, apresentamos a transformada de Laplace da derivada fracionária segundo Caputo. A demonstração dessa propriedade será deixada como exercício ao leitor e pode ser encontrada na referência [15[15] A.A. Kilbas, H.M. Srivastava e J.H. Trujillo, Theory and Application of Fractional Differential Equations (Elsevier, Reino Unido, 2006)., 16[16] K.M. Owolabi, A. Atangana, Numerical Methods for Fractional Differentiation (Springer, Singapura, 2019).]. Assim, temos
em que representa a derivada de ordem inteira k da função no ponto t=0.
3. Equações Diferenciais Fracionárias
Nesta seção apresentaremos exemplos de equações diferenciais fracionárias no escopo de Caputo que podemos solucionar analiticamente. Não abordaremos teoremas que garantem a existência e unicidade da solução, pois fogem ao escopo desta apresentação, contudo o leitor interessado pode encontrar esse assunto nas referências [15[15] A.A. Kilbas, H.M. Srivastava e J.H. Trujillo, Theory and Application of Fractional Differential Equations (Elsevier, Reino Unido, 2006)., 20[20] K. Diethelm, The Analysis of Fractional Differential Equations: An Application-Oriented Exposition Using Differential Operators of Caputo Type (Springer, Singapura, 2010).].
As equações diferenciais fracionárias que serão estudadas neste trabalho possuem a forma geral
em que , é uma função suficientemente regular de t e , e representa a derivada de ordem inteira k da função no instante t=0.
Mostraremos a seguir o método da trasformada de Laplace para tratar tais problemas. A viabilidade da solução analítica através desse tipo de estratégia dependerá da forma da função . Observamos que nas aplicações que serão apresentadas na seção 5 5. Aplicações Nesta seção aplicaremos o método descrito na seção anterior para analisar brevemente três modelos baseados em equações diferenciais fracionárias. O primeiro deles é o oscilador harmônico fracionário submetido a uma força externa, outro é a versão fracionária da equação logística e, por fim, trataremos de um modelo de viscoelasticidade para tecidos moles. Esperamos que o tratamento aqui apresentado contribua de forma significativa nas pesquisas sobre cálculo fracionário e sirva de referência literária introdutória na área. 5.1. O oscilador harmônico fracionário submetido a uma força externa senoidal O oscilador harmônico é um sistema físico constituído por um objeto de massa m preso a uma mola de constante k. A equação diferencial que governa o movimento do oscilador, isto é, a equação que fornece a posição do oscilador em função do tempo é dada por (34) d 2 x ( t ) d t 2 + ω o 2 x ( t ) = 0 , em que ωo2=km. A Eq. (34) é usada para analisar o oscilador livre. Se por acaso o oscilador estiver submetido a alguma força externa de módulo Fext, podemos escrever (35) d 2 x ( t ) d t 2 + ω o 2 x ( t ) = F e x t . É relevante aqui recordarmos que o oscilador harmônico é um modelo que descreve bem muitos sistemas físicos reais, como por exemplo sistemas que oscilam em torno de uma posição de equilíbrio. Ademais, se considerarmos um oscilador com dissipação, isto é, se as forças de atrito forem consideradas, obtemos os casos especiais do oscilador amortecido, superamortecido e subamortecido, cada qual com a sua aplicação prática. Um fato bastante peculiar no caso do oscilador harmônico fracionário é a possibilidade de se descrever dissipação sem a necessidade de introduzir o termo associado à força de atrito na equação do movimento; bastando para isso variar a ordem da derivada fracionária. Além disso, podemos generalizar o oscilador para um tratamento quântico, submetendo a equação de Schrödinger ao potencial do oscilador harmônico, o que nos permite descrever por exemplo o comportamento de moléculas diatômicas com bastante acurácia. Porém, o problema do oscilador harmônico no âmbito da mecânica quântica ainda é um problema muito pouco explorado. O correspondente fracionário do modelo físico dado na Eq. (35) quando submetido a uma força externa senoidal é dado por [3, 37] (36) D t α x ( t ) + ω o α x ( t ) = F o sin ( ω t ) , em que α∈(1,2) e Fo é a amplitude da força externa. Observe que quando 1<α<2 duas condições iniciais são necessárias para a solução do problema, a saber: x(0) e x′(0). Observa-se que uma interpretação física para o oscilador harmônico fracionário foi dada por Stanislavsky [23]. Ele sugeriu que um oscilador harmônico fracionário deve ser interpretado como uma média de osciladores harmônicos comuns (de derivada ordinária) governado por uma seta do tempo estocástica. Além disso, do ponto de vista matemático os osciladores harmônicos fracionários também são úteis para determinar a solução da equação de Schrödinger fracionária livre [7, 24]. Ressalta-se ainda que aspectos fisicamente significativos do sistema dinâmico proporcionado pela equação do oscilador harmônico fracionário, como energia total e resposta dinâmica do sistema, também já foram estudados na literatura [25, 26, 27]. A solução da Eq. (36) para o caso livre, isto é, F0 = 0, pode ser obtida de maneira analítica diretamente a partir da transformada de Laplace, considerando as Eq. (11) e Eq. (14), e é dada por (37) x ( t ) = x ( 0 ) E α ,1 ( - ω 0 α t α ) + x ′ ( 0 ) t E α ,2 ( - ω 0 α t α ) . No entanto, o mesmo não se pode dizer sobre a solução analítica para o caso em que F0≠0. Ainda é possível encontrar solução analítica para casos particulares da Eq. (36) utilizando a transformada de Laplace e o teorema dos resíduos, mas a solução está longe de ser encontrada trivialmente. Por exemplo, para α=3/2, ω0=ω=1, F0 = 1 e x(0)=x′(0)=1 é possível mostrar que a solução analítica da Eq. (36) é (38) x ( t ) = E 3 2 ,1 ( - t 3 / 2 ) + t E 3 2 ,2 ( - t 3 / 2 ) + 4 3 e - t / 2 cos ( 3 2 t ) + sen ( t ) 2 - cos ( t ) 2 2 - 2 . Nesse sentido, trataremos aqui a Eq. (36) utilizando o método numérico desenvolvido na seção anterior. Nosso objetivo aqui é duplo: no primeiro momento estudaremos o que ocorre no modelo fracionário quando ω≠ωo, e em seguida, veremos como o sistema fracionário se comporta no caso em que a frequência natural coincide com a da força externa, isto é, ω=ωo. Primeiramente, vamos avaliar a qualidade da aproximação numérica a comparando à solução analítica dada pela Eq. (38). Na Figura 3 podemos ver tal comparação e qualitativamente é possível observar que a solução obtida com o método numérico descrito se aproxima satisfatoriamente da solução analítica. Figura 3 Solução analítica e numérica da equação Dt3/2=sen(t)-y no intervalo [0,100] com condições iniciais y(0)=1 e y′(0)=1. Para este problema a solução analítica é dada pela Eq. (38). Além disso, para a solução numérica foi utilizado o passo h=10-1. Na Figura 4 podemos ver o caso do oscilador harmônico para ω0=(0,5)2/3 e ω=1 para diferentes valores de α (1,2; 1,4; 1,6 e 1,8). Percebemos que quanto menor é a ordem da derivada fracionária, menor é a amplitude da oscilação. E ainda, percebemos também que neste caso há uma espécie de amortecimento da oscilação, mesmo não tendo considerado ações de forças dissipativas. Figura 4 Solução numérica do oscilador harmônico fracionário Dtα=sen(t)-12y no intervalo [0,100] com condições iniciais y(0)=1 e y′(0)=-1 para diferentes valores de α∈(1,2). Além disso, para a solução numérica foi utilizado o passo h=10-1. Finalmente, na Figura 5 podemos ver o caso do oscilador harmônico para ω0=ω=1 para diferentes valores de α (1,2; 1,4; 1,6 e 1,8). Conforme já conhecido da literatura, quando as frequência da fonte externa é igual à frequência natural de oscilação do sistema, no cálculo de ordem inteira, temos o fenômeno da ressonância, ou seja, a amplitude vai crescendo indefinidamente. Contudo, notamos no gráfico e também no exemplo dado pela equação (38), que enquanto a ordem da derivada fracionária diminui, a oscilação se estabiliza, isto é, além da amplitude do movimento oscilatório ser menor, a ressonância deixa de acontecer. No caso de oscilações forçadas de um sistema mecânico, infere-se deste resultado que mesmo uma força externa com mesma frequência não levará o sistema a se romper. Figura 5 Solução numérica do oscilador harmônico fracionário Dtα=sen(t)-y no intervalo [0,100] com condições iniciais y(0)=1 e y′(0)=-1 para diferentes valores de α∈(1,2). Além disso, para a solução numérica foi utilizado o passo h=10-1. 5.2. Equação logística fracionária A equação logística (39) d y ( t ) d t = k y ( t ) ( 1 - y ( t ) ) , em que k é uma constante, foi proposta originalmente por Verhulst em 1838 [28]. Em seu trabalho, Verhulst tentou aprimorar o modelo de Malthus para a estatística populacional [29]. Desde então, a equação logística passou a ser aplicada em diversos contextos, desde modelos biológicos de proliferação de bactérias até predição de sistemas econômicos [30, 31]. Uma hipótese básica relacionada à Eq. (39) é que, para o processo de variação da população, o tempo é uma quantidade passiva homogênea que indexa os eventos de mudança da quantidade por meio do tempo regularmente medido por um relógio. No entanto, este não é necessariamente o caso para fenômenos complexos e, assim, a derivada temporal é continuada analiticamente para uma derivada fracionária para produzir a equação logística de ordem fracionária [32, 30] (40) D t α y ( t ) = k α y ( t ) ( 1 - y ( t ) ) , em que k é uma constante, 0<α<1 e com uma condição inicial y(0)=y0. Cabe ressaltar que devido aos efeitos de memória e também a não-localidade das derivadas de Caputo, diversos modelos populacionais da literatura, e não apenas a equação logística, recentemente foram estendidos para suas respectivas versões fracionárias com o intuito de avaliar melhor a evolução de doenças [33, 34, 35]. Figura 6 Solução numérica da equação logística fracionária Dtα=y(1-y) no intervalo [0,10] para diferentes valores de α∈(0,1). (a) Solução com condição inicial y0=1/2; (b) solução com condição inicial y0=3/2. Além disso, para a solução numérica foi utilizado o passo (a) h=10-2; e (b) h=10-3. Note que diferentemente do oscilador harmônico a Eq. (40) é não-linear e portanto os métodos envolvendo a transformada de Laplace dada pela Eq. (14) não podem ser aplicados. A resolução desta equação pode aparecer na literatura de algumas formas distintas: o tratamento de uma versão linearizada [31], solução analítica através de expansão em algum tipo de série [32], solução analítica através de alguma hipótese específica sobre a forma da solução [36] e abordagem via métodos numéricos [30]. A nossa abordagem aqui será realizada com o uso do método construído na seção precedente. Na Figura 6(a) podemos ver as soluções da equação logística para diferentes valores de α (0,2; 0,4; 0,6; 0,8 e 1,0) com a condição inicial y0=1/2 e h=10-2. Vemos assim que quanto menor a ordem da derivada fracionária, a solução da equação logística fracionária cresce de forma mais suave e se estabiliza num valor mais baixo, isto é, o valor suporte é menor. Se tomarmos como exemplo o crescimento de uma população, se o problema puder ser modelado via cálculo fracionário, a população que cresce mais devagar e se estabiliza num valor menor é modelada por uma derivada de ordem fracionária mais baixa. Por outro lado, na Figura 6(b) vemos as soluções da equação logística para diferentes valores de α (0,2; 0,4; 0,6; 0,8 e 1,0) com a condição inicial y0=3/2 e h=10-3. Note que neste caso o passo é menor que o exemplo anterior para garantir uma solução suave para o caso α=0,2. Novamente podemos ver que há uma estabilização mais suave da solução quando a ordem da derivada fracionária é menor. 5.3. Modelo de viscoelasticidade de ordem fracionária para tecidos moles Os tecidos biológicos moles podem ser considerados como materiais viscoelásticos, isto é, apresentam propriedades elásticas e viscosas ao mesmo tempo [38]. A lei elástica de Fung [39] é o modelo padrão-ouro na literatura para modelar a resposta elástica de tecidos moles [38, 11, 12]. Ela descreve o comportamento não-linear da tensão-deformação de um tecido biológico mole e foi deduzida a partir de observações experimentais. A lei de Fung é descrita pela equação diferencial ordinária (41) d s ( λ ) d λ = E + β s ( λ ) , s ( 1 ) = 0 , em que s é a tensão Lagrangiana, λ o alongamento, E é uma constante com unidades de tensão e β é uma constante adimensional. Ambas as constantes são características do material envolvido. De forma não-trivial pode-se estender a lei de Fung para o modelo [11] (42) ( 1 + τ d d t ) s ( t ) = E β ( e β ϵ ( t ) - 1 ) + ρ ( E + β s ( t ) ) d d t ϵ ( t ) que descreve as características viscoelásticas dos tecidos moles. Aqui os parâmetros materiais são constantes não-negativas tais que E tem unidades de tensão, β é adimensional, ρ e τ tem unidades de tempo com ρ>τ. Além disso, a função ϵ(t) representa a deformação e é supostamente conhecida. No entanto, este modelo apresentado peca ao modelar com precisão dados experimentais [11]. Já foi mostrado que as propriedades viscoelásticas desses tecidos são melhores tratadas com modelos viscoelásticos de ordem fracionária [40, 41]. Desta forma, o modelo dado pela Eq. (42) é continuado analiticamente para a equação diferencial fracionária [11] (43) ( 1 + τ α D t α ) s ( t ) = E β ( e β ϵ ( t ) - 1 ) + ρ α ( E + β s ( t ) ) D t α ϵ ( t ) , em que 0<α<1. Aqui vamos replicar o experimento do teste de tensão da referência [11]. Com isso, iremos simular o comportamento do colágeno para uma função de deformação ϵ(t) conhecida. Desta forma, a Eq. (43) será resolvida numericamente com os seguintes parâmetros: α=0,33; β=12,4; E=1 MPa, τ=10μs e ρ=10 ms. O teste de tensão será feito utilizando uma deformação linear ϵ(t)=kt, com taxa de deformação k constante, para um material em estado virgem, isto é, s(0)=0. Além disso, serão considerados três valores distintos de k, a saber: k=0,01, k=0,01 e k=0,001. A simulação será feita no intervalo de tempo 0≤t≤tmax, com tmax escolhido tal que ϵ(tmax)=0,15. Finalmente, para a função de deformação escolhida temos, pela Eq. (13), que Dtαϵ(t)=kt1-αΓ(2-α). Com isso, o problema a ser resolvido numericamente é dado pela equação fracionária linear (44) D t α s ( t ) = τ - α ( ρ α β k t 1 - α Γ ( 2 - α ) - 1 ) s ( t ) + τ - α E ( 1 β ( e β k t - 1 ) + ρ α k t 1 - α Γ ( 2 - α ) ) . Na Figura 7 podemos ver a solução para o teste de tensão com gráfico apresentado na forma tensão-deformação. Os resultados apresentados são compatíveis com resultados experimentais e correspondem aos mesmos apresentados por Freed et al. [11]. Novamente, percebe-se a adequação do método numérico ao problema de viscoelasticidade do colágeno e a sua importância para a solução de problemas teóricos e reais. Figura 7 Resultado do teste de tensão para a simulação do colágeno modelado pela Eq. (44). Aqui foram utilizados os seguintes dados: deformação linear ϵ(t)=kt, com k=0,1, k=0,01 e k=0,001. Cada uma das simulações foi realizada no intervalo de tempo [0,tmax], de tal forma que para cada valor de k a igualdade ϵ(tmax)=0,15 deve ser respeitada. apenas o oscilador harmônico permitirá essa análise analítica, pois é a única em que a função f é linear na variável y. No caso da equação logística veremos que a solução analítica através da transformada de Laplace é inviável, considerando que neste caso a que é não linear.
Nesse sentido, como exemplo ilustrativo da utilização da transformada de Laplace na solução de equações fracionárias, vamos resolver o problema de valor inicial fracionário dado por
em que . Para esse fim, tomemos a transformada de Laplace da Eq. (16), isto é
Denominando e utilizando os resultados das seções anteriores obtemos
ou seja
Usando a Eq. (11), podemos calcular a transformada inversa da última equação, obtendo
Ou seja, a solução é dada em termos da função de Mittag-Leffler de um parâmetro. O gráfico apresentado na Figura 1 mostra o comportamento da solução dada na Eq. (17) para diferentes valores de (; ; ; e ). Observamos que à medida que aumentamos o valor da ordem da derivada fracionária, o decrescimento da função y se torna mais rápido, lembrando que quando temos que .
Gráfico da solução do problema de valor inicial Eq. (15) para diferentes valores de e y0=1.
4. Um Método Numérico para Equações Diferenciais Fracionárias
Diferentemente das derivadas comuns de ordem inteira, que são funcionais locais, as derivadas fracionárias, por definição, são funcionais que possuem uma memória total de estados passados. Isso se deve ao fato das derivadas fracionárias serem definidas a partir de uma convolução entre um núcleo e uma função. Essa observação nos leva a crer que os métodos do tipo Runge-Kutta não serão os mais eficientes para tratar as equações fracionárias, exatamente por levarem em consideração apenas o comportamento local da função. Notamos aqui que métodos dessa forma podem ser facilmente estendidos para os casos de equações diferenciais fracionárias através da série de Taylor generalizada e são encontrados na literatura apresentando bons resultados para casos particulares [17[17] Z.M. Odibat e N.T. Shawagfeh, Applied Mathematics and Computation 186, 286 (2007)., 18[18] S. ÖLrekçi, Advances in Mathematical Physics, 507970 (2015)., 19[19] M.S. Arshad, D. Baleanu, M.B. Riaz e M. Abbas, Discrete Dynamics in Nature and Society, 1020472 (2020).]. Apesar disso, nosso interesse é trabalhar com um método numérico que não tenha as limitações locais dos métodos do tipo Runge-Kutta. Portanto, iremos descrever nessa seção um método em duas etapas do tipo preditor-corretor, que pode ser visto como uma variação do método de Adams-Bashforth, para lidar com problemas do tipo dado pela Eq. (15). O método aqui descrito é devido a Diethelm et al. [1[1] K. Diethelm, N.J. Ford e A.D. Freed, Nonlinear Dynamics 29, 3 (2002)., 20[20] K. Diethelm, The Analysis of Fractional Differential Equations: An Application-Oriented Exposition Using Differential Operators of Caputo Type (Springer, Singapura, 2010).] e maiores detalhes analíticos sobre o método, como análise de erro e convergência, podem ser encontrados nas referências [1[1] K. Diethelm, N.J. Ford e A.D. Freed, Nonlinear Dynamics 29, 3 (2002)., 20[20] K. Diethelm, The Analysis of Fractional Differential Equations: An Application-Oriented Exposition Using Differential Operators of Caputo Type (Springer, Singapura, 2010)., 21[21] K. Diethelm, N.J. Ford e A.D. Freed, Numerical Algorithms 36, 31 (2004)., 22[22] K. Diethelm, N.J. Ford, A.D. Freed e Y. Luchko, Computer Methods in Applied Mechanics and Engineering, 194, 743 (2005).]. A descrição do método aqui será feita de tal forma que o leitor possa desenvolver métodos numéricos análogos que possam satisfazer suas necessidades.
Aqui queremos tratar do problema de valor inicial dado pela equação diferencial fracionária
em que e , com a condição inicial
Observamos que é possível descrever o método para e a escolha por é apenas por motivos didáticos.
Começamos a construção do método numérico definindo o operador integral
conhecido como integral fracionária de Riemann-Liouville [16[16] K.M. Owolabi, A. Atangana, Numerical Methods for Fractional Differentiation (Springer, Singapura, 2019)., 20[20] K. Diethelm, The Analysis of Fractional Differential Equations: An Application-Oriented Exposition Using Differential Operators of Caputo Type (Springer, Singapura, 2010).]. Aplicando o operador em temos que
Fazendo uma mudança de ordem de integração é possível reescrever a Eq. (19) na forma
Agora com a mudança de variável podemos verificar que a integral em Eq. (20) é dada por
Novamente realizando uma mudança de variável fazendo temos que a integral em Eq. (21) será reescrita na forma
Notamos, pela identidade dada na Eq. (5), que
e finalmente temos o efeito da aplicação do operador na derivada fracionária :
Aplicando o operador na equação diferencial fracionária dada pela Eq. (18) encontramos a seguinte equação integral que será a base para o método aqui descrito
em que .
Considere uma discretização equiespaçada do domínio com , para todo . Assim, no ponto temos
em que para todo . O que faremos agora é tentar aproximar as integrais que aparecem na Eq. (22) usando a regra dos trapézios, considerando o termo como uma função peso. Logo, sendo a interpolação linear da função no intervalo temos que
Assim, as integrais no lado direito da Eq. (22) podem ser aproximadas por
em que os coeficientes são calculadas através das seguintes integrais
Observamos que as integrais dadas pelas Eqs. (23)–(25) são imediatas e calculando-as obtemos os coeficientes desejados
Finalmente, podemos escrever que a solução aproximada no instante de tempo é dada por
Note que na fórmula encontrada na Eq. (29), para a aproximação da solução da equação diferencial, o valor da função que queremos calcular no instante , , aparece naturalmente também no lado direito da expressão, como argumento da função f. Como esse valor não é conhecido precisamos estimá-lo e ao valor desta estimativa denominaremos . Para tal utilizaremos um procedimento análogo ao já descrito anteriormente modificando apenas a quadratura realizada nas integrais da Eq. (22), ao invés de as aproximarmos utilizando a regra dos trapézios iremos considerar que a função dentro de cada subintervalo é constante e igual a . Isto é,
com para todo , sendo os coeficientes dados por
Finalmente, temos a partir da Eq. (22) que a estimativa para pode ser escrita como
Portanto, a solução numérica a ser obtida pelo método de Adam-Bashforth fracionário é construída a partir da predição de , utilizando a Eq. (31), seguida da correção dada pela Eq. (29) e que fornece a aproximação desejada correspondente à .
A fim de verificarmos a qualidade do método aqui descrito, temos no gráfico apresentado na Figura 2 uma comparação entre a solução analítica e numérica para um problema de valor inicial. Aqui foi determinada a solução do problema de valor inicial dado pela equação , no intervalo , com condição inicial . Conforme mostramos na Eq. (17) a solução analítica deste problema é dada pela função de Mittag-Leffler . No caso da solução numérica observamos que foi utilizado o algoritmo de Adam-Bashforth fracionário com passo . É imediato verificar qualitativamente a adequação da solução numérica em relação a solução analítica.
Finalmente, se faz importante observar que para problemas de valor inicial para na forma
é imediato mostrar que as Eqs. (29)–(31) são ligeiramente modificadas e assumem a forma
Este formato será importante na aplicação do oscilador harmônico fracionário.
Solução analítica e numérica da equação no intervalo com condição inicial . Para este problema a solução analítica é dada pela função de Mittag-Leffler . Além disso, para a solução numérica foi utilizado o passo .
5. Aplicações
Nesta seção aplicaremos o método descrito na seção anterior para analisar brevemente três modelos baseados em equações diferenciais fracionárias. O primeiro deles é o oscilador harmônico fracionário submetido a uma força externa, outro é a versão fracionária da equação logística e, por fim, trataremos de um modelo de viscoelasticidade para tecidos moles. Esperamos que o tratamento aqui apresentado contribua de forma significativa nas pesquisas sobre cálculo fracionário e sirva de referência literária introdutória na área.
5.1. O oscilador harmônico fracionário submetido a uma força externa senoidal
O oscilador harmônico é um sistema físico constituído por um objeto de massa m preso a uma mola de constante k. A equação diferencial que governa o movimento do oscilador, isto é, a equação que fornece a posição do oscilador em função do tempo é dada por
em que . A Eq. (34) é usada para analisar o oscilador livre. Se por acaso o oscilador estiver submetido a alguma força externa de módulo Fext, podemos escrever
É relevante aqui recordarmos que o oscilador harmônico é um modelo que descreve bem muitos sistemas físicos reais, como por exemplo sistemas que oscilam em torno de uma posição de equilíbrio. Ademais, se considerarmos um oscilador com dissipação, isto é, se as forças de atrito forem consideradas, obtemos os casos especiais do oscilador amortecido, superamortecido e subamortecido, cada qual com a sua aplicação prática. Um fato bastante peculiar no caso do oscilador harmônico fracionário é a possibilidade de se descrever dissipação sem a necessidade de introduzir o termo associado à força de atrito na equação do movimento; bastando para isso variar a ordem da derivada fracionária. Além disso, podemos generalizar o oscilador para um tratamento quântico, submetendo a equação de Schrödinger ao potencial do oscilador harmônico, o que nos permite descrever por exemplo o comportamento de moléculas diatômicas com bastante acurácia. Porém, o problema do oscilador harmônico no âmbito da mecânica quântica ainda é um problema muito pouco explorado.
O correspondente fracionário do modelo físico dado na Eq. (35) quando submetido a uma força externa senoidal é dado por [3[3] R.F. Camargo e E.C. Oliveira, Cálculo fracionário (Livraria da Física, São Paulo, 2015)., 37[37] G.S. Teodoro, D.S. Olievira, C. Oliveira, Rev. Bras. de Ensino de Física 40, e2307 (2018).]
em que e Fo é a amplitude da força externa. Observe que quando duas condições iniciais são necessárias para a solução do problema, a saber: e . Observa-se que uma interpretação física para o oscilador harmônico fracionário foi dada por Stanislavsky [23[23] A.A. Stanislavsky, Physical Review E 70, 051103 (2004).]. Ele sugeriu que um oscilador harmônico fracionário deve ser interpretado como uma média de osciladores harmônicos comuns (de derivada ordinária) governado por uma seta do tempo estocástica. Além disso, do ponto de vista matemático os osciladores harmônicos fracionários também são úteis para determinar a solução da equação de Schrödinger fracionária livre [7[7] R. Herrmann, Fractional Calculus: An Introduction to Physicists (World Scientific, Singapura (2011)., 24[24] X. Guo e M. Xu, Journal of Mathematical Physics 47, 082104 (2006).]. Ressalta-se ainda que aspectos fisicamente significativos do sistema dinâmico proporcionado pela equação do oscilador harmônico fracionário, como energia total e resposta dinâmica do sistema, também já foram estudados na literatura [25[25] B.N.N. Achar, J.W. Hanneken, T. Enck e T. Clarke, Physica A 309, 275 (2002)., 26[26] B.N.N. Achar, J.W. Hanneken, T. Enck e T. Clarke, Physica A 297, 361 (2001)., 27[27] A. Al-rabtah, V.S. Ertürk e S. Momani, Computers & Mathematics with Applications 59, 1356 (2010).].
A solução da Eq. (36) para o caso livre, isto é, F0 = 0, pode ser obtida de maneira analítica diretamente a partir da transformada de Laplace, considerando as Eq. (11) e Eq. (14), e é dada por
No entanto, o mesmo não se pode dizer sobre a solução analítica para o caso em que . Ainda é possível encontrar solução analítica para casos particulares da Eq. (36) utilizando a transformada de Laplace e o teorema dos resíduos, mas a solução está longe de ser encontrada trivialmente. Por exemplo, para , , F0 = 1 e é possível mostrar que a solução analítica da Eq. (36) é
Nesse sentido, trataremos aqui a Eq. (36) utilizando o método numérico desenvolvido na seção anterior. Nosso objetivo aqui é duplo: no primeiro momento estudaremos o que ocorre no modelo fracionário quando , e em seguida, veremos como o sistema fracionário se comporta no caso em que a frequência natural coincide com a da força externa, isto é, .
Primeiramente, vamos avaliar a qualidade da aproximação numérica a comparando à solução analítica dada pela Eq. (38). Na Figura 3 podemos ver tal comparação e qualitativamente é possível observar que a solução obtida com o método numérico descrito se aproxima satisfatoriamente da solução analítica.
Solução analítica e numérica da equação no intervalo com condições iniciais e . Para este problema a solução analítica é dada pela Eq. (38). Além disso, para a solução numérica foi utilizado o passo .
Na Figura 4 podemos ver o caso do oscilador harmônico para e para diferentes valores de (; ; e ). Percebemos que quanto menor é a ordem da derivada fracionária, menor é a amplitude da oscilação. E ainda, percebemos também que neste caso há uma espécie de amortecimento da oscilação, mesmo não tendo considerado ações de forças dissipativas.
Solução numérica do oscilador harmônico fracionário no intervalo com condições iniciais e para diferentes valores de . Além disso, para a solução numérica foi utilizado o passo .
Finalmente, na Figura 5 podemos ver o caso do oscilador harmônico para para diferentes valores de (; ; e ). Conforme já conhecido da literatura, quando as frequência da fonte externa é igual à frequência natural de oscilação do sistema, no cálculo de ordem inteira, temos o fenômeno da ressonância, ou seja, a amplitude vai crescendo indefinidamente. Contudo, notamos no gráfico e também no exemplo dado pela equação (38), que enquanto a ordem da derivada fracionária diminui, a oscilação se estabiliza, isto é, além da amplitude do movimento oscilatório ser menor, a ressonância deixa de acontecer. No caso de oscilações forçadas de um sistema mecânico, infere-se deste resultado que mesmo uma força externa com mesma frequência não levará o sistema a se romper.
Solução numérica do oscilador harmônico fracionário no intervalo com condições iniciais e para diferentes valores de . Além disso, para a solução numérica foi utilizado o passo .
5.2. Equação logística fracionária
A equação logística
em que k é uma constante, foi proposta originalmente por Verhulst em 1838 [28[28] P.F. Verhulst, Correspondance mathématique et physique 59, 113 (1838).]. Em seu trabalho, Verhulst tentou aprimorar o modelo de Malthus para a estatística populacional [29[29] R.C. Bassanezzi, Ciência e Natura 36, 97 (2014).]. Desde então, a equação logística passou a ser aplicada em diversos contextos, desde modelos biológicos de proliferação de bactérias até predição de sistemas econômicos [30[30] A.M.A. El-Sayed, A.E.M. El-Mesiry e H.A.A. El-Saka, Applied Mathematics Letters 20, 817 (2007)., 31[31] R.F. Camargo, M.M. Theodoro, Revista Eletrônica Paulista de Matemática 17, 61 (2020).].
Uma hipótese básica relacionada à Eq. (39) é que, para o processo de variação da população, o tempo é uma quantidade passiva homogênea que indexa os eventos de mudança da quantidade por meio do tempo regularmente medido por um relógio. No entanto, este não é necessariamente o caso para fenômenos complexos e, assim, a derivada temporal é continuada analiticamente para uma derivada fracionária para produzir a equação logística de ordem fracionária [32[32] B.J. West, Physica A: Statistical Mechanics and its Applications 429, 103 (2015)., 30[30] A.M.A. El-Sayed, A.E.M. El-Mesiry e H.A.A. El-Saka, Applied Mathematics Letters 20, 817 (2007).]
em que k é uma constante, e com uma condição inicial . Cabe ressaltar que devido aos efeitos de memória e também a não-localidade das derivadas de Caputo, diversos modelos populacionais da literatura, e não apenas a equação logística, recentemente foram estendidos para suas respectivas versões fracionárias com o intuito de avaliar melhor a evolução de doenças [33[33] I. Ahmed, I.A. Baba, A. Yusuf, P. Kumam, W. Kumam, Advances Differences Equations, 394, 1 (2020)., 34[34] M.A.A. Oud, A. Ali, H. Alrabaiah, S. Ullah, M.A. Khan, A.A. Islam, Advances Differences Equations, 106, 1 (2021)., 35[35] F. Alcántara-López, C. Fuentes, C. Chávez, F. Brambila-Paz, A. Quevedo, Mathematics 9 1915 (2021).].
Solução numérica da equação logística fracionária no intervalo para diferentes valores de . (a) Solução com condição inicial ; (b) solução com condição inicial . Além disso, para a solução numérica foi utilizado o passo (a) ; e (b) .
Note que diferentemente do oscilador harmônico a Eq. (40) é não-linear e portanto os métodos envolvendo a transformada de Laplace dada pela Eq. (14) não podem ser aplicados. A resolução desta equação pode aparecer na literatura de algumas formas distintas: o tratamento de uma versão linearizada [31[31] R.F. Camargo, M.M. Theodoro, Revista Eletrônica Paulista de Matemática 17, 61 (2020).], solução analítica através de expansão em algum tipo de série [32[32] B.J. West, Physica A: Statistical Mechanics and its Applications 429, 103 (2015).], solução analítica através de alguma hipótese específica sobre a forma da solução [36[36] V.E. Tarasov, Mathematics 8, 2231 (2020).] e abordagem via métodos numéricos [30[30] A.M.A. El-Sayed, A.E.M. El-Mesiry e H.A.A. El-Saka, Applied Mathematics Letters 20, 817 (2007).]. A nossa abordagem aqui será realizada com o uso do método construído na seção precedente.
Na Figura 6(a) podemos ver as soluções da equação logística para diferentes valores de (; ; ; e ) com a condição inicial e . Vemos assim que quanto menor a ordem da derivada fracionária, a solução da equação logística fracionária cresce de forma mais suave e se estabiliza num valor mais baixo, isto é, o valor suporte é menor. Se tomarmos como exemplo o crescimento de uma população, se o problema puder ser modelado via cálculo fracionário, a população que cresce mais devagar e se estabiliza num valor menor é modelada por uma derivada de ordem fracionária mais baixa.
Por outro lado, na Figura 6(b) vemos as soluções da equação logística para diferentes valores de (; ; ; e ) com a condição inicial e . Note que neste caso o passo é menor que o exemplo anterior para garantir uma solução suave para o caso . Novamente podemos ver que há uma estabilização mais suave da solução quando a ordem da derivada fracionária é menor.
5.3. Modelo de viscoelasticidade de ordem fracionária para tecidos moles
Os tecidos biológicos moles podem ser considerados como materiais viscoelásticos, isto é, apresentam propriedades elásticas e viscosas ao mesmo tempo [38[38] Y.C. Fung, Biomechanics: mechanical properties of living tissues (Springer Verlag, New York, 1993).].
A lei elástica de Fung [39[39] Y.C. Fung, American Journal of Physiology 28, 1532 (1967).] é o modelo padrão-ouro na literatura para modelar a resposta elástica de tecidos moles [38[38] Y.C. Fung, Biomechanics: mechanical properties of living tissues (Springer Verlag, New York, 1993)., 11[11] A.D. Freed e K. Diethelm, Fractional Calculus and Applied Analysis 10, 219 (2007)., 12[12] T.W.J. Almeida, Sistema para análise viscoelástica de tecidos moles por ondas de cisalhamento usando excitação magnética e medida ultrassônica. Tese de Doutorado, Universidade de São Paulo, Ribeirão Preto (2015).]. Ela descreve o comportamento não-linear da tensão-deformação de um tecido biológico mole e foi deduzida a partir de observações experimentais.
A lei de Fung é descrita pela equação diferencial ordinária
em que s é a tensão Lagrangiana, o alongamento, E é uma constante com unidades de tensão e é uma constante adimensional. Ambas as constantes são características do material envolvido.
De forma não-trivial pode-se estender a lei de Fung para o modelo [11[11] A.D. Freed e K. Diethelm, Fractional Calculus and Applied Analysis 10, 219 (2007).]
que descreve as características viscoelásticas dos tecidos moles. Aqui os parâmetros materiais são constantes não-negativas tais que E tem unidades de tensão, é adimensional, e tem unidades de tempo com . Além disso, a função representa a deformação e é supostamente conhecida. No entanto, este modelo apresentado peca ao modelar com precisão dados experimentais [11[11] A.D. Freed e K. Diethelm, Fractional Calculus and Applied Analysis 10, 219 (2007).].
Já foi mostrado que as propriedades viscoelásticas desses tecidos são melhores tratadas com modelos viscoelásticos de ordem fracionária [40[40] T.C. Doehring, A.D. Freed, E.O. Carew, I. Vesely, Journal of Biomechanical Engineering 127, 700 (2005)., 41[41] A.D. Freed, K. Diethelm, Biomechanics and Modeling in Mechanobiology 5, 203 (2006).]. Desta forma, o modelo dado pela Eq. (42) é continuado analiticamente para a equação diferencial fracionária [11[11] A.D. Freed e K. Diethelm, Fractional Calculus and Applied Analysis 10, 219 (2007).]
em que .
Aqui vamos replicar o experimento do teste de tensão da referência [11[11] A.D. Freed e K. Diethelm, Fractional Calculus and Applied Analysis 10, 219 (2007).]. Com isso, iremos simular o comportamento do colágeno para uma função de deformação conhecida. Desta forma, a Eq. (43) será resolvida numericamente com os seguintes parâmetros: ; ; E=1 MPa, s e ms. O teste de tensão será feito utilizando uma deformação linear , com taxa de deformação k constante, para um material em estado virgem, isto é, . Além disso, serão considerados três valores distintos de k, a saber: , e . A simulação será feita no intervalo de tempo , com tmax escolhido tal que . Finalmente, para a função de deformação escolhida temos, pela Eq. (13), que . Com isso, o problema a ser resolvido numericamente é dado pela equação fracionária linear
Na Figura 7 podemos ver a solução para o teste de tensão com gráfico apresentado na forma tensão-deformação. Os resultados apresentados são compatíveis com resultados experimentais e correspondem aos mesmos apresentados por Freed et al. [11[11] A.D. Freed e K. Diethelm, Fractional Calculus and Applied Analysis 10, 219 (2007).]. Novamente, percebe-se a adequação do método numérico ao problema de viscoelasticidade do colágeno e a sua importância para a solução de problemas teóricos e reais.
Resultado do teste de tensão para a simulação do colágeno modelado pela Eq. (44). Aqui foram utilizados os seguintes dados: deformação linear , com , e . Cada uma das simulações foi realizada no intervalo de tempo , de tal forma que para cada valor de k a igualdade deve ser respeitada.
6. Considerações Finais e Perspectivas
Neste trabalho, apresentam-se os conceitos das equações diferenciais fracionárias no escopo de Caputo enfatizando sua resposta analítica que é deduzida pelo método da Transformada de Laplace cuja solução geral é dada em termos da função de Mittag-Leffer de um parâmetro.
Além disso, enfatizou-se um método do tipo preditor-corretor para a resolução numérica destas equações diferenciais fracionárias, na qual é uma variação do método de Adam-Bashforth, que caracterizou-se pela sua eficácia em relação aos métodos do tipo Runge-Kutta usualmente abordados na literatura para a solução numérica destas equações.
As três aplicações expostas: o oscilador harmônico fracionário submetido a uma força externa senoidal, a equação logística e o modelo de viscoelasticidade, exemplificam o funcionamento do método numérico descrito, na qual os resultados numéricos obtidos mostraram valores com convergência rápida do valor exato e com pequenos erros.
Espera-se com este trabalho contribuia para uma melhor compreensão das equações diferenciais fracionárias no escopo de Caputo e uma reflexão sobre sua importância aplicação nos fenômenos físicos, tanto de teoria (que traz uma novo conceito teórico) quanto na solução (que traz uma abordagem eficaz para seu desenvolvimento numérico).
Referências
-
[1]K. Diethelm, N.J. Ford e A.D. Freed, Nonlinear Dynamics 29, 3 (2002).
-
[2]B. Ross, Historia Mathematica 4, 75 (1977).
-
[3]R.F. Camargo e E.C. Oliveira, Cálculo fracionário (Livraria da Física, São Paulo, 2015).
-
[4]M. Caputo, Geophysical Journal of the Royal Astronomical Society 13, 529 (1967).
-
[5]K.B. Oldhan e J. Spanier, The Fractional Calculus: Theory and Aplications, Differentiation and Integration to Arbitrary Order (Academic, New York, 1974).
-
[6]N. Varalta, A.V. Gomes e R.F. Camargo, TEMA 15, 211 (2014).
-
[7]R. Herrmann, Fractional Calculus: An Introduction to Physicists (World Scientific, Singapura (2011).
-
[8]P.J. Torvik e R.L. Bagley, Journal of Applied Mechanics 51, 94 (1984).
-
[9]W. Labecca, O. Guimarães e J.R.C. Piqueira, Mathematical Problems in Engineering, 591715 (2015).
-
[10]J. Vaz Jr e E.C. Oliveira, Mathematics in Engineering 4, 1 (2022).
-
[11]A.D. Freed e K. Diethelm, Fractional Calculus and Applied Analysis 10, 219 (2007).
-
[12]T.W.J. Almeida, Sistema para análise viscoelástica de tecidos moles por ondas de cisalhamento usando excitação magnética e medida ultrassônica Tese de Doutorado, Universidade de São Paulo, Ribeirão Preto (2015).
-
[13]J.P.S. Ramírez, Função Gama, disponível em: https://www.ime.unicamp.br/~ftorres/ENSINO/MONOGRAFIAS/JP_VC2_2015.pdf (2015).
» https://www.ime.unicamp.br/~ftorres/ENSINO/MONOGRAFIAS/JP_VC2_2015.pdf -
[14]D.S. Oliveira, Derivada Fracionária e as Funções de Mittag-Leffler Dissertação de Mestrado, Universidade Estadual de Campinas, Campinas (2014).
-
[15]A.A. Kilbas, H.M. Srivastava e J.H. Trujillo, Theory and Application of Fractional Differential Equations (Elsevier, Reino Unido, 2006).
-
[16]K.M. Owolabi, A. Atangana, Numerical Methods for Fractional Differentiation (Springer, Singapura, 2019).
-
[17]Z.M. Odibat e N.T. Shawagfeh, Applied Mathematics and Computation 186, 286 (2007).
-
[18]S. ÖLrekçi, Advances in Mathematical Physics, 507970 (2015).
-
[19]M.S. Arshad, D. Baleanu, M.B. Riaz e M. Abbas, Discrete Dynamics in Nature and Society, 1020472 (2020).
-
[20]K. Diethelm, The Analysis of Fractional Differential Equations: An Application-Oriented Exposition Using Differential Operators of Caputo Type (Springer, Singapura, 2010).
-
[21]K. Diethelm, N.J. Ford e A.D. Freed, Numerical Algorithms 36, 31 (2004).
-
[22]K. Diethelm, N.J. Ford, A.D. Freed e Y. Luchko, Computer Methods in Applied Mechanics and Engineering, 194, 743 (2005).
-
[23]A.A. Stanislavsky, Physical Review E 70, 051103 (2004).
-
[24]X. Guo e M. Xu, Journal of Mathematical Physics 47, 082104 (2006).
-
[25]B.N.N. Achar, J.W. Hanneken, T. Enck e T. Clarke, Physica A 309, 275 (2002).
-
[26]B.N.N. Achar, J.W. Hanneken, T. Enck e T. Clarke, Physica A 297, 361 (2001).
-
[27]A. Al-rabtah, V.S. Ertürk e S. Momani, Computers & Mathematics with Applications 59, 1356 (2010).
-
[28]P.F. Verhulst, Correspondance mathématique et physique 59, 113 (1838).
-
[29]R.C. Bassanezzi, Ciência e Natura 36, 97 (2014).
-
[30]A.M.A. El-Sayed, A.E.M. El-Mesiry e H.A.A. El-Saka, Applied Mathematics Letters 20, 817 (2007).
-
[31]R.F. Camargo, M.M. Theodoro, Revista Eletrônica Paulista de Matemática 17, 61 (2020).
-
[32]B.J. West, Physica A: Statistical Mechanics and its Applications 429, 103 (2015).
-
[33]I. Ahmed, I.A. Baba, A. Yusuf, P. Kumam, W. Kumam, Advances Differences Equations, 394, 1 (2020).
-
[34]M.A.A. Oud, A. Ali, H. Alrabaiah, S. Ullah, M.A. Khan, A.A. Islam, Advances Differences Equations, 106, 1 (2021).
-
[35]F. Alcántara-López, C. Fuentes, C. Chávez, F. Brambila-Paz, A. Quevedo, Mathematics 9 1915 (2021).
-
[36]V.E. Tarasov, Mathematics 8, 2231 (2020).
-
[37]G.S. Teodoro, D.S. Olievira, C. Oliveira, Rev. Bras. de Ensino de Física 40, e2307 (2018).
-
[38]Y.C. Fung, Biomechanics: mechanical properties of living tissues (Springer Verlag, New York, 1993).
-
[39]Y.C. Fung, American Journal of Physiology 28, 1532 (1967).
-
[40]T.C. Doehring, A.D. Freed, E.O. Carew, I. Vesely, Journal of Biomechanical Engineering 127, 700 (2005).
-
[41]A.D. Freed, K. Diethelm, Biomechanics and Modeling in Mechanobiology 5, 203 (2006).
Datas de Publicação
-
Publicação nesta coleção
17 Jan 2022 -
Data do Fascículo
2022
Histórico
-
Recebido
09 Dez 2021 -
Aceito
15 Dez 2021