Harjoitus 8 >> load data-har8 >> >> >> da0 da0 = 0 0.4444 0.8889 1.3333 1.7778 2.2222 2.6667 3.1111 3.5556 4.0000 0 0.2066 -0.1324 -0.4176 0.5898 0.6377 -1.8433 -0.3695 4.7206 -2.1273 >> da1 da1 = 1.4950 3.0525 3.5100 3.5647 3.9181 5.1705 5.6749 6.0375 6.1241 6.3616 -3.4546 0.3094 -4.8565 -4.4540 -0.0274 6.1305 -3.0787 1.4855 -0.9023 4.4694 >> da2 da2 = -12.0000 -9.3333 -6.6667 -4.0000 -1.3333 1.3333 4.0000 6.6667 9.3333 12.0000 0.0000 0.0000 0.0000 0.0000 0.0912 0.9088 1.0000 1.0000 1.0000 1.0000 Tehtävä 1 da0: spline: derivaatan arvot alku- ja loppupisteissä -2 ja 0 >> xx=0:.01:4; >> y=[-2,da0(2,:),0] y = -2.0000 0 0.2066 -0.1324 -0.4176 0.5898 0.6377 -1.8433 -0.3695 4.7206 -2.1273 0 >> yy=spline(da0(1,:),y,xx); >> plot(da0(1,:),da0(2,:),'ro',xx,yy) >> title('H8T1, spline: derivaatan arvot -2 ja 0') >> legend('da0','splini') csape: not-a-knot >> pp=csape(da0(1,:),da0(2,:),'not-a-knot'); >> pp pp = form: 'pp' breaks: [0 0.4444 0.8889 1.3333 1.7778 2.2222 2.6667 3.1111 3.5556 4] coefs: [9x4 double] pieces: 9 order: 4 dim: 1 >> [solmut,kertoimet,l,k,d]=unmkpp(pp) solmut = 0 0.4444 0.8889 1.3333 1.7778 2.2222 2.6667 3.1111 3.5556 4.0000 kertoimet = 0.3986 -1.9124 1.2360 0 0.3986 -1.3809 -0.2277 0.2066 4.8342 -0.8495 -1.2190 -0.1324 -5.6259 5.5961 0.8906 -0.4176 -7.9818 -1.9050 2.5311 0.5898 19.6760 -12.5475 -3.8922 0.6377 3.1317 13.6873 -3.3856 -1.8433 -36.0607 17.8629 10.6366 -0.3695 -36.0607 -30.2180 5.1455 4.7206 l = 9 k = 4 d = 1 >> xx=0:.01:4; >> yy=ppval(pp,xx); >> plot(da0(1,:),da0(2,:),'ro',xx,yy) >> legend('da0','splini') >> title('H8T1, csape: not-a-knot') da1: csape >> pp=csape(da1(1,:),da1(2,:)); >> pp pp = form: 'pp' breaks: [1.4950 3.0525 3.5100 3.5647 3.9181 5.1705 5.6749 6.0375 6.1241 6.3616] coefs: [9x4 double] pieces: 9 order: 4 dim: 1 >> [solmut,kertoimet,l,k,d]=unmkpp(pp) solmut = 1.4950 3.0525 3.5100 3.5647 3.9181 5.1705 5.6749 6.0375 6.1241 6.3616 kertoimet = 1.0e+003 * 0.0208 -0.0812 0.0785 -0.0035 0.0225 0.0159 -0.0233 0.0003 -0.2034 0.0468 0.0054 -0.0049 -0.0072 0.0134 0.0087 -0.0045 -0.0113 0.0057 0.0155 -0.0000 0.0926 -0.0367 -0.0233 0.0061 -0.2686 0.1035 0.0104 -0.0031 1.2396 -0.1887 -0.0205 0.0015 0.2876 0.1335 -0.0253 -0.0009 l = 9 k = 4 d = 1 >> xx=linspace(da1(1,1),da1(1,10),400); >> yy=ppval(pp,xx); >> plot(da1(1,:),da1(2,:),'ro',xx,yy) >> legend('da1','splini') >> title('H8T1, da1, csape') csape: clamped -1 ja 0 >> pp=csape(da1(1,:),[-1,da1(2,:),0],'clamped'); >> pp pp = form: 'pp' breaks: [1.4950 3.0525 3.5100 3.5647 3.9181 5.1705 5.6749 6.0375 6.1241 6.3616] coefs: [9x4 double] pieces: 9 order: 4 dim: 1 >> [solmut,kertoimet,l,k,d]=unmkpp(pp) solmut = 1.4950 3.0525 3.5100 3.5647 3.9181 5.1705 5.6749 6.0375 6.1241 6.3616 kertoimet = 1.0e+003 * -0.0082 0.0149 -0.0010 -0.0035 0.0639 -0.0233 -0.0140 0.0003 -0.3190 0.0644 0.0048 -0.0049 -0.0054 0.0120 0.0090 -0.0045 -0.0117 0.0063 0.0154 -0.0000 0.0972 -0.0377 -0.0240 0.0061 -0.2983 0.1093 0.0122 -0.0031 2.3038 -0.2152 -0.0262 0.0015 -1.0083 0.3836 -0.0116 -0.0009 l = 9 k = 4 d = 1 >> xx=linspace(da1(1,1),da1(1,10),400); >> yy=ppval(pp,xx); >> plot(da1(1,:),da1(2,:),'ro',xx,yy) >> title('H8T1, da1, csape: clamped -1 ja 0') >> legend('da1','splini') Tehtävä 2: da1: csape: clamped -1 ja 0 (a) >> fnplt(pp) >> hold on >> plot(da1(1,:),da1(2,:),'ro') >> legend('splini','da1') >> title('H8T2(a), da1, csape: clamped') (b) >> I=fnint(pp) I = form: 'pp' breaks: [1.4950 3.0525 3.5100 3.5647 3.9181 5.1705 5.6749 6.0375 6.1241 6.3616] coefs: [9x5 double] pieces: 9 order: 5 dim: 1 >> [solmut,kertoimet,l,k,d]=unmkpp(I) solmut = 1.4950 3.0525 3.5100 3.5647 3.9181 5.1705 5.6749 6.0375 6.1241 6.3616 kertoimet = -2.0446 4.9772 -0.5000 -3.4546 0 15.9639 -7.7606 -7.0027 0.3094 0.1797 -79.7603 21.4525 2.3931 -4.8565 -1.1881 -1.3400 3.9899 4.4820 -4.4540 -1.4440 -2.9284 2.0955 7.7081 -0.0274 -2.3031 24.2952 -12.5745 -11.9775 6.1305 6.6648 -74.5858 36.4448 6.0833 -3.0787 6.6686 575.9444 -71.7269 -13.1053 1.4855 6.8002 -252.0818 127.8743 -5.8083 -0.9023 6.8163 l = 9 k = 5 d = 1 >> xx=linspace(da1(1,1),da1(1,10),400); >> yy=ppval(I,xx); >> plot(da1(1,:),da1(2,:),'ro',xx,yy) >> legend('da1','fnint') >> title('H8T2(b), da1, csape: clamped, fnint') (c) >> D=fnder(pp) D = form: 'pp' breaks: [1.4950 3.0525 3.5100 3.5647 3.9181 5.1705 5.6749 6.0375 6.1241 6.3616] coefs: [9x3 double] pieces: 9 order: 3 dim: 1 >> [solmut,kertoimet,l,k,d]=unmkpp(D) solmut = 1.4950 3.0525 3.5100 3.5647 3.9181 5.1705 5.6749 6.0375 6.1241 6.3616 kertoimet = 1.0e+003 * -0.0245 0.0299 -0.0010 0.1916 -0.0466 -0.0140 -0.9571 0.1287 0.0048 -0.0161 0.0239 0.0090 -0.0351 0.0126 0.0154 0.2915 -0.0754 -0.0240 -0.8950 0.2187 0.0122 6.9113 -0.4304 -0.0262 -3.0250 0.7672 -0.0116 l = 9 k = 3 d = 1 >> xx=linspace(da1(1,1),da1(1,10),400); >> yy=ppval(D,xx); >> plot(da1(1,:),da1(2,:),'ro',xx,yy) >> title('H8T2(c), da1, csape: clamped, fnder') >> legend('da1','fnder') (d) >> A=fnzeros(pp) A = 2.1502 3.0738 3.9199 5.3931 5.8169 6.0851 6.2005 2.1502 3.0738 3.9199 5.3931 5.8169 6.0851 6.2005 >> xx=linspace(da1(1,1),da1(1,10),400); >> yy=ppval(pp,xx); >> plot(da1(1,:),da1(2,:),'ro',xx,yy) >> hold on >> plot(A(1,:),zeros(length(A)),'g*') >> legend('da1','splini','nollakohdat') >> title('H8T2(d), da1, csape: clamped, nollakohdat') >> D=fnder(pp) D = form: 'pp' breaks: [1.4950 3.0525 3.5100 3.5647 3.9181 5.1705 5.6749 6.0375 6.1241 6.3616] coefs: [9x3 double] pieces: 9 order: 3 dim: 1 >> B=fnzeros(D) B = 1.5294 2.6777 3.4704 4.7831 5.6144 5.9659 6.1403 1.5294 2.6777 3.4704 4.7831 5.6144 5.9659 6.1403 >> xx=linspace(da1(1,1),da1(1,10),400); >> yy=ppval(pp,xx); >> yyy=ppval(D,xx); >> plot(da1(1,:),da1(2,:),'ro',xx,yy) >> hold on >> plot(xx,yyy,'k') >> plot(B(1,:),zeros(length(B)),'g*') >> legend('da1','splini','fnder','derivaatan nollakohdat') >> title('H8T2(d), da1, csape: clamped, derivaatan nollakohdat') (e) >> min=fnmin(pp) min = -4.9490 >> xx=linspace(da1(1,1),da1(1,10),400); >> yy=ppval(pp,xx); >> plot(da1(1,:),da1(2,:),'ro',xx,yy) >> hold on >> plot(xx,min,'k') >> max=-fnmin(fncmb(pp,-1)) max = 10.4302 >> plot(xx,max,'k') >> legend('da1','splini','min=-4.9490','max=10.4302') >> title('H8T2(e), da1, csape: clamped, min ja max') Tehtävä 3: >> pp=csape(da2(1,:),da2(2,:),'periodic'); >> [solmut,kertoimet,l,k,d]=unmkpp(pp) solmut = -12.0000 -9.3333 -6.6667 -4.0000 -1.3333 1.3333 4.0000 6.6667 9.3333 12.0000 kertoimet = -0.0002 0.0000 0.0011 0.0000 0.0008 -0.0012 -0.0022 0.0000 -0.0030 0.0050 0.0078 0.0000 0.0159 -0.0187 -0.0289 0.0000 -0.0271 0.1084 0.2102 0.0912 0.0159 -0.1084 0.2102 0.9088 -0.0030 0.0187 -0.0289 1.0000 0.0008 -0.0050 0.0078 1.0000 -0.0002 0.0012 -0.0022 1.0000 l = 9 k = 4 d = 1 >> xx=linspace(-12,12,1000); >> yy=ppval(pp,xx); >> plot(da2(1,:),da2(2,:),'ro',xx,yy) >> hold on >> PP=pchip(da2(1,:),da2(2,:)) PP = form: 'pp' breaks: [-12 -9.3333 -6.6667 -4 -1.3333 1.3333 4 6.6667 9.3333 12] coefs: [9x4 double] pieces: 9 order: 4 dim: 1 >> [Solmut,Kertoimet,L,K,dim]=unmkpp(PP) Solmut = -12.0000 -9.3333 -6.6667 -4.0000 -1.3333 1.3333 4.0000 6.6667 9.3333 12.0000 Kertoimet = -0.0000 0.0000 0 0.0000 -0.0000 0.0000 0.0000 0.0000 -0.0000 0.0000 0.0000 0.0000 -0.0010 0.0154 0.0000 0.0000 -0.0689 0.2757 0.0615 0.0912 -0.0010 -0.0077 0.0615 0.9088 -0.0000 -0.0000 0.0000 1.0000 -0.0000 -0.0000 0.0000 1.0000 0 0 0 1.0000 L = 9 K = 4 dim = 1 >> yyy=ppval(PP,xx); >> plot(xx,yyy,'k') >> legend('da2','splini','pchip') >> title('H8T3, da2, csape: periodic, pchip-splini') Tehtävä 4: Ellipsin yhtälö on x^2/a^2+y^2/b^2=1. Nyt x^2+2y^2=3, eli x^2/3+2y^2/3=1, joten x^2/(sqrt(3))^2+y^2/((sqrt(3/2))^2=1. Siis a=sqrt(3) ja b=sqrt(3/2). Ellipsin parametriesitys on x=a*cos(t) y=b*sin(t), missä t \in [0,2*pi]. Siis nyt ellipsin parametriesitys on x=sqrt(3)*cos(t) y=sqrt(3/2)*sin(t), missä t \in [0,2*pi]. >> t=linspace(0,2*pi,100); >> x=sqrt(3)*cos(t); >> y=sqrt(3/2)*sin(t); >> plot(x,y) >> legend('ellipsi x^{2}+2y^{2}=3') >> title('H8T4 ellipsi') >> T=linspace(0,2*pi,10); >> X=sqrt(3)*cos(T); >> pp=csape(T,X,'periodic') pp = form: 'pp' breaks: [0 0.6981 1.3963 2.0944 2.7925 3.4907 4.1888 4.8869 5.5851 6.2832] coefs: [9x4 double] pieces: 9 order: 4 dim: 1 >> [solmut,kertoimet,l,k,d]=unmkpp(pp) solmut = 0 0.6981 1.3963 2.0944 2.7925 3.4907 4.1888 4.8869 5.5851 6.2832 kertoimet = 0.1007 -0.9017 0.0000 1.7321 0.2551 -0.6908 -1.1118 1.3268 0.2900 -0.1566 -1.7034 0.3008 0.1893 0.4509 -1.4979 -0.8660 0.0000 0.8474 -0.5916 -1.6276 -0.1893 0.8474 0.5916 -1.6276 -0.2900 0.4509 1.4979 -0.8660 -0.2551 -0.1566 1.7034 0.3008 -0.1007 -0.6908 1.1118 1.3268 l = 9 k = 4 d = 1 >> Y=sqrt(3/2)*sin(T); >> pp2=csape(T,Y,'periodic') pp2 = form: 'pp' breaks: [0 0.6981 1.3963 2.0944 2.7925 3.4907 4.1888 4.8869 5.5851 6.2832] coefs: [9x4 double] pieces: 9 order: 4 dim: 1 >> [Solmut,Kertoimet,L,K,D]=unmkpp(pp2) Solmut = 0 0.6981 1.3963 2.0944 2.7925 3.4907 4.1888 4.8869 5.5851 6.2832 Kertoimet = -0.1957 0 1.2230 0 -0.1041 -0.4099 0.9369 0.7873 0.0362 -0.6279 0.2124 1.2061 0.1595 -0.5522 -0.6115 1.0607 0.2083 -0.2181 -1.1493 0.4189 0.1595 0.2181 -1.1493 -0.4189 0.0362 0.5522 -0.6115 -1.0607 -0.1041 0.6279 0.2124 -1.2061 -0.1957 0.4099 0.9369 -0.7873 L = 9 K = 4 D = 1 >> t=linspace(0,2*pi,100); >> plot(ppval(pp,t),ppval(pp2,t)) >> hold on >> plot(X,Y,'ro') >> legend('splini','(X,Y)') >> title('H8T4, interpoloiva splini') Tehtävä 5: Tehtävä 5 p = 0 -1.0000 1.0000 1.5000 2.0000 1.5000 1.0000 -2.0000 2.0000 3.0000 0 1.0000 2.0000 2.5000 2.0000 0 0.5000 1.0000 2.0000 1.0000 >> A=p' A = 0 0 -1.0000 1.0000 1.0000 2.0000 1.5000 2.5000 2.0000 2.0000 1.5000 0 1.0000 0.5000 -2.0000 1.0000 2.0000 2.0000 3.0000 1.0000 >> bezier(A) ans = 175.0175 >> legend('Bezier-käyrä') >> title('H8T5, Bezier-käyrä')