Laboratorijske vaje Numeriˇ cne metode
2. Vaja
B. Jurˇciˇc Zlobec
Numeriˇcne metode FE, Ljubljana, 15. november 2013
Funkcije v Octave
Definiraj novo funkcijoprvaf() s pomoˇcjo funkcijske datoteke.
Funkcija je podana z naslednjim predpisom:
prvaf(x) =
πex, x <0
π+ sin(πx), 0≤x < 12 π+ 1−(x−12)2, 12 ≤x < 12+√
π+ 1
0, drugod
(1)
Izraˇcunaj funkcijske vrednosti na vektorjux=[1.2,3.4].
Funkcijska datoteka
Definicijo funkcije zapiˇsemo na datotekoprvaf.m.
Funkcijska datoteka se zaˇcne s stavkom function.
Vsebina datotekeprvaf.m:
function y=prvaf(x);
%y=prvaf(x)
%vektor x je vrstica, rezultat je vrstica y.
y=[];
for xi=x,
if xi<0, yi=pi*exp(xi);
elseif xi<1/2, yi=pi+sin(pi*xi);
elseif xi<1/2+sqrt(pi+1), yi=pi+1-(xi-1/2)^2;
else yi=0; end;
y=[y,yi];
end;
Preiskus
Uporaba ukazahelp help prvaf
y=prvaf(x)
vektor x je vrstica, rezultat je vrstica y.
Raˇcunanje funkcijskih vrednosti x=[1.2,2.4];
y=prvaf(x);
printf(’[%0.5f,%0.2f]=...
prvaf([%0.5f,%0.2f])\n’,y,x);
[3.65159,0.53]=prvaf[1.20000,2.40]
Graf funkcije
Uporabilinspacein plotter nariˇsi graf na intervalu [−3,3].
x=linspace(-3,3); y=prvaf(x); plot(x,y);
-3 -2 -1 1 2 3
1 2 3 4
Eulerjeva transformacija
S pomoˇcjo Eulerjeve transformacije doloˇci pribliˇzno vrednost vsote slabo konvergentne alternirajoˇce vrste, katere ˇcleni po absolutni vrednosti enakomerno padajo proti 0.
S =
∞
X
i=0
(−1)i
√i + 1 (2)
Ce vzamemoˇ n-to delne vsoto, je napaka razmeroma velika.
Sn=
n
X
i=0
(−1)i
√i + 1, n= 0, . . . ,N∈N.
|S−Sn|< 1
√n+ 2
Povpreˇ cenje delnih vsot
Natanˇcnejˇsi rezultat dobimo s povpreˇcenjem delnih vsot (Eulerjeva transformacija).
Sn(0)=
n
X
i=0
(−1)i
√i+ 1, n = 0, . . . ,N ∈N
Sn(1)= (Sn+1(0) +Sn(0))/2, n = 0, . . . ,N−1,
Sn(2)= (Sn+1(1) +Sn(1))/2, n = 0, . . . ,N−2, . . .
Sn(k)= (Sn+1(k−1)+Sn(k−1))/2, n = 0, . . . ,N−k.
Izbrali bomoN= 40 in k = 30, ter primerjali z rezultatom, ki ga dobimo s programomMathematicaR.
Program
Uporabili bomo funkcijocumsum(),kumulativna vsotain kontekstni parameterend.
Pri risanju grafa smo uporabili ˇseones()in size().
n=0:40;
s=(-1).^n./sqrt(n+1);
ss=cumsum(s); si=ss;
for i=1:30
si=(si(1:end-1)+si(2:end))/2;
end;
printf(’Delna vsota vrste S40=\%0.16f\n’,ss(end));
printf(’Vsota vrste S=\%0.16f\n’,si(end));
plot(n,ss,’b’,n,si(end)*ones(size(n)),’r’);
Delna vsota vrste S40=0.682509473280156 Vsota vrste S=0.604898643421630
Graf delnih vsot
10 20 30 40 50 60
0.3 0.4 0.5 0.6 0.7 0.8 0.9
Modroso delne vsote,zelenoje prvo povpreˇcenje inrdeˇceje limitna vrednost.
S pomoˇcjo programaMathematicaR dobimo naslednji rezultat:
S = 0.6048986434216303702472659142359554997598
Vektorske funkcije
Funkcijesum(),max(),min()in norm()in dot().
x=[4,2,3,6,5,4,7,6];
y=[2,3,-4,5,2,-4,-3,1];
printf(’max(x)=%d, min(x)=%d\’,max(x),min(x));
[max,maxi]=max(x); [mix,mixi]=min(x);
printf(’max index=%d, min index=%d\n’,maxi,mixi);
printf(’sum(x)=%d, dot(x,y)=%d sum(x.*y)=%d\n’,...
sum(x),dot(x,y),sum(x.*y));
max(x)=7, min(x)=3
max index=7, min index=3
sum(x)=37, dot(x,y)=13 sum(x.*y)=13
Vektorske norme
1 Druga ali Evklidska norma kxk2 = q
x12+x22+· · ·+xn2.
2 Prva ali Manhattanska norma kxk1 =|x1|+|x2|+· · ·+|xn|.
3 Neskonˇcna normakxk∞= max
1≤i≤n{|xi|}.
printf(’norm(x,2)=%0.4f, sqrt(sum(x.*x))=%0.4f, ’,...
norm(x,2),sqrt(sum(x.*x)));
printf(’norm(y,1)=%0.4f, sum(abs(y)=%0.4f\n’,...
norm(y,1),sum(abs(y))));
printf(’norm(y,inf)=%0.2f, max(abs(y))=%0.2f\n’,...
norm(y,inf),max(abs(y)));
norm(x,2)=13.67, sqrt(sum(x.*x))=13.67 norm(y,1)=24.00, sum(abs(y)=24.00 norm(y,inf)=5.00, max(abs(y))=5.00
Enotske krogle vektorskih norm in Mahattanska razdalja
Modulska aritmetika
Operatorjamod() inrem(). Uporabili smo ˇse funkcijo length().
x=[12,4,17,5,2,7,9,13]
formatn=[]; for i=1:length(x),...
formatn=[formatn,’%d ’]; end;
format=[’mod(x,%d)=[ ’,formatn,’]\n’];
n=7; printf(format,n, mod(x,n));
n=4; printf(format,n, mod(x,n));
n=-4;printf(format,n, mod(x,n));
format=[’rem(x,%d)=[ ’,formatn,’]\n’];
n=-4;printf(format,n, rem(x,n));
mod(x,7)=[ 5 4 3 5 2 0 2 6 ] mod(x,4)=[ 0 0 1 1 2 3 1 1 ]
mod(x,-4)=[ 0 0 -3 -3 -2 -1 -3 -3 ] rem(x,-4)=[ 0 0 1 1 2 3 1 1 ]
Zaokroˇ zevanje
Funkcijefloor(),ceil(),round()in fix().
x=[-2.6, -1.4, 2.8, 3.2]
formatn=[]; for i=1:length(x),...
formatn=[formatn,’%d ’]; end;
printf([’floor(x)=[ ’,formatn,’]\n’],floor(x));
printf([’ceil(x)=[ ’,formatn,’]\n’], ceil(x));
printf([’round(x)=[ ’,formatn,’]\n’], round(x));
printf([’fix(x)=[ ’,formatn,’]\n’], fix(x));
floor(x)=[ -3 -2 2 3 ] ceil(x)=[ -2 -1 3 4 ] round(x)=[ -3 -1 3 3 ] fix(x)=[ -2 -1 2 3 ]
Logiˇ cni in relacijski operatorji
Logiˇcni operatorji:
or()infiksno |,and() infiksno&, not()prefiksno! ali~in xor().
Logiˇcni konstanti:
true->1in false->0.
Relacijski operatorji:
enako==, veˇcje >, manjˇse <,
veˇcje ali enako >=, manjˇse ali enako<= in ni enako!= oziroma~=.
Naj omenimo na tem mestu ˇse naslednje funkcije:
find(),all()in any().
Primer
Komponente vektorjax so ˇcleni zaporedja:
y0 = 0.5, a= 3.9, in n= 10.
yn=ayn−1(1−yn−1), xn=bn ync, n = 1, . . . ,n.
1 Poiˇsˇci najveˇcjo in najmanjˇso komponento vektorjax.
2 Doloˇci indeks najveˇcjih in najmanjˇsih komponent vektorja.
3 Koliko komponent je veˇcjih ali enakih 3 in manjˇsih od 7?
4 Koliko komponent je deljivih s 3?
5 Kje se nahajajo komponente, ki so deljive z dva?
6 Koliko komponent je enakih 9?
7 Koliko je lihih komponent?
8 Kje se nahajajo komponente, ki so enake 9 ali 8?
Reˇsitve
x=[]; y=1/2; a=3.9;
for i=1:10, y=a*y*(1-y); x=[x,floor(10*y)]; end; x printf(’Naloga 1. min(x)=%d, max(x)=%d\n’,...
min(x),max(x));
printf(’Naloga 2.\n’);
find(x==min(x)), find(x==max(x))
printf(’Naloga 3. stevilo=%d\n’, sum(3<=x & x < 7));
printf(’Naloga 4. stevilo=%d\n’, sum(mod(x,3)==0));
printf(’Naloga 5.\n’);
find(mod(x,2)==0)
printf(’Naloga 6. stevilo=%d\n’, sum(x==9));
printf(’Naloga 7. stevilo=%d\n’, sum(mod(x,2)!=0));
printf(’Naloga 8.\n’);
find(x==9 | x==8)
Izpisi
x = 9 0 3 8 4 9 1 4 9 1
Naloga 1. min(x)=0, max(x)=9 Naloga 2.
ans = 2
ans = 1 6 9 Naloga 3. stevilo=3 Naloga 4. stevilo=5 Naloga 5.
ans = 2 4 5 8 Naloga 6. stevilo=3 Naloga 7. stevilo=6 Naloga 8.
ans = 1 4 6 9