Instruções de uso dos subprogramas para resolução de problemas de valor inicial

DMPA/IM/UFRGS
PROF. RUDNEI DIAS DA CUNHA

Os subprogramas que implementam os métodos de passo simples (Euler, Heum, Runge-Kutta de 4ª ordem e Runge-Kutta-Fehlberg) e de passo múltiplo (Adams-Bashforth de 4ª ordem e Adams-Moulton de 4ª ordem) em Scilab foram escritos para resolver problemas de valor inicial escritos na forma autônoma,
{𝑋′=𝐹(𝑋) ; 𝑋(𝑥0)=𝑋0
 
onde 𝑋′={𝑓0=1 ; 𝑓1(𝑥0,𝑥1,…,𝑥𝑛) ⋮ 𝑓𝑛(𝑥0,𝑥1,…,𝑥𝑛), 𝑋=[𝑥0≡𝑡𝑥1⋮𝑥𝑛]

o que implica na necessidade do usuário adicionar, explicitamente, a equação diferencial para a variável 𝑡 (mais especificamente, 𝑓0=1).

A interface dos subprogramas é definida como segue:

a) Euler: Y=euler(F,x,h,t1)
b) Heum: Y=heum(F,x,h,t1)
c) Runge-Kutta de 4ª ordem: Y=rk4(F,x,h,t1)
d) Runge-Kutta-Fehlberg: Y=rkf45(F,x,h,t1,delta)
e) Adams-Bashforth de 4ª ordem: Y=ab4(F,x,h,t1)
f) Adams-Moulton de 4ª ordem: Y=ab4_am4(F,x,h,t1)
 
onde F especifica as funções 𝑓0,𝑓1,…,𝑓𝑛; x é um vetor contendo os valores iniciais do PVI (incluindo o valor inicial de 𝑡); h é o passo de integração; t1 é o valor final; e delta é a tolerância para o erro no método RKF.

A especificação do conjunto de funções 𝑓0,𝑓1,…,𝑓𝑛 para utilização com os subprogramas pode ser feito de duas formas:

a) Especificando-as como expressões aritméticas em “strings”, armazenadas como elementos de um vetor coluna;
b) Especificando-as num subprograma externo (com a extensão .sci), definido pelo usuário. Esse subprograma deverá receber um vetor coluna como argumento de entrada (o qual contém os valores das variáveis de interesse) e deverá retornar um outro vetor coluna (o qual conterá os valores das funções 𝑓0,𝑓1,…,𝑓𝑛 avaliados nas variáveis de interesse). Note que ambos os vetores deverão ter a mesma quantidade de elementos.

Cabe ressaltar que, em ambos os casos, as variáveis de interesse devem ser referenciadas através dos elementos de um arranjo, como se pode notar no exemplo a seguir.

Exemplo: seja o sistema de equações diferenciais ordinárias

{𝑥′=2𝑥−0,02𝑥𝑦 ; 𝑦′=0,0005𝑥𝑦−0,8𝑦

sujeito às condições iniciais 𝑥(0)=3000 e 𝑦(0)=120. Quais os valores de 𝑥 e 𝑦 em 𝑡=10, usando o método de Runge-Kutta de 4ª ordem (rk4) com ℎ=0,25.

Resolução com arranjo de “strings”:

Como o sistema deve ser escrito na forma autônoma, devemos considerar 𝑡 como uma variável do problema, e adicionar ao sistema original a equação 𝑓0=1. Assumindo que as variáveis 𝑡, 𝑥 e 𝑦 serão armazenadas nos elementos 1, 2 e 3 do vetor coluna x, podemos escrever a chamada ao subprograma rk4 como

F = [ '1'; '2*x(2)-0.02*x(2)*x(3)'; '0.0005*x(2)*x(3)-0.8*x(3)' ]
X0 = [ 0; 3000; 120 ]
Y = rk4(F, X0, 0.25, 10.)
O resultado fornecido pelo subprograma rk4 é uma matriz Y que conterá tantas linhas quantas foram as iterações realizadas pelo método; a quantidade de colunas será igual a quantidade de variáveis de interesse, incluindo a variável 𝑡.

Resolução com subprograma externo:

Novamente assumindo que as variáveis 𝑡, 𝑥 e 𝑦 serão armazenadas nos elementos 1, 2 e 3 do vetor coluna x, pode-se escrever um subprograma denominado presa_predador, como segue:

function F=presa_predador(var)
 // Para facilitar a escrita das equações,
 // utiliza os nomes das variáveis
 x = var(2)
 y = var(3)
 // Escreve as equações;
 // F(1)=1 sempre, pois o PVI deve ser escrito na
 // forma autônoma
 F(1) = 1
 F(2) = 2*x-0.02*x*y
 F(3) = 0.0005*x*y-0.8*y
endfunction
Uma vez tendo carregado o subprograma presa_predador no Scilab, a chamada ao subprograma rk4 pode ser feita, então, como

X0 = [ 0; 3000; 120 ]
Y = rk4(presa_predador, X0, 0.25, 10.)