Forum Ogólnotematyczne / Bez rejestracji

mat · Gość

DOWOLNE WĘZŁY

clear;
funcprot(0);
clf; a=-1;
disp("a = "+string(a));
b=1;
disp("b = "+string(b));
n=11;
disp("n = "+string(n));
m=n+1;
disp("m = "+string(m));
bt1=-0.32
al1=1.6
function tri=tr_inv(t);
tri=((b-a)*t+(a+b))/2;
endfunction;
xi=zeros(m,1)
xi(1)=a;
for i=0:(n-2)
xi(i+2)=tr_inv(-cos(((2*i+1)*%pi)/(2*(n- 1))));
end;
xi(m)=b;
disp(xi,"xi=");
function f=f(x);
f=1/(1+25*(x^2));
endfunction;
yi=zeros(m,1);
for i=1:m yi(i)=f(xi(i));
end;
disp(yi,"yi=");
hi=zeros(n);
for i=1:n
hi(i)=xi(i+1)-xi(i);
end;
disp(hi,"hi=");
di=zeros(n,1);
for i=2:(n-1)

di(i)=(6/(hi(i)+hi(i+1)))*((yi(i+2)-yi(i+1))/hi(i+1)-(yi(i+1)-yi(i))/hi(i));
end;
di(n)=(6/(hi(n)))*(bt1-((yi(n)-yi(n-1))/hi(n))); di(1)=6/hi(1)*((yi(2)-yi(1))/hi(i+1)-al1)

disp(di,"di=");
lai=zeros(n);
for i=1:(n-1)
lai(i)=hi(i+1)/(hi(i)+hi(i+1));
end;
lai(n)=hi(1)/(hi(1)+hi(n));
disp(lai,"lai=");
mii=zeros(n);
for i=1:n
mii(i)=1-lai(i);
end;
disp(mii,"mii=");
mm=zeros(m,m);
mm(1,1)=2;
mm(1,n+1)=1;
for i=1:(n-1)
mm(i+1,i)=mii(i);
mm(i+1,i+1)=2;
mm(i+1,i+2)=lai(i);
end;

mm(n+1,n)=1;
mm(n+1,n+1)=2;
disp(mm,"mm=");
fv=zeros(m,1);
fv(2:m,1)=di;
disp(fv,"fv=");
M=mm\fv;
disp(M,"M = ");
function si=si(i,x);
si=M(i)*((xi(i+1)-x)^3)/(6*hi(i))+M(i+1)*((x-xi(i))^3)/(6*hi(i))+((yi(i+1)-yi(i))/hi(i)-(hi(i)/6)*(M(i+1)-M(i)))*(x-xi(i))+yi(i)-M(i)*(hi(i)^2)/6;
endfunction;

function s=s(x);
while x<a>=b x=x-b+a;
end;
for i=1:n if x<xi(i+1) then s=si(i,x);
break;
end;
end;
endfunction;
d=linspace(a,b,10000);
plot(d,f);
plot(xi,yi,"ok");
plot(d,s,"r");
clear;

mat · Gość

phi^3_3
clear;
funcprot(0);
clf;
a=-1;
b=1;
n=11;
h=(b-a)/n;
xi=zeros(n+1,1);
for i=0:n
xi(i+1)=a+i*h;
end;
disp(xi,"xi=");
function f=f(x);
f=1/(1+25*(x^2));
endfunction;
yi=zeros(n+3,1);
yi(1)=0;
for i=0:n
yi(i+2)=f(xi(i+1));
end;
yi(n+3)=0;
disp(yi,"yi=");
M=zeros(n+3,n+3);
M(1,1)=1;
M(1,3)=-1;
M(1,n+1)=-1;
M(1,n+3)=1
for i=2:(n+2)
for j=1:(n+3)
if i==j then M(i,j)=4;
end;
if i+1==j|i==j+1 then M(i,j)=1;
end;
end;
end;
M(n+3,1)=1
M(n+3,2)=-2
M(n+3,3)=1
M(n+3,n+1)=-1;
M(n+3,n+2)=2
M(n+3,n+3)=-1
disp(M,"M=");
c=M\yi;
disp(c,"c = ");
function y=phi3(i,x);
xip=zeros(5);
for k=1:5
xip(k)=a+(i-3+k)*h;
end;
if x<xip(1) then y=0;
else if x<xip(2) then y=((x-xip(1))^3)/(h^3);
else if x<xip(3) then y=1+3*(x-xip(2))/h+3*((x-xip(2))^2)/(h^2)-3*((x-xip(2))^3)/(h^3);
else if x<xip(4) then y=1+3*(xip(4)-x)/h+3*((xip(4)-x)^2)/(h^2)-3*((xip(4)-x)^3)/(h^3);
else if x<xip(5) then y=((xip(5)-x)^3)/(h^3);
else y=0;
end;
end;
end;
end;
end;
endfunction;
function fs=fspline(x);
fs=0;
for i=1:(n+3)
fs=fs+c(i)*phi3(i-2,x);
end;
endfunction;
d=linspace(a,b,100);
plot(d,f);
plot(xi,yi(2:(n+2)),"ok");
plot(d,fspline,"r");
clear;

mat · Gość

DYSKRETNA TRANSFORMATA SINUSOWA

clear;
funcprot(0);
clf;
function f=f(x);
f=x^(sin(x));
endfunction;
a=0;
disp("a = "+string(a));
b=2;
disp("b = "+string(b));
c=1;
disp("c = "+string(c));
N=40;
disp("N = "+string(N));
m=3;
disp("m = "+string(m));

function f=f(x);
if x<0 then f=0; else
if x<0.5 then f=1; else
if x<=1 then f=3; else
if x<1.5 then f=-2*x+5; else
if x==1.5 then f=0; else
if x<2 then f=2; else
if x==2 then f=1; else
f=0;
end;
end;
end;
end;
end;
end;
end;
endfunction;
if f(a)==f(b) then c=f(a);
end;
function tr=tr(x);
tr=(%pi/(b-a))*(x-a);
endfunction;
function tri=tr_inv(t);
tri=((b-a)/%pi)*t+a;
endfunction;
function g=g(t);
g=f(tr_inv(t))-c;
endfunction;
xk=zeros(N);
gk=zeros(N);
Gn=zeros(N);
bn=zeros(N);
for i=0: (N-1)
xk(i+1)=(%pi*i)/N;
gk(i+1)=g(xk(i+1));
end;
for j=0:(N-1)
gks=zeros(N-1);
for i=1:(N-1)
gks(i)=gk(i+1)*sin((%pi*i*j)/N);
end;
Gn(j+1)=sum(gks)
bn(j+1)=(2/N)*Gn(j+1)
end;

function G=G(t);
bns=zeros(N-1);
for j=1:(N-1);
bns(j)=bn(j+1)*sin(j*t)
end;
G=sum(bns);
endfunction;
function F=F(x);
F=G(tr(x))+c;
endfunction;
function E=E(x);
E=abs(f(x)-F(x));
endfunction
h=b-a;
d=linspace(a,b,1000);
dl=linspace(a-2*h,b+2*h,5000);
d1=linspace(a,0.5-%eps,250);
d2=linspace(0.5,b-%eps,750);
clf(0);
plot(d1,f);
plot(d2,f);
plot({0,0.5,1,1.5,2},{1,3,3,0,1},".");
plot({0.5,1.5,2},{1,2,2},"o");
plot(d,F,"r");
clf(1);
plot(d1,E)
plot(d2,E)
plot({0.5,1.5,2},{abs(1-F(0.5)),abs(2-f(1.5)),abs(2-F(2))},"o");
plot({0.5,1.5,2},{E(0.5),E(1.5),E(2)},".");

clear;

Gość

clear;
funcprot(0);
clf;
a=-1;
disp("a = "+string(a));
b=1;
disp("b = "+string(b));
n=11;
disp("n = "+string(n));
al1=1.634;
disp("alpha_1 = "+string(al1));
bt1=-4.673;
disp("beta_1 = "+string(bt1));
h=(b-a)/n;
disp("h = "+string(h));
xi=zeros(n+1,1);
for i=0:n
xi(i+1)=a+i*h;
end;
disp(xi,"xi=");
function f=f(x);
f=1/(1+25*(x^2));
endfunction;
yi=zeros(n+3,1);
yi(1)=-(h*al1)/3;
for i=0:n
yi(i+2)=f(xi(i+1));
end;
yi(n+3)=-(h*bt1)/3;
disp(yi,"yi=");
M=zeros(n+3,n+3);
M(1,1)=1;
M(1,3)=-1;
for i=2:(n+2)
for j=1:(n+3)
if i==j then M(i,j)=4;
end;
if i+1==j|i==j+1 then M(i,j)=1;
end;
end;
end;
M(n+3,n+1)=1;
M(n+3,n+3)=-1
disp(M,"M=");
c=M\yi;
disp(c,"c = ");
function y=phi3(i,x);
xip=zeros(5);
for k=1:5
xip(k)=a+(i-3+k)*h;
end;
if x<xip(1) then
y=0;
else
if x<xip(2) then
y=((x-xip(1))^3)/(h^3);
else
if x<xip(3) then
y=1+3*(x-xip(2))/h+3*((x-xip(2))^2)/(h^2)-3*((x-xip(2))^3)/(h^3);
else
if x<xip(4) then
y=1+3*(xip(4)-x)/h+3*((xip(4)-x)^2)/(h^2)-3*((xip(4)-x)^3)/(h^3);
else
if x<xip(5) then
y=((xip(5)-x)^3)/(h^3);
else
y=0;
end;
end;
end;
end;
end;
endfunction;
function fs=fspline(x);
fs=0;
for i=1:(n+3)
fs=fs+c(i)*phi3(i-2,x);
end;
endfunction;
d=linspace(a,b,100);
plot(d,f);
plot(xi,yi(2:(n+2)),"ok");
plot(d,fspline,"r");
clear;

mat · Gość

clear;
funcprot(0);
clf;
a=-1;
disp("a = "+string(a));
b=1;
disp("b = "+string(b));
n=14;
disp("n = "+string(n));
al2=0.1;
disp("alpha_2 = "+string(al2));
bt2=0;
disp("beta_2 = "+string(bt2));
h=(b-a)/n;
disp("h = "+string(h));
xi=zeros(n+1,1);
for i=0:n
xi(i+1)=a+i*h;
end;
disp(xi,"xi=");
function f=f(x);
f=1/(1+25*(x^2));
endfunction;
yi=zeros(n+3,1);
yi(1)=((h^2)*al2)/6;
for i=0:n
yi(i+2)=f(xi(i+1));
end;
yi(n+3)=((h^2)*bt2)/6;
disp(yi,"yi=");
M=zeros(n+3,n+3);
M(1,1)=1;
M(1,2)=-2;
M(1,3)=1;
for i=2:(n+2)
for j=1:(n+3)
if i==j then M(i,j)=4;
end;
if i+1==j|i==j+1 then M(i,j)=1;
end;
end;
end;
M(n+3,n+1)=1;
M(n+3,n+2)=-2;
M(n+3,n+3)=1;
disp(M,"M=");
c=M\yi;
disp(c,"c = ");
function y=phi3(i,x);
xip=zeros(5);
for k=1:5
xip(k)=a+(i-3+k)*h;
end;
if x<xip(1) then
y=0;
else
if x<xip(2) then
y=((x-xip(1))^3)/(h^3);
else
if x<xip(3) then
y=1+3*(x-xip(2))/h+3*((x-xip(2))^2)/(h^2)-3*((x-xip(2))^3)/(h^3);
else
if x<xip(4) then
y=1+3*(xip(4)-x)/h+3*((xip(4)-x)^2)/(h^2)-3*((xip(4)-x)^3)/(h^3);
else
if x<xip(5) then
y=((xip(5)-x)^3)/(h^3);
else
y=0;
end;
end;
end;
end;
end;
endfunction;
function fs=fspline(x);
fs=0;
for i=1:(n+3)
fs=fs+c(i)*phi3(i-2,x);
end;
endfunction;
d=linspace(a,b,100);
plot(d,f);
plot(xi,yi(2:(n+2)),"ok");
plot(d,fspline,"r");
clear;

mat · Gość

clear;
funcprot(0);
clf;
function f=f(x);
f=x^(sin(x));
endfunction;
a=0;
disp("a = "+string(a));
b=2;
disp("b = "+string(b));
N=5;
disp("N = "+string(N));
m=3;
disp("m = "+string(m));
function tr=tr(x);
tr=(2*x-(a+b))/(b-a)
endfunction;
function tri=tr_inv(t);
tri=((b-a)*t+(a+b))/2;
endfunction;
function g=g(t);
g=f(tr_inv(t))
endfunction;
c=zeros(m+1);
l=nearfloat("pred",1);
tc=linspace(-l,l,N);
lct=zeros(N);
mct=zeros(N);
for j=0:m
for i=1:N
lct(i)=g(tc(i))*legendre(j,0,tc(i));
mct(i)=legendre(j,0,tc(i))*legendre(j,0,tc(i));
end;
lc=inttrap(tc,lct);
mc=inttrap(tc,mct);
c(j+1)=lc/mc;
end;
disp(c,"c=");
function G=G(t);
cp=zeros(m+1);
if t<1>=1 then t=nearfloat("pred",1);
end
for j=0:m
cp(j+1)=c(j+1)*legendre(j,0,t);
end;
G=sum(cp);
endfunction;
function F=F(x);
F=G(tr(x));
endfunction;
d=linspace(a,b,1000);
plot(d,f);
plot(d,F,"r");
clear;