Octave: Směrové pole ODR

Směrové pole

Definuji vektorovou funkci:

f = @(x,y) 2*x/(3*x.^2 + y);

Určím hrubost zobrazování vektorů h a vytvořím mřížku (ta se vytiskne až po zadání příkazu quiver):

h = 0.2;
[x, y] = meshgrid(0:h:2, -1:h:1);

Numericky derivuji funkci f(x,y) (Aproximuji malý úsek křivky na nějak nakloněnou (=> nějak dlouhou) úsečku. Pro posun po ose x o krok 1 přiřazuji funkční hodnotu s. Celková délka úsečky D = velikost vektoru (1,s):

s = f(x,y);
D = sqrt(1 + s.^2);

Příkaz quiver (do češtiny lze přeložit jako toulec - tj. obsahuje uvnitř šípy/šipky) slouží k vykreslení pole (šipek) vektorů (quiver(X,Y,U,V, velikost-sipky), axis tight přilne osy přímo k vektorovému poli.

quiver(x,y, 1./D,s./D, 0.3),axis tight;

Celý kód:

f = @(x,y) -(0.5*x-1)./y.^2;
h = 0.2;
[x,y] = meshgrid(0:h:2, -1:h:1);
s = f(x,y);
D = sqrt(1 +s.^2);
quiver(x,y,1./D,s./D, 0.3);
axis tight;
xlabel ' osa x', ylabel 'osa y';
title 'Směrové pole pro dy/dx = -(0.5*x-1)/y^2)';

image.png

Pro zobrazení chybějících hodnot mimo definiční obor:

f = @(x,y) -x./y;
h = 0.2;
[x,y] = meshgrid(0:h:2, -1:h:1);
s = f(x,y);
D = sqrt(1 +s.^2);
hold on;

quiver(x,y,1./D,s./D, 0.3, color='g');
X = D==Inf;
A = D*0;
A(:) = 150;
G = X*A;

axis tight;
hold on;
quiver(x,y,G.*0.01, G, 0.1, Color="r");
hold on;
xlabel ' osa x', ylabel 'osa y';
title 'Směrové pole pro dy/dx = -x/y)';

image.png