ode scilab

Esta é uma pré-visualização de arquivo. Entre para ver o arquivo original

ode
Solucionador de equações diferenciais ordinárias
Seqüência de Chamamento
y = ode(y0,t0,t,f)
[y,w,iw] = ode([type],y0,t0,t [,rtol [,atol]],f [,jac] [,w,iw])
[y,rd,w,iw] = ode("root",y0,t0,t [,rtol [,atol]],f [,jac],ng,g [,w,iw])
y = ode("discrete",y0,k0,kvect,f)
Parâmetros
y0
matriz ou vetor de reais (condições iniciais).
t0
escalar real (tempo inicial).
t
vetor de reais (tempos nos quais a solução é computada).
f
função externa (função, lista ou string).
type
um dos strings seguintes: "adams", "stiff", "rk", "rkf", "fix", "discrete", "roots".
rtol, atol
constantes reais ou vetores com o mesmo tamanho que y.
jac
função externa (função, lista ou string).
ng
inteiro
g
função externa (função, lista ou string).
k0
inteiro (tempo inicial).
kvect
vetor de inteiros.
y
a real vector or matrix. The solution.
rd
a real vector
w, iw
vetores de reais. Ver ode() optional output
Descrição
ode é a função padrão para se resolver sistemas de EDO explícitos definidos por: dy/dt=f(t,y) , y(t0)=y0. É uma interface entre vários
solucionadores, em particular a ODEPACK. O tipo de problema resolvido e o método a ser utilizado dependem do primeiro argumento
opcional type que pode ser um dos strings seguintes:
<not given>:
O solucionador lsoda do pacote ODEPACK é chamado por padrão. Ele escolhe automaticamente entre o método preditor-corretor
não-rígido de Adams e a Fórmula de Diferenciação Retroativa (FDR) rígida. Ele utiliza o método não rígido inicialmente e monitora
ode - Solucionador de equações diferenciais ordinárias https://help.scilab.org/docs/6.0.0/pt_BR/ode.html
1 de 4 28/11/2017 12:35
os dados para decidir qual método utilizar.
"adams":
Este é para problemas não rígidos. O solucionador lsode do pacote ODEPACK é chamado e utiliza o método de Adams.
"stiff":
Este é para problemas rígidos. O solucionador lsode do pacote ODEPACK é chamado e é utilizado o método FDR.
"rk":
Método adaptativo de Runge-Kutta de ordem 4 (RK4).
"rkf":
O programa de Shampine e Watts baseado no par Runge-Kutta de Fehlberg de ordem 4 e 5 (RKF45) é utilizado. Este é utilizado
para problemas não-rígidos e mediamente rígidos quando as computações de derivação não são custosas. Este método não é
recomendado ser utilizado quando o usuário requer uma maior precisão.
"fix":
Mesmo solucionador que "rkf", mas a interface do usuário é bem simples, i.e. apenas os parâmetros rtol e atol podem ser
passados ao solucionador. Este é o método mais simples a se tentar.
"root":
Solucionador de EDOs com capacidade de encontrar raízes. O solucionador lsodar do pacote ODEPACK é utilizado. É uma
variante do solucionador lsoda onde se acha raízes de uma dada função vetorial. Ver ajuda em ode_roots para mais detalhes.
"discrete":
Simulação de tempo discreto. Ver ajuda em ode_discrete para mais detalhes.
Nesta ajuda podemos apenas descrever o uso de ode para sistemas padrões explícitos EDOs .
A chamada mais simples a ode é: y=ode(y0,t0,t,f) onde y0 é o vetor de condições iniciais, t0 é o tempo inicial, t é o vetor
de tempos onde a solução y é computada e y é a matriz de vetores soluções y=[y(t(1)),y(t(2)),...].
O argumento de entrada f define o lado direito da equação diferencial de primeira ordem: dy/dt=f(t,y). É uma função externa, isto
é, uma função com sintaxe especificada, ou o nome de uma subrotina Fortran ou uma subfunção C (string) com seqüência de
chamamento especificada, ou uma lista:
Se f for uma função do Scilab, a sua sintaxe deve ser ydot = f(t,y), onde t é um escalar real (tempo) e y é um vetor
de reais (estado). ydot é um vetor de reais (dy/dt)
Se f é um string ele se refere ao nome de uma subrotina Fortran ou uma subfunção C, i.e. Se ode(y0,t0,t,"fex") for
o comando, então a subrotina fex será chamada.
A rotina Fortran deve ter a seguinte seqüência de chamamento: fex(n,t,y,ydot), com n um inteiro, t um escalar de
dupla precisão, y e ydot vetores de dupla precisão.
A função C deve ter o seguinte protótipo: fex(int *n,double *t,double *y,double *ydot)
t é o tempo, y o estado e ydot a derivada do estado (dy/dt)
Esta função externa pode ser construída em um SO de modo indpendente através de ilib_for_link e ligada
dinamicamente através da função do Scilab link.
O argumento f também pode ser uma lista com a seguinte estrutura: lst=list(realf,u1,u2,...un) onde realf é
uma função do Scilab com a sintaxe: ydot = f(t,y,u1,u2,...,un)
Esta sintaxe permite utilizar parâmetros como argumentos de realf.
A função f pode retornar uma matriz p x q ao invés de um vetor. Com esta notação de matriz, nos resolvemos o sistema n=p+q
de EDOs dY/dt=F(t,Y) onde Y é uma matriz p x q matrix. Então, as condições iniciais, Y0, também devem ser uma matriz p x
q eo resultado de ode é a matriz p x q(T+1) [Y(t_0),Y(t_1),...,Y(t_T)].
Argumentos de entrada opcionais podem ser fornecidos para o erro da solução: rtol e atol são limiares para os erros relativos e
absolutos estimados. O erro estimado em y(i) é: rtol(i)*abs(y(i))+atol(i)
e uma intergração é feita enquanto este erro é pequeno para todos os componentes do estado. Se rtol e/ou atol for uma
constante rtol(i) e/ou atol(i) são ajustados para esta constante. Os valores padrões para rtol e atol são respectivamente
rtol=1.d-5 e atol=1.d-7 para a maior parte dos solucionadores rtol=1.d-3 e atol=1.d-4 para "rfk" e "fix".
Para problemas rígidos, é melhor fornecer o Jacobiano da função do lado direito da equação como o argumento opcional jac. É
uma função externa, istoé i.e. uma função com sintaxe especificada, ou o nome de uma subrotina Fortran ou uma subfunção C
(string) com seqüência de chamamento especificada, ou uma lista.
Se jac for uma função, a sintaxe deve ser J=jac(t,y)
ode - Solucionador de equações diferenciais ordinárias https://help.scilab.org/docs/6.0.0/pt_BR/ode.html
2 de 4 28/11/2017 12:35
onde t é um escalar real (tempo), y é um vetor de reais (estado). A matriz resultante J deve fornecer df/dx i.e. J(k,i) =
dfk/dxi com fk = k-ésimo componente de f.
Se jac for um sting, ele se refere a uma subrotina Fortran, ou uma subfunção C, com as seguinte seqüência de chamamento:
No caso Fortran:
subroutine fex(n,t,y,ml,mu,J,nrpd)
integer n,ml,mu,nrpd
double precision t,y(*),J(*)
No caso C:
void fex(int *n,double *t,double *y,int *ml,int *mu,double *J,int *nrpd,)
jac(n,t,y,ml,mu,J,nrpd). Na maior parte dos casos, você não tem que se referir a ml, mu e nrpd.
Se jac for uma lista, as mesmas convenções que para f se aplicam.
Os argumentos opcionais w e iw e vetores para armazenamento de informações podem ser retornados pela rotina de integração
(ver ode_optional_output para detalhes). Queando esses vetores são fornecidos no lado direito de ode, a integração reinicia com
os mesmos parâmetros da parada anteiror.
Mais opções podem ser fornecidas a solucionadores ODEPACK utilizando-se a variável %ODEOPTIONS. Ver odeoptions.
Exemplos
// ---------- EDO simples de uma dimensão (função externa do Scilab)
// dy/dt=y^2-y sin(t)+cos(t), y(0)=0
function ydot=f(t, y),ydot=y^2-y*sin(t)+cos(t),endfunction
y0=0;t0=0;t=0:0.1:%pi;
y=ode(y0,t0,t,f)
plot(t,y)
// ---------- EDO simples de uma dimensão (função externa codificada em C)
ccode=['#include <math.h>'
'void myode(int *n,double *t,double *y,double *ydot)'
'{'
' ydot[0]=y[0]*y[0]-y[0]*sin(*t)+cos(*t);'
'}']
mputl(ccode,TMPDIR+'/myode.c') //criando o arquivo C
ilib_for_link('myode','myode.c',[],'c',[],TMPDIR+'/loader.sce');//compilando
exec(TMPDIR+'/loader.sce') //linking incremental
y0=0;t0=0;t=0:0.1:%pi;
y=ode(y0,t0,t,'myode');
// ---------- Simulação de dx/dt = A x(t) + B u(t) com u(t)=sin(omega*t),
// x0=[1;0]
// solução x(t) desejada em t=0.1, 0.2, 0.5 ,1.
// A e a função u são passados para a função do lado direito em uma lista.
// B e omega são passados como
variáveis globais
function xdot=linear(t, x, A, u),xdot=A*x+B*u(t),endfunction
function ut=u(t),ut=sin(omega*t),endfunction
A=[1 1;0 2];B=[1;1];omega=5;
ode([1;0],0,[0.1,0.2,0.5,1],list(linear,A,u))
// ---------- Integração com notação matricial da equação diferencial de Ricatti
// Xdot=A'*X + X*A - X'*B*X + C , X(0)=Identity
// Solução em t=[1,2]
function Xdot=ric(t, X),Xdot=A'*X+X*A-X'*B*X+C,endfunction
A=[1,1;0,2]; B=[1,0;0,1]; C=[1,0;0,1];
t0=0;t=0:0.1:%pi;
X=ode(eye(A),0,t,ric)
// ---------- Notação matricial, computação de exp(A)
A=[1,1;0,2];
function xdot=f(t, x),xdot=A*x;,endfunction
ode(eye(A),0,1,f)
ode("adams",eye(A),0,1,f)
// ---------- Notação de matriz, computação de exp(A) com matriz rígida, Jacobian fornecida
A=[10,0;0,-1];
function xdot=f(t, x),xdot=A*x,endfunction
function J=Jacobian(t, y),J=A,endfunction
ode("stiff",[0;1],0,1,f,Jacobian)
Ver Também
odeoptions — ajusta opções para solucionadores de EDO
ode_optional_output — descrição de saídas opcionais de solucionadores de EDOs
ode - Solucionador de equações diferenciais ordinárias https://help.scilab.org/docs/6.0.0/pt_BR/ode.html
3 de 4 28/11/2017 12:35
ode_roots — solucionador de EDOs com busca de raízes
ode_discrete — solucionador de equações diferenciais ordinárias, simulação de tempo discreto
dassl — Equação diferencial algébrica
impl — equações diferenciais algébricas
odedc — solucionador de EDOs contÃnuas/discretas
csim — simulation (time response) of linear system
ltitr — discrete time response (state space)
rtitr — discrete time response (transfer matrix)
intg — integral definida
Bibliografia
Alan C. Hindmarsh, lsode and lsodi, two new initial value ordinary differential equation solvers, acm-signum newsletter, vol. 15, no. 4 (1980),
pp. 10-11.
Funções Utilizadas
As rotinas associadas podem ser encontradas no diretório SCI/modules/differential_equations/src/fortran:
lsode.f lsoda.f lsodar.f
ode - Solucionador de equações diferenciais ordinárias https://help.scilab.org/docs/6.0.0/pt_BR/ode.html
4 de 4 28/11/2017 12:35