push
authorNils Forssén <forssennils@gmail.com>
Wed, 12 Oct 2022 14:05:17 +0000 (16:05 +0200)
committerNils Forssén <forssennils@gmail.com>
Wed, 12 Oct 2022 14:05:17 +0000 (16:05 +0200)
15 files changed:
LorenzAttractor.m [new file with mode: 0644]
Project/N_20.fig [new file with mode: 0644]
Project/N_20.svg [new file with mode: 0644]
Project/N_5.fig [new file with mode: 0644]
Project/N_5.svg [new file with mode: 0644]
Project/lufact.m [new file with mode: 0644]
Project/script.m [new file with mode: 0644]
Project/thomas.asv [new file with mode: 0644]
Project/thomas.m [new file with mode: 0644]
RungeKutta4.m [new file with mode: 0644]
SimpsonsRule.m [new file with mode: 0644]
adaptativeSimpson.m [new file with mode: 0644]
forwardEuler.m [new file with mode: 0644]
lab3_k.m [new file with mode: 0644]
plot_composite_quad.m [new file with mode: 0644]

diff --git a/LorenzAttractor.m b/LorenzAttractor.m
new file mode 100644 (file)
index 0000000..bc7ee53
--- /dev/null
@@ -0,0 +1,63 @@
+function LorenzAttractor(T)
+%%
+%  Nonlinear ODE example from nonlinear dynamics to solve the
+%  Lorenz system using a Runge-Kutta 4 time integrator
+%
+%      d(y1)/dt = 10*(y2-y1)
+%      d(y2)/dt = 28*y1 - y2 - y1*y3
+%      d(y3)/dt = y1*y2 - 8/3*By3
+%
+%  INPUT:
+%
+%     T - final time
+%
+%  OUTPUT:
+%
+%     NONE - this makes a movie
+%
+%%
+%  adhoc tactic to set the time step
+   num_steps = 2500;
+   dt = T/num_steps;
+%%
+% set the initial condition and allocate memory for the solution vector
+   y   = zeros(3,num_steps);
+   ini = 5*ones(3,1);
+%%
+% save the initial condition and enter the Runge-Kutta loop to solve the ODE
+   y(:,1) = ini;
+   for i = 2:num_steps
+      k1 = RHS(y(:,i-1));
+      k2 = RHS(y(:,i-1) + 0.5*dt.*k1);
+      k3 = RHS(y(:,i-1) + 0.5*dt.*k2);
+      k4 = RHS(y(:,i-1) + dt.*k3);
+      y(:,i) = y(:,i-1) + (dt/6).*(k1 + 2.*k2 + 2.*k3 + k4);
+   end
+%%
+% visualize the solution progessively to make a movie
+   clf;
+   hold on;
+   xlim([-20 20]);
+   ylim([-30 30]);
+   zlim([0 50]);
+   grid on;
+   for n = 1:(length(y)/7)
+      plot3(y(1,1:7*n),y(2,1:7*n),y(3,1:7*n),'b');
+      drawnow
+      view(-37.5-n, 24)
+      pause(0.02)
+   end
+end
+
+%%%
+%  auxiliary function to get the right hand side
+%%%
+ function dy = RHS(y)
+%%
+% evaluate the right hand side of the nonlinear ODE for the Lorenz
+% attractor. The input y is a vector with three components.
+   dy    = zeros(3,1);
+   dy(1) = 
+   dy(2) = 
+   dy(3) = 
+ end
\ No newline at end of file
diff --git a/Project/N_20.fig b/Project/N_20.fig
new file mode 100644 (file)
index 0000000..5eb85ad
Binary files /dev/null and b/Project/N_20.fig differ
diff --git a/Project/N_20.svg b/Project/N_20.svg
new file mode 100644 (file)
index 0000000..64b7199
--- /dev/null
@@ -0,0 +1,180 @@
+<?xml version="1.0"?>
+<!DOCTYPE svg PUBLIC '-//W3C//DTD SVG 1.0//EN'
+          'http://www.w3.org/TR/2001/REC-SVG-20010904/DTD/svg10.dtd'>
+<svg xmlns:xlink="http://www.w3.org/1999/xlink" style="fill-opacity:1; color-rendering:auto; color-interpolation:auto; text-rendering:auto; stroke:black; stroke-linecap:square; stroke-miterlimit:10; shape-rendering:auto; stroke-opacity:1; fill:black; stroke-dasharray:none; font-weight:normal; stroke-width:1; font-family:'Dialog'; font-style:normal; stroke-linejoin:miter; font-size:12px; stroke-dashoffset:0; image-rendering:auto;" width="700" height="525" xmlns="http://www.w3.org/2000/svg"
+><!--Generated by the Batik Graphics2D SVG Generator--><defs id="genericDefs"
+  /><g
+  ><defs id="defs1"
+    ><clipPath clipPathUnits="userSpaceOnUse" id="clipPath1"
+      ><path d="M0 0 L700 0 L700 525 L0 525 L0 0 Z"
+      /></clipPath
+      ><font horiz-adv-x="75.0" id="font1"
+      ><font-face ascent="100.53711" descent="21.972656" units-per-em="100" style="font-style:normal; font-family:Dialog; font-weight:normal;"
+        /><missing-glyph horiz-adv-x="75.0" d="M12.5 0 L12.5 62.5 L62.5 62.5 L62.5 0 L12.5 0 ZM14.0625 1.5625 L60.9375 1.5625 L60.9375 60.9375 L14.0625 60.9375 L14.0625 1.5625 Z"
+        /><glyph unicode="1" horiz-adv-x="55.615234" d="M37.25 0 L28.4688 0 L28.4688 56 Q25.2969 52.9844 20.1406 49.9531 Q14.9844 46.9219 10.8906 45.4062 L10.8906 53.9062 Q18.2656 57.375 23.7812 62.3047 Q29.2969 67.2344 31.5938 71.875 L37.25 71.875 L37.25 0 Z"
+        /><glyph unicode="." horiz-adv-x="27.783203" d="M9.0781 0 L9.0781 10.0156 L19.0938 10.0156 L19.0938 0 L9.0781 0 Z"
+        /><glyph unicode="2" horiz-adv-x="55.615234" d="M50.3438 8.4531 L50.3438 0 L3.0312 0 Q2.9375 3.1719 4.0469 6.1094 Q5.8594 10.9375 9.8359 15.625 Q13.8125 20.3125 21.3438 26.4688 Q33.0156 36.0312 37.1172 41.625 Q41.2188 47.2188 41.2188 52.2031 Q41.2188 57.4219 37.4766 61.0078 Q33.7344 64.5938 27.7344 64.5938 Q21.3906 64.5938 17.5781 60.7891 Q13.7656 56.9844 13.7188 50.25 L4.6875 51.1719 Q5.6094 61.2812 11.6641 66.5781 Q17.7188 71.875 27.9375 71.875 Q38.2344 71.875 44.2422 66.1641 Q50.25 60.4531 50.25 52 Q50.25 47.7031 48.4922 43.5547 Q46.7344 39.4062 42.6562 34.8125 Q38.5781 30.2188 29.1094 22.2188 Q21.1875 15.5781 18.9453 13.2109 Q16.7031 10.8438 15.2344 8.4531 L50.3438 8.4531 Z"
+        /><glyph unicode="3" horiz-adv-x="55.615234" d="M4.2031 18.8906 L12.9844 20.0625 Q14.5 12.5938 18.1406 9.2969 Q21.7812 6 27 6 Q33.2031 6 37.4766 10.2969 Q41.75 14.5938 41.75 20.9531 Q41.75 27 37.7969 30.9297 Q33.8438 34.8594 27.7344 34.8594 Q25.25 34.8594 21.5312 33.8906 L22.5156 41.6094 Q23.3906 41.5 23.9219 41.5 Q29.5469 41.5 34.0391 44.4297 Q38.5312 47.3594 38.5312 53.4688 Q38.5312 58.2969 35.2578 61.4766 Q31.9844 64.6562 26.8125 64.6562 Q21.6875 64.6562 18.2656 61.4297 Q14.8438 58.2031 13.875 51.7656 L5.0781 53.3281 Q6.6875 62.1562 12.3984 67.0156 Q18.1094 71.875 26.6094 71.875 Q32.4688 71.875 37.3984 69.3594 Q42.3281 66.8438 44.9453 62.5 Q47.5625 58.1562 47.5625 53.2656 Q47.5625 48.6406 45.0703 44.8281 Q42.5781 41.0156 37.7031 38.7656 Q44.0469 37.3125 47.5625 32.6953 Q51.0781 28.0781 51.0781 21.1406 Q51.0781 11.7656 44.2422 5.25 Q37.4062 -1.2656 26.9531 -1.2656 Q17.5312 -1.2656 11.3047 4.3516 Q5.0781 9.9688 4.2031 18.8906 Z"
+        /><glyph unicode="4" horiz-adv-x="55.615234" d="M32.3281 0 L32.3281 17.1406 L1.2656 17.1406 L1.2656 25.2031 L33.9375 71.5781 L41.1094 71.5781 L41.1094 25.2031 L50.7812 25.2031 L50.7812 17.1406 L41.1094 17.1406 L41.1094 0 L32.3281 0 ZM32.3281 25.2031 L32.3281 57.4688 L9.9062 25.2031 L32.3281 25.2031 Z"
+        /><glyph unicode="5" horiz-adv-x="55.615234" d="M4.1562 18.75 L13.375 19.5312 Q14.4062 12.7969 18.1406 9.3984 Q21.875 6 27.1562 6 Q33.5 6 37.8906 10.7891 Q42.2812 15.5781 42.2812 23.4844 Q42.2812 31 38.0625 35.3516 Q33.8438 39.7031 27 39.7031 Q22.75 39.7031 19.3359 37.7734 Q15.9219 35.8438 13.9688 32.7656 L5.7188 33.8438 L12.6406 70.6094 L48.25 70.6094 L48.25 62.2031 L19.6719 62.2031 L15.8281 42.9688 Q22.2656 47.4688 29.3438 47.4688 Q38.7188 47.4688 45.1641 40.9688 Q51.6094 34.4688 51.6094 24.2656 Q51.6094 14.5469 45.9531 7.4688 Q39.0625 -1.2188 27.1562 -1.2188 Q17.3906 -1.2188 11.2109 4.25 Q5.0312 9.7188 4.1562 18.75 Z"
+        /><glyph unicode="6" horiz-adv-x="55.615234" d="M49.75 54.0469 L41.0156 53.375 Q39.8438 58.5469 37.7031 60.8906 Q34.125 64.6562 28.9062 64.6562 Q24.7031 64.6562 21.5312 62.3125 Q17.3906 59.2812 14.9922 53.4688 Q12.5938 47.6562 12.5 36.9219 Q15.6719 41.75 20.2656 44.0938 Q24.8594 46.4375 29.8906 46.4375 Q38.6719 46.4375 44.8516 39.9688 Q51.0312 33.5 51.0312 23.25 Q51.0312 16.5 48.125 10.7188 Q45.2188 4.9375 40.1406 1.8594 Q35.0625 -1.2188 28.6094 -1.2188 Q17.625 -1.2188 10.6953 6.8594 Q3.7656 14.9375 3.7656 33.5 Q3.7656 54.25 11.4219 63.6719 Q18.1094 71.875 29.4375 71.875 Q37.8906 71.875 43.2891 67.1406 Q48.6875 62.4062 49.75 54.0469 ZM13.875 23.1875 Q13.875 18.6562 15.7969 14.5078 Q17.7188 10.3594 21.1875 8.1797 Q24.6562 6 28.4688 6 Q34.0312 6 38.0391 10.4922 Q42.0469 14.9844 42.0469 22.7031 Q42.0469 30.125 38.0859 34.3984 Q34.125 38.6719 28.125 38.6719 Q22.1719 38.6719 18.0234 34.3984 Q13.875 30.125 13.875 23.1875 Z"
+        /><glyph unicode="7" horiz-adv-x="55.615234" d="M4.7344 62.2031 L4.7344 70.6562 L51.0781 70.6562 L51.0781 63.8125 Q44.2344 56.5469 37.5234 44.4844 Q30.8125 32.4219 27.1562 19.6719 Q24.5156 10.6875 23.7812 0 L14.75 0 Q14.8906 8.4531 18.0625 20.4141 Q21.2344 32.375 27.1719 43.4844 Q33.1094 54.5938 39.7969 62.2031 L4.7344 62.2031 Z"
+        /><glyph unicode="8" horiz-adv-x="55.615234" d="M17.6719 38.8125 Q12.2031 40.8281 9.5703 44.5391 Q6.9375 48.25 6.9375 53.4219 Q6.9375 61.2344 12.5547 66.5547 Q18.1719 71.875 27.4844 71.875 Q36.8594 71.875 42.5781 66.4297 Q48.2969 60.9844 48.2969 53.1719 Q48.2969 48.1875 45.6797 44.5078 Q43.0625 40.8281 37.75 38.8125 Q44.3438 36.6719 47.7812 31.8828 Q51.2188 27.0938 51.2188 20.4531 Q51.2188 11.2812 44.7266 5.0312 Q38.2344 -1.2188 27.6406 -1.2188 Q17.0469 -1.2188 10.5469 5.0547 Q4.0469 11.3281 4.0469 20.7031 Q4.0469 27.6875 7.5938 32.3984 Q11.1406 37.1094 17.6719 38.8125 ZM15.9219 53.7188 Q15.9219 48.6406 19.1953 45.4141 Q22.4688 42.1875 27.6875 42.1875 Q32.7656 42.1875 36.0156 45.3828 Q39.2656 48.5781 39.2656 53.2188 Q39.2656 58.0625 35.9141 61.3594 Q32.5625 64.6562 27.5938 64.6562 Q22.5625 64.6562 19.2422 61.4297 Q15.9219 58.2031 15.9219 53.7188 ZM13.0938 20.6562 Q13.0938 16.8906 14.875 13.375 Q16.6562 9.8594 20.1719 7.9297 Q23.6875 6 27.7344 6 Q34.0312 6 38.1328 10.0547 Q42.2344 14.1094 42.2344 20.3594 Q42.2344 26.7031 38.0156 30.8594 Q33.7969 35.0156 27.4375 35.0156 Q21.2344 35.0156 17.1641 30.9141 Q13.0938 26.8125 13.0938 20.6562 Z"
+        /><glyph unicode="9" horiz-adv-x="55.615234" d="M5.4688 16.5469 L13.9219 17.3281 Q14.9844 11.375 18.0156 8.6875 Q21.0469 6 25.7812 6 Q29.8281 6 32.8828 7.8594 Q35.9375 9.7188 37.8906 12.8203 Q39.8438 15.9219 41.1641 21.1953 Q42.4844 26.4688 42.4844 31.9375 Q42.4844 32.5156 42.4375 33.6875 Q39.7969 29.5 35.2344 26.8828 Q30.6719 24.2656 25.3438 24.2656 Q16.4531 24.2656 10.3047 30.7109 Q4.1562 37.1562 4.1562 47.7031 Q4.1562 58.5938 10.5781 65.2344 Q17 71.875 26.6562 71.875 Q33.6406 71.875 39.4297 68.1172 Q45.2188 64.3594 48.2188 57.3984 Q51.2188 50.4375 51.2188 37.25 Q51.2188 23.5312 48.2422 15.4062 Q45.2656 7.2812 39.3828 3.0312 Q33.5 -1.2188 25.5938 -1.2188 Q17.1875 -1.2188 11.8672 3.4453 Q6.5469 8.1094 5.4688 16.5469 ZM41.4531 48.1406 Q41.4531 55.7188 37.4297 60.1562 Q33.4062 64.5938 27.7344 64.5938 Q21.875 64.5938 17.5312 59.8125 Q13.1875 55.0312 13.1875 47.4062 Q13.1875 40.5781 17.3125 36.3047 Q21.4375 32.0312 27.4844 32.0312 Q33.5938 32.0312 37.5234 36.3047 Q41.4531 40.5781 41.4531 48.1406 Z"
+        /><glyph unicode="0" horiz-adv-x="55.615234" d="M4.1562 35.2969 Q4.1562 48 6.7656 55.7422 Q9.375 63.4844 14.5234 67.6797 Q19.6719 71.875 27.4844 71.875 Q33.25 71.875 37.5938 69.5547 Q41.9375 67.2344 44.7734 62.8672 Q47.6094 58.5 49.2188 52.2266 Q50.8281 45.9531 50.8281 35.2969 Q50.8281 22.7031 48.2422 14.9688 Q45.6562 7.2344 40.5078 3.0078 Q35.3594 -1.2188 27.4844 -1.2188 Q17.1406 -1.2188 11.2344 6.2031 Q4.1562 15.1406 4.1562 35.2969 ZM13.1875 35.2969 Q13.1875 17.6719 17.3125 11.8359 Q21.4375 6 27.4844 6 Q33.5469 6 37.6719 11.8594 Q41.7969 17.7188 41.7969 35.2969 Q41.7969 52.9844 37.6719 58.7891 Q33.5469 64.5938 27.3906 64.5938 Q21.3438 64.5938 17.7188 59.4688 Q13.1875 52.9375 13.1875 35.2969 Z"
+      /></font
+    ></defs
+    ><g style="fill:white; stroke:white;"
+    ><rect x="0" y="0" width="700" style="clip-path:url(#clipPath1); stroke:none;" height="525"
+    /></g
+    ><g style="fill:white; text-rendering:optimizeSpeed; color-rendering:optimizeSpeed; image-rendering:optimizeSpeed; shape-rendering:crispEdges; stroke:white; color-interpolation:sRGB;"
+    ><rect x="0" width="700" height="525" y="0" style="stroke:none;"
+      /><path style="stroke:none;" d="M91 467 L634 467 L634 39 L91 39 Z"
+    /></g
+    ><g style="fill:rgb(38,38,38); text-rendering:geometricPrecision; image-rendering:optimizeQuality; color-rendering:optimizeQuality; stroke-linejoin:round; stroke:rgb(38,38,38); color-interpolation:linearRGB; stroke-width:0.8333;"
+    ><line y2="467" style="fill:none;" x1="91" x2="634" y1="467"
+      /><line y2="39" style="fill:none;" x1="91" x2="634" y1="39"
+      /><line y2="461.57" style="fill:none;" x1="91" x2="91" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="145.3" x2="145.3" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="199.6" x2="199.6" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="253.9" x2="253.9" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="308.2" x2="308.2" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="362.5" x2="362.5" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="416.8" x2="416.8" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="471.1" x2="471.1" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="525.4" x2="525.4" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="579.7" x2="579.7" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="634" x2="634" y1="467"
+      /><line y2="44.43" style="fill:none;" x1="91" x2="91" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="145.3" x2="145.3" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="199.6" x2="199.6" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="253.9" x2="253.9" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="308.2" x2="308.2" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="362.5" x2="362.5" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="416.8" x2="416.8" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="471.1" x2="471.1" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="525.4" x2="525.4" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="579.7" x2="579.7" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="634" x2="634" y1="39"
+    /></g
+    ><g transform="translate(91,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-5" xml:space="preserve" y="17" style="stroke:none;"
+      >1</text
+    ></g
+    ><g transform="translate(145.3,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.1</text
+    ></g
+    ><g transform="translate(199.6,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.2</text
+    ></g
+    ><g transform="translate(253.9,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.3</text
+    ></g
+    ><g transform="translate(308.2,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.4</text
+    ></g
+    ><g transform="translate(362.5,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.5</text
+    ></g
+    ><g transform="translate(416.8,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.6</text
+    ></g
+    ><g transform="translate(471.1,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.7</text
+    ></g
+    ><g transform="translate(525.4,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.8</text
+    ></g
+    ><g transform="translate(579.7,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.9</text
+    ></g
+    ><g transform="translate(634,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-5" xml:space="preserve" y="17" style="stroke:none;"
+      >2</text
+    ></g
+    ><g style="fill:rgb(38,38,38); text-rendering:geometricPrecision; image-rendering:optimizeQuality; color-rendering:optimizeQuality; stroke-linejoin:round; stroke:rgb(38,38,38); color-interpolation:linearRGB; stroke-width:0.8333;"
+    ><line y2="39" style="fill:none;" x1="91" x2="91" y1="467"
+      /><line y2="39" style="fill:none;" x1="634" x2="634" y1="467"
+      /><line y2="467" style="fill:none;" x1="91" x2="96.43" y1="467"
+      /><line y2="424.2" style="fill:none;" x1="91" x2="96.43" y1="424.2"
+      /><line y2="381.4" style="fill:none;" x1="91" x2="96.43" y1="381.4"
+      /><line y2="338.6" style="fill:none;" x1="91" x2="96.43" y1="338.6"
+      /><line y2="295.8" style="fill:none;" x1="91" x2="96.43" y1="295.8"
+      /><line y2="253" style="fill:none;" x1="91" x2="96.43" y1="253"
+      /><line y2="210.2" style="fill:none;" x1="91" x2="96.43" y1="210.2"
+      /><line y2="167.4" style="fill:none;" x1="91" x2="96.43" y1="167.4"
+      /><line y2="124.5999" style="fill:none;" x1="91" x2="96.43" y1="124.5999"
+      /><line y2="81.8" style="fill:none;" x1="91" x2="96.43" y1="81.8"
+      /><line y2="39" style="fill:none;" x1="91" x2="96.43" y1="39"
+      /><line y2="467" style="fill:none;" x1="634" x2="628.57" y1="467"
+      /><line y2="424.2" style="fill:none;" x1="634" x2="628.57" y1="424.2"
+      /><line y2="381.4" style="fill:none;" x1="634" x2="628.57" y1="381.4"
+      /><line y2="338.6" style="fill:none;" x1="634" x2="628.57" y1="338.6"
+      /><line y2="295.8" style="fill:none;" x1="634" x2="628.57" y1="295.8"
+      /><line y2="253" style="fill:none;" x1="634" x2="628.57" y1="253"
+      /><line y2="210.2" style="fill:none;" x1="634" x2="628.57" y1="210.2"
+      /><line y2="167.4" style="fill:none;" x1="634" x2="628.57" y1="167.4"
+      /><line y2="124.5999" style="fill:none;" x1="634" x2="628.57" y1="124.5999"
+      /><line y2="81.8" style="fill:none;" x1="634" x2="628.57" y1="81.8"
+      /><line y2="39" style="fill:none;" x1="634" x2="628.57" y1="39"
+    /></g
+    ><g transform="translate(84.3333,467)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-24" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.5</text
+    ></g
+    ><g transform="translate(84.3333,424.2)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.52</text
+    ></g
+    ><g transform="translate(84.3333,381.4)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.54</text
+    ></g
+    ><g transform="translate(84.3333,338.6)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.56</text
+    ></g
+    ><g transform="translate(84.3333,295.8)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.58</text
+    ></g
+    ><g transform="translate(84.3333,253)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-24" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.6</text
+    ></g
+    ><g transform="translate(84.3333,210.2)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.62</text
+    ></g
+    ><g transform="translate(84.3333,167.4)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.64</text
+    ></g
+    ><g transform="translate(84.3333,124.5999)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.66</text
+    ></g
+    ><g transform="translate(84.3333,81.8)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.68</text
+    ></g
+    ><g transform="translate(84.3333,39)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-24" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.7</text
+    ></g
+    ><g style="stroke-linecap:butt; fill:rgb(0,114,189); text-rendering:geometricPrecision; image-rendering:optimizeQuality; color-rendering:optimizeQuality; stroke-linejoin:round; stroke:rgb(0,114,189); color-interpolation:linearRGB; stroke-width:0.8333;"
+    ><path d="M91 467 L119.579 373.4814 L148.1579 296.6363 L176.7368 237.4048 L205.3158 191.7797 L233.8947 156.7113 L262.4737 129.8598 L291.0526 109.4154 L319.6316 93.9697 L348.2105 82.4194 L376.7895 73.8966 L405.3684 67.7155 L433.9474 63.3316 L462.5263 60.312 L491.1053 58.3114 L519.6842 57.0543 L548.2632 56.3201 L576.8421 55.9331 L605.421 55.7526 L634 53.665" style="fill:none; fill-rule:evenodd;"
+      /><path d="M91 467 L119.579 367.9323 L148.1579 291.6423 L176.7368 232.8554 L205.3158 187.5936 L233.8947 152.8281 L262.4737 126.2331 L291.0526 106.0095 L319.6316 90.7562 L348.2105 79.3757 L376.7895 71.0043 L405.3684 64.9593 L433.9474 60.6987 L462.5263 57.7918 L491.1053 55.8947 L519.6842 54.7329 L548.2632 54.0872 L576.8421 53.7825 L605.421 53.6788 L634 53.665" style="fill:none; fill-rule:evenodd; stroke:rgb(217,83,25);"
+    /></g
+  ></g
+></svg
+>
diff --git a/Project/N_5.fig b/Project/N_5.fig
new file mode 100644 (file)
index 0000000..d622a05
Binary files /dev/null and b/Project/N_5.fig differ
diff --git a/Project/N_5.svg b/Project/N_5.svg
new file mode 100644 (file)
index 0000000..f5e0e04
--- /dev/null
@@ -0,0 +1,180 @@
+<?xml version="1.0"?>
+<!DOCTYPE svg PUBLIC '-//W3C//DTD SVG 1.0//EN'
+          'http://www.w3.org/TR/2001/REC-SVG-20010904/DTD/svg10.dtd'>
+<svg xmlns:xlink="http://www.w3.org/1999/xlink" style="fill-opacity:1; color-rendering:auto; color-interpolation:auto; text-rendering:auto; stroke:black; stroke-linecap:square; stroke-miterlimit:10; shape-rendering:auto; stroke-opacity:1; fill:black; stroke-dasharray:none; font-weight:normal; stroke-width:1; font-family:'Dialog'; font-style:normal; stroke-linejoin:miter; font-size:12px; stroke-dashoffset:0; image-rendering:auto;" width="700" height="525" xmlns="http://www.w3.org/2000/svg"
+><!--Generated by the Batik Graphics2D SVG Generator--><defs id="genericDefs"
+  /><g
+  ><defs id="defs1"
+    ><clipPath clipPathUnits="userSpaceOnUse" id="clipPath1"
+      ><path d="M0 0 L700 0 L700 525 L0 525 L0 0 Z"
+      /></clipPath
+      ><font horiz-adv-x="75.0" id="font1"
+      ><font-face ascent="100.53711" descent="21.972656" units-per-em="100" style="font-style:normal; font-family:Dialog; font-weight:normal;"
+        /><missing-glyph horiz-adv-x="75.0" d="M12.5 0 L12.5 62.5 L62.5 62.5 L62.5 0 L12.5 0 ZM14.0625 1.5625 L60.9375 1.5625 L60.9375 60.9375 L14.0625 60.9375 L14.0625 1.5625 Z"
+        /><glyph unicode="1" horiz-adv-x="55.615234" d="M37.25 0 L28.4688 0 L28.4688 56 Q25.2969 52.9844 20.1406 49.9531 Q14.9844 46.9219 10.8906 45.4062 L10.8906 53.9062 Q18.2656 57.375 23.7812 62.3047 Q29.2969 67.2344 31.5938 71.875 L37.25 71.875 L37.25 0 Z"
+        /><glyph unicode="." horiz-adv-x="27.783203" d="M9.0781 0 L9.0781 10.0156 L19.0938 10.0156 L19.0938 0 L9.0781 0 Z"
+        /><glyph unicode="2" horiz-adv-x="55.615234" d="M50.3438 8.4531 L50.3438 0 L3.0312 0 Q2.9375 3.1719 4.0469 6.1094 Q5.8594 10.9375 9.8359 15.625 Q13.8125 20.3125 21.3438 26.4688 Q33.0156 36.0312 37.1172 41.625 Q41.2188 47.2188 41.2188 52.2031 Q41.2188 57.4219 37.4766 61.0078 Q33.7344 64.5938 27.7344 64.5938 Q21.3906 64.5938 17.5781 60.7891 Q13.7656 56.9844 13.7188 50.25 L4.6875 51.1719 Q5.6094 61.2812 11.6641 66.5781 Q17.7188 71.875 27.9375 71.875 Q38.2344 71.875 44.2422 66.1641 Q50.25 60.4531 50.25 52 Q50.25 47.7031 48.4922 43.5547 Q46.7344 39.4062 42.6562 34.8125 Q38.5781 30.2188 29.1094 22.2188 Q21.1875 15.5781 18.9453 13.2109 Q16.7031 10.8438 15.2344 8.4531 L50.3438 8.4531 Z"
+        /><glyph unicode="3" horiz-adv-x="55.615234" d="M4.2031 18.8906 L12.9844 20.0625 Q14.5 12.5938 18.1406 9.2969 Q21.7812 6 27 6 Q33.2031 6 37.4766 10.2969 Q41.75 14.5938 41.75 20.9531 Q41.75 27 37.7969 30.9297 Q33.8438 34.8594 27.7344 34.8594 Q25.25 34.8594 21.5312 33.8906 L22.5156 41.6094 Q23.3906 41.5 23.9219 41.5 Q29.5469 41.5 34.0391 44.4297 Q38.5312 47.3594 38.5312 53.4688 Q38.5312 58.2969 35.2578 61.4766 Q31.9844 64.6562 26.8125 64.6562 Q21.6875 64.6562 18.2656 61.4297 Q14.8438 58.2031 13.875 51.7656 L5.0781 53.3281 Q6.6875 62.1562 12.3984 67.0156 Q18.1094 71.875 26.6094 71.875 Q32.4688 71.875 37.3984 69.3594 Q42.3281 66.8438 44.9453 62.5 Q47.5625 58.1562 47.5625 53.2656 Q47.5625 48.6406 45.0703 44.8281 Q42.5781 41.0156 37.7031 38.7656 Q44.0469 37.3125 47.5625 32.6953 Q51.0781 28.0781 51.0781 21.1406 Q51.0781 11.7656 44.2422 5.25 Q37.4062 -1.2656 26.9531 -1.2656 Q17.5312 -1.2656 11.3047 4.3516 Q5.0781 9.9688 4.2031 18.8906 Z"
+        /><glyph unicode="4" horiz-adv-x="55.615234" d="M32.3281 0 L32.3281 17.1406 L1.2656 17.1406 L1.2656 25.2031 L33.9375 71.5781 L41.1094 71.5781 L41.1094 25.2031 L50.7812 25.2031 L50.7812 17.1406 L41.1094 17.1406 L41.1094 0 L32.3281 0 ZM32.3281 25.2031 L32.3281 57.4688 L9.9062 25.2031 L32.3281 25.2031 Z"
+        /><glyph unicode="5" horiz-adv-x="55.615234" d="M4.1562 18.75 L13.375 19.5312 Q14.4062 12.7969 18.1406 9.3984 Q21.875 6 27.1562 6 Q33.5 6 37.8906 10.7891 Q42.2812 15.5781 42.2812 23.4844 Q42.2812 31 38.0625 35.3516 Q33.8438 39.7031 27 39.7031 Q22.75 39.7031 19.3359 37.7734 Q15.9219 35.8438 13.9688 32.7656 L5.7188 33.8438 L12.6406 70.6094 L48.25 70.6094 L48.25 62.2031 L19.6719 62.2031 L15.8281 42.9688 Q22.2656 47.4688 29.3438 47.4688 Q38.7188 47.4688 45.1641 40.9688 Q51.6094 34.4688 51.6094 24.2656 Q51.6094 14.5469 45.9531 7.4688 Q39.0625 -1.2188 27.1562 -1.2188 Q17.3906 -1.2188 11.2109 4.25 Q5.0312 9.7188 4.1562 18.75 Z"
+        /><glyph unicode="6" horiz-adv-x="55.615234" d="M49.75 54.0469 L41.0156 53.375 Q39.8438 58.5469 37.7031 60.8906 Q34.125 64.6562 28.9062 64.6562 Q24.7031 64.6562 21.5312 62.3125 Q17.3906 59.2812 14.9922 53.4688 Q12.5938 47.6562 12.5 36.9219 Q15.6719 41.75 20.2656 44.0938 Q24.8594 46.4375 29.8906 46.4375 Q38.6719 46.4375 44.8516 39.9688 Q51.0312 33.5 51.0312 23.25 Q51.0312 16.5 48.125 10.7188 Q45.2188 4.9375 40.1406 1.8594 Q35.0625 -1.2188 28.6094 -1.2188 Q17.625 -1.2188 10.6953 6.8594 Q3.7656 14.9375 3.7656 33.5 Q3.7656 54.25 11.4219 63.6719 Q18.1094 71.875 29.4375 71.875 Q37.8906 71.875 43.2891 67.1406 Q48.6875 62.4062 49.75 54.0469 ZM13.875 23.1875 Q13.875 18.6562 15.7969 14.5078 Q17.7188 10.3594 21.1875 8.1797 Q24.6562 6 28.4688 6 Q34.0312 6 38.0391 10.4922 Q42.0469 14.9844 42.0469 22.7031 Q42.0469 30.125 38.0859 34.3984 Q34.125 38.6719 28.125 38.6719 Q22.1719 38.6719 18.0234 34.3984 Q13.875 30.125 13.875 23.1875 Z"
+        /><glyph unicode="7" horiz-adv-x="55.615234" d="M4.7344 62.2031 L4.7344 70.6562 L51.0781 70.6562 L51.0781 63.8125 Q44.2344 56.5469 37.5234 44.4844 Q30.8125 32.4219 27.1562 19.6719 Q24.5156 10.6875 23.7812 0 L14.75 0 Q14.8906 8.4531 18.0625 20.4141 Q21.2344 32.375 27.1719 43.4844 Q33.1094 54.5938 39.7969 62.2031 L4.7344 62.2031 Z"
+        /><glyph unicode="8" horiz-adv-x="55.615234" d="M17.6719 38.8125 Q12.2031 40.8281 9.5703 44.5391 Q6.9375 48.25 6.9375 53.4219 Q6.9375 61.2344 12.5547 66.5547 Q18.1719 71.875 27.4844 71.875 Q36.8594 71.875 42.5781 66.4297 Q48.2969 60.9844 48.2969 53.1719 Q48.2969 48.1875 45.6797 44.5078 Q43.0625 40.8281 37.75 38.8125 Q44.3438 36.6719 47.7812 31.8828 Q51.2188 27.0938 51.2188 20.4531 Q51.2188 11.2812 44.7266 5.0312 Q38.2344 -1.2188 27.6406 -1.2188 Q17.0469 -1.2188 10.5469 5.0547 Q4.0469 11.3281 4.0469 20.7031 Q4.0469 27.6875 7.5938 32.3984 Q11.1406 37.1094 17.6719 38.8125 ZM15.9219 53.7188 Q15.9219 48.6406 19.1953 45.4141 Q22.4688 42.1875 27.6875 42.1875 Q32.7656 42.1875 36.0156 45.3828 Q39.2656 48.5781 39.2656 53.2188 Q39.2656 58.0625 35.9141 61.3594 Q32.5625 64.6562 27.5938 64.6562 Q22.5625 64.6562 19.2422 61.4297 Q15.9219 58.2031 15.9219 53.7188 ZM13.0938 20.6562 Q13.0938 16.8906 14.875 13.375 Q16.6562 9.8594 20.1719 7.9297 Q23.6875 6 27.7344 6 Q34.0312 6 38.1328 10.0547 Q42.2344 14.1094 42.2344 20.3594 Q42.2344 26.7031 38.0156 30.8594 Q33.7969 35.0156 27.4375 35.0156 Q21.2344 35.0156 17.1641 30.9141 Q13.0938 26.8125 13.0938 20.6562 Z"
+        /><glyph unicode="9" horiz-adv-x="55.615234" d="M5.4688 16.5469 L13.9219 17.3281 Q14.9844 11.375 18.0156 8.6875 Q21.0469 6 25.7812 6 Q29.8281 6 32.8828 7.8594 Q35.9375 9.7188 37.8906 12.8203 Q39.8438 15.9219 41.1641 21.1953 Q42.4844 26.4688 42.4844 31.9375 Q42.4844 32.5156 42.4375 33.6875 Q39.7969 29.5 35.2344 26.8828 Q30.6719 24.2656 25.3438 24.2656 Q16.4531 24.2656 10.3047 30.7109 Q4.1562 37.1562 4.1562 47.7031 Q4.1562 58.5938 10.5781 65.2344 Q17 71.875 26.6562 71.875 Q33.6406 71.875 39.4297 68.1172 Q45.2188 64.3594 48.2188 57.3984 Q51.2188 50.4375 51.2188 37.25 Q51.2188 23.5312 48.2422 15.4062 Q45.2656 7.2812 39.3828 3.0312 Q33.5 -1.2188 25.5938 -1.2188 Q17.1875 -1.2188 11.8672 3.4453 Q6.5469 8.1094 5.4688 16.5469 ZM41.4531 48.1406 Q41.4531 55.7188 37.4297 60.1562 Q33.4062 64.5938 27.7344 64.5938 Q21.875 64.5938 17.5312 59.8125 Q13.1875 55.0312 13.1875 47.4062 Q13.1875 40.5781 17.3125 36.3047 Q21.4375 32.0312 27.4844 32.0312 Q33.5938 32.0312 37.5234 36.3047 Q41.4531 40.5781 41.4531 48.1406 Z"
+        /><glyph unicode="0" horiz-adv-x="55.615234" d="M4.1562 35.2969 Q4.1562 48 6.7656 55.7422 Q9.375 63.4844 14.5234 67.6797 Q19.6719 71.875 27.4844 71.875 Q33.25 71.875 37.5938 69.5547 Q41.9375 67.2344 44.7734 62.8672 Q47.6094 58.5 49.2188 52.2266 Q50.8281 45.9531 50.8281 35.2969 Q50.8281 22.7031 48.2422 14.9688 Q45.6562 7.2344 40.5078 3.0078 Q35.3594 -1.2188 27.4844 -1.2188 Q17.1406 -1.2188 11.2344 6.2031 Q4.1562 15.1406 4.1562 35.2969 ZM13.1875 35.2969 Q13.1875 17.6719 17.3125 11.8359 Q21.4375 6 27.4844 6 Q33.5469 6 37.6719 11.8594 Q41.7969 17.7188 41.7969 35.2969 Q41.7969 52.9844 37.6719 58.7891 Q33.5469 64.5938 27.3906 64.5938 Q21.3438 64.5938 17.7188 59.4688 Q13.1875 52.9375 13.1875 35.2969 Z"
+      /></font
+    ></defs
+    ><g style="fill:white; stroke:white;"
+    ><rect x="0" y="0" width="700" style="clip-path:url(#clipPath1); stroke:none;" height="525"
+    /></g
+    ><g style="fill:white; text-rendering:optimizeSpeed; color-rendering:optimizeSpeed; image-rendering:optimizeSpeed; shape-rendering:crispEdges; stroke:white; color-interpolation:sRGB;"
+    ><rect x="0" width="700" height="525" y="0" style="stroke:none;"
+      /><path style="stroke:none;" d="M91 467 L634 467 L634 39 L91 39 Z"
+    /></g
+    ><g style="fill:rgb(38,38,38); text-rendering:geometricPrecision; image-rendering:optimizeQuality; color-rendering:optimizeQuality; stroke-linejoin:round; stroke:rgb(38,38,38); color-interpolation:linearRGB; stroke-width:0.8333;"
+    ><line y2="467" style="fill:none;" x1="91" x2="634" y1="467"
+      /><line y2="39" style="fill:none;" x1="91" x2="634" y1="39"
+      /><line y2="461.57" style="fill:none;" x1="91" x2="91" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="145.3" x2="145.3" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="199.6" x2="199.6" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="253.9" x2="253.9" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="308.2" x2="308.2" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="362.5" x2="362.5" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="416.8" x2="416.8" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="471.1" x2="471.1" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="525.4" x2="525.4" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="579.7" x2="579.7" y1="467"
+      /><line y2="461.57" style="fill:none;" x1="634" x2="634" y1="467"
+      /><line y2="44.43" style="fill:none;" x1="91" x2="91" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="145.3" x2="145.3" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="199.6" x2="199.6" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="253.9" x2="253.9" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="308.2" x2="308.2" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="362.5" x2="362.5" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="416.8" x2="416.8" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="471.1" x2="471.1" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="525.4" x2="525.4" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="579.7" x2="579.7" y1="39"
+      /><line y2="44.43" style="fill:none;" x1="634" x2="634" y1="39"
+    /></g
+    ><g transform="translate(91,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-5" xml:space="preserve" y="17" style="stroke:none;"
+      >1</text
+    ></g
+    ><g transform="translate(145.3,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.1</text
+    ></g
+    ><g transform="translate(199.6,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.2</text
+    ></g
+    ><g transform="translate(253.9,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.3</text
+    ></g
+    ><g transform="translate(308.2,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.4</text
+    ></g
+    ><g transform="translate(362.5,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.5</text
+    ></g
+    ><g transform="translate(416.8,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.6</text
+    ></g
+    ><g transform="translate(471.1,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.7</text
+    ></g
+    ><g transform="translate(525.4,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.8</text
+    ></g
+    ><g transform="translate(579.7,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-12" xml:space="preserve" y="17" style="stroke:none;"
+      >1.9</text
+    ></g
+    ><g transform="translate(634,473.6667)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-5" xml:space="preserve" y="17" style="stroke:none;"
+      >2</text
+    ></g
+    ><g style="fill:rgb(38,38,38); text-rendering:geometricPrecision; image-rendering:optimizeQuality; color-rendering:optimizeQuality; stroke-linejoin:round; stroke:rgb(38,38,38); color-interpolation:linearRGB; stroke-width:0.8333;"
+    ><line y2="39" style="fill:none;" x1="91" x2="91" y1="467"
+      /><line y2="39" style="fill:none;" x1="634" x2="634" y1="467"
+      /><line y2="467" style="fill:none;" x1="91" x2="96.43" y1="467"
+      /><line y2="424.2" style="fill:none;" x1="91" x2="96.43" y1="424.2"
+      /><line y2="381.4" style="fill:none;" x1="91" x2="96.43" y1="381.4"
+      /><line y2="338.6" style="fill:none;" x1="91" x2="96.43" y1="338.6"
+      /><line y2="295.8" style="fill:none;" x1="91" x2="96.43" y1="295.8"
+      /><line y2="253" style="fill:none;" x1="91" x2="96.43" y1="253"
+      /><line y2="210.2" style="fill:none;" x1="91" x2="96.43" y1="210.2"
+      /><line y2="167.4" style="fill:none;" x1="91" x2="96.43" y1="167.4"
+      /><line y2="124.5999" style="fill:none;" x1="91" x2="96.43" y1="124.5999"
+      /><line y2="81.8" style="fill:none;" x1="91" x2="96.43" y1="81.8"
+      /><line y2="39" style="fill:none;" x1="91" x2="96.43" y1="39"
+      /><line y2="467" style="fill:none;" x1="634" x2="628.57" y1="467"
+      /><line y2="424.2" style="fill:none;" x1="634" x2="628.57" y1="424.2"
+      /><line y2="381.4" style="fill:none;" x1="634" x2="628.57" y1="381.4"
+      /><line y2="338.6" style="fill:none;" x1="634" x2="628.57" y1="338.6"
+      /><line y2="295.8" style="fill:none;" x1="634" x2="628.57" y1="295.8"
+      /><line y2="253" style="fill:none;" x1="634" x2="628.57" y1="253"
+      /><line y2="210.2" style="fill:none;" x1="634" x2="628.57" y1="210.2"
+      /><line y2="167.4" style="fill:none;" x1="634" x2="628.57" y1="167.4"
+      /><line y2="124.5999" style="fill:none;" x1="634" x2="628.57" y1="124.5999"
+      /><line y2="81.8" style="fill:none;" x1="634" x2="628.57" y1="81.8"
+      /><line y2="39" style="fill:none;" x1="634" x2="628.57" y1="39"
+    /></g
+    ><g transform="translate(84.3333,467)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-24" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.5</text
+    ></g
+    ><g transform="translate(84.3333,424.2)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.52</text
+    ></g
+    ><g transform="translate(84.3333,381.4)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.54</text
+    ></g
+    ><g transform="translate(84.3333,338.6)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.56</text
+    ></g
+    ><g transform="translate(84.3333,295.8)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.58</text
+    ></g
+    ><g transform="translate(84.3333,253)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-24" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.6</text
+    ></g
+    ><g transform="translate(84.3333,210.2)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.62</text
+    ></g
+    ><g transform="translate(84.3333,167.4)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.64</text
+    ></g
+    ><g transform="translate(84.3333,124.5999)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.66</text
+    ></g
+    ><g transform="translate(84.3333,81.8)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-33" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.68</text
+    ></g
+    ><g transform="translate(84.3333,39)" style="font-size:16.6667px; fill:rgb(38,38,38); text-rendering:geometricPrecision; color-rendering:optimizeQuality; image-rendering:optimizeQuality; stroke:rgb(38,38,38); color-interpolation:linearRGB;"
+    ><text x="-24" xml:space="preserve" y="6.5" style="stroke:none;"
+      >0.7</text
+    ></g
+    ><g style="stroke-linecap:butt; fill:rgb(0,114,189); text-rendering:geometricPrecision; image-rendering:optimizeQuality; color-rendering:optimizeQuality; stroke-linejoin:round; stroke:rgb(0,114,189); color-interpolation:linearRGB; stroke-width:0.8333;"
+    ><path d="M91 467 L226.75 263.5816 L362.5 145.4532 L498.25 108.2161 L634 53.665" style="fill:none; fill-rule:evenodd;"
+      /><path d="M91 467 L226.75 160.6728 L362.5 74.8603 L498.25 55.5447 L634 53.665" style="fill:none; fill-rule:evenodd; stroke:rgb(217,83,25);"
+    /></g
+  ></g
+></svg
+>
diff --git a/Project/lufact.m b/Project/lufact.m
new file mode 100644 (file)
index 0000000..5ca3886
--- /dev/null
@@ -0,0 +1,49 @@
+function [L,U,P] = lufact(A)
+%%
+%  Computes the LU factorization of the matrix A. The factorization is done
+%  in-place and then the L and U matrices are extracted at the end for
+%  output. The factorization is PA = LU
+%
+%  INPUT: 
+%    A - n by n square, non-singular matrix
+%
+%  OUTPUT:
+%    L - lower triangular matrix with ones along the main diagonal
+%    U - upper triangular matrix
+%    P - permutation matrix for pivoting
+%%
+%  get the size of the system
+   n  = length(A);
+%%
+%  initialize the pivoting matrix
+   P = eye(n);
+%%
+%  Vectorized in-place LU factorization (with row pivoting) that keeps
+%  track of the total permutations by scrambleing the matrix P
+   for j = 1:n-1
+        A
+        [v i] = max(abs(A(j:n, j)));
+        A([j i+j-1],:) = A([i+j-1 j],:);
+        P([j i+j-1],:) = P([i+j-1 j],:);
+        A
+%%
+% Find the index of the largest pivot element in the current column
+
+
+%%
+% Swap the rows within the in-place array as well as the permutation matrix P
+    
+
+
+%%
+% Perform the in-place elimination and save the new column of L
+      i = j+1:n; % indices for the "active" matrix portion
+      A(i,j) = A(i,j)/A(j,j);
+      A(i,i) = A(i,i) - A(i,j)*A(j,i);
+      %A, return
+   end
+%%
+% Extract L and U from the in-place form
+   U = triu(A);
+   L = eye(n) + tril(A,-1);
+end
diff --git a/Project/script.m b/Project/script.m
new file mode 100644 (file)
index 0000000..1790f1f
--- /dev/null
@@ -0,0 +1,59 @@
+%% Section 1
+
+N = 5;
+
+% create three random vectors
+main = N*ones(N,1) + rand(N,1);
+upper = rand(N-1,1);
+lower = rand(N-1,1);
+% use the ’diag’ command to create tridiagonal matrix
+A = diag(main,0) + diag(upper,1) - diag(lower,-1);
+
+known_x = ones(N,1); % make known solution vector of 1s
+manuf_b = A*known_x; % manufacture the right-hand-side
+
+A
+
+%% Section 2
+
+N = 5;
+x = linspace(1,2, N);
+h = x(2) - x(1); 
+
+% Build tridiag matrix A
+main = ones(N-2,1);
+for i=2:N-1
+    main(i-1,1) = (2/(x(i)^2)) - (2/(h^2));
+end
+lower = ones(N-3,1);
+for i=2:N-2
+    lower(i-1,1) = (1/(h^2)) - (2/(x(i+1)*h));
+end
+upper = ones(N-3,1);
+for i=2:N-2
+    upper(i-1,1) = (1/(h^2)) + (2/(x(i)*h));
+end
+
+A = diag(main,0) + diag(upper,1) + diag(lower,-1);
+
+% Build constant matrix b
+b = zeros(N-2, 1);
+for i=2:N-1
+    b(i-1) = (2*log(x(i)))/(x(i)^2);
+end
+b(1) = b(1) + -((1/(h^2)) - (2/(x(1)*h)))*0.5;
+b(end) = b(end) + -((1/(h^2)) + (2/(x(end)*h)))*log(2);
+
+answer = thomas(A, b);
+answer = [0.5;answer;log(2)];
+
+plot(x, answer);
+hold on;
+real_values = zeros(N, 1);
+for i=1:N
+    real_values(i) = (4/x(i)) - (2/(x(i)^2)) + log(x(i)) - (3/2);
+end
+plot(x, real_values);
+
+
+norm(answer - real_values, inf)
\ No newline at end of file
diff --git a/Project/thomas.asv b/Project/thomas.asv
new file mode 100644 (file)
index 0000000..f4c60ae
--- /dev/null
@@ -0,0 +1,21 @@
+function [x] = thomas(A,b)
+    % from wikipedia
+    n = length(b);
+    c_prime = zeros(n-1, 1);
+    d_prime = zeros(n, 1);
+    c_prime(1) = A(1,2) / A(1,1);
+    d_prime(1) = b(1) / A(1,1);
+    for i=2:n-1
+        c_prime(i) = A(i,1 + i) / (A(i,i) - (A(i, i-1)*c_prime(i-1)));
+        d_prime(i) = (b(i) - (A(i, i-1)*d_prime(i-1))) / (A(i, i) - (A(i, i-1)*c_prime(i-1)));
+    end
+    d_prime(n) = (b(n) - (A(n, n-1)*d_prime(n-1))) / (A(n, n) - (A(n, n-1)*c_prime(n-1)));
+    
+    % backward sub
+    x = zeros(n, 1);
+    x(n) = d_prime(n);
+    for j=n-1:-1:1
+        x(j) = d_prime(j) - (c_prime(j)*x(j+1));
+    end
+
+end
\ No newline at end of file
diff --git a/Project/thomas.m b/Project/thomas.m
new file mode 100644 (file)
index 0000000..f4c60ae
--- /dev/null
@@ -0,0 +1,21 @@
+function [x] = thomas(A,b)
+    % from wikipedia
+    n = length(b);
+    c_prime = zeros(n-1, 1);
+    d_prime = zeros(n, 1);
+    c_prime(1) = A(1,2) / A(1,1);
+    d_prime(1) = b(1) / A(1,1);
+    for i=2:n-1
+        c_prime(i) = A(i,1 + i) / (A(i,i) - (A(i, i-1)*c_prime(i-1)));
+        d_prime(i) = (b(i) - (A(i, i-1)*d_prime(i-1))) / (A(i, i) - (A(i, i-1)*c_prime(i-1)));
+    end
+    d_prime(n) = (b(n) - (A(n, n-1)*d_prime(n-1))) / (A(n, n) - (A(n, n-1)*c_prime(n-1)));
+    
+    % backward sub
+    x = zeros(n, 1);
+    x(n) = d_prime(n);
+    for j=n-1:-1:1
+        x(j) = d_prime(j) - (c_prime(j)*x(j+1));
+    end
+
+end
\ No newline at end of file
diff --git a/RungeKutta4.m b/RungeKutta4.m
new file mode 100644 (file)
index 0000000..c7cecff
--- /dev/null
@@ -0,0 +1,44 @@
+function [t,y] = RungeKutta4(f,t0,T,y0,h)
+%%
+%  Classical four-stage Runge-Kutta time integration technique
+%  to solve the initial value problem
+%
+%     y' = f(t,y),  y(t0) = y0
+%
+%  INPUT:
+%
+%     f  - right-hand-side function that can depend on
+%          the time t and the function y
+%     t0 - intial time
+%     T  - final time
+%     y0 - initial solution value
+%     h  - fixed time step size
+%
+%  OUTPUT:
+%
+%     t - vector of time values
+%     y - solution vector of approximate values
+%
+%%
+%  discretize the time variable and set the problem size
+   t = t0:h:T;
+   y = zeros(length(t),1);
+%%
+%  save the intial values
+   y(1) = y0;
+   t(1) = t0;
+%%
+%  use Runge-Kutta 4 stage with the given timestep h to solve the problem and
+%  save the information for later analysis
+   for j = 2:length(t)
+        k_1 = f(t(j-1), y(j-1));
+        k_2 = f(t(j-1) + (h/2), y(j-1) + ((h/2)*k_1));
+        k_3 = f(t(j-1) + (h/2), y(j-1) + ((h/2)*k_2));
+        k_4 = f(t(j-1) + h, y(j-1) + ((h*k_3)));
+        y(j) = y(j-1) + ((h/6)*(k_1 + (2*k_2) + (2*k_3) + k_4));
+
+
+
+
+   end
+end
diff --git a/SimpsonsRule.m b/SimpsonsRule.m
new file mode 100644 (file)
index 0000000..cecf9ad
--- /dev/null
@@ -0,0 +1,72 @@
+function [S,x] = SimpsonsRule(f,a,b,n)
+%%
+%  Implementation of the (possibly) composite Simpsons's quadrature
+%  rule to approximate the definite integral
+%
+%           b
+%          /\
+%          |
+%          | f(x) dx
+%          |
+%         \/
+%         a
+%
+%  INPUT:
+%
+%     f    - anonymous function integrand
+%     a    - lower integration bound
+%     b    - upper integration bound
+%     n    - number of intervals (must be an even integer)
+%
+%  OUTPUT:
+%
+%     In - approximation to the definite integral on the given interval [a,b]
+%     x  - set of n+1 quadrature nodes
+%
+%%
+%  Simpson's rule needs at least two intervals
+   if n == 1
+      error('Simpsons rule requires n >= 2');
+   end
+%%
+%  check if the number of intervals is even
+   if mod(n,2) ~= 0
+      error('The number of intervals n must be even');
+   end
+%%
+%  Find the interval spacing h
+   h = (b-a)/n;
+
+%%
+%  build the set of uniformly spaced quadrature nodes
+   x = zeros(n+1,1);
+   x(1) = a;
+   for i=2:n+1
+       x(i) = x(i-1) + h;
+
+%%
+%  Initialize the integral value to zero
+   S = f(x(1)) + f(x(end));
+   for i=2:n
+       if mod(i, 2) == 0
+            S = S + 4*f(x(i));
+       else
+            S = S + 2*f(x(i));
+       end
+   end
+   
+   S = (h/3) * S;
+
+%%
+%  Approximate the integral with (composite) Simpson's rule
+
+
+
+
+
+
+
+
+
+
+end
diff --git a/adaptativeSimpson.m b/adaptativeSimpson.m
new file mode 100644 (file)
index 0000000..5b96887
--- /dev/null
@@ -0,0 +1,67 @@
+function [In,t] = adaptativeSimpson(f,a,b,tol)
+%%
+%  Adaptive Simpson's rule to approximate a definite integral. This uses
+%  a posterioi error analysis to recursively apply Simpson's rule where it
+%  is needed
+%
+%  INPUT:
+%
+%     f   - integrand function
+%     a   - lower bound of interval
+%     b   - upper bound of the interval
+%     tol - error tolerance for adaptivity
+%
+%  OUTPUT:
+%
+%     In - adaptive approximation of the definite integral to
+%          within the error tolerance
+%     t  - vector of the adapted quadrature nodes
+%
+%%
+%%
+%  use recursive bisection with error estimation to compute the integral approximation
+   m      = 0.5*(b+a); % find the current midpoint
+   [In,t] = do_integral(a,f(a),b,f(b),m,f(m),tol);
+%%
+%  this is a recursive function within Matlab so it calls itself
+      function [In,t] = do_integral(a,fa,b,fb,m,fm,tol)
+      %%
+      %  need the two midpoints of the sub-intervals and the function
+      %  evalutions for the recursion
+         xL  = 0.5*(a+m);
+         fxL = f(xL);
+         xR = 0.5*(b+m);
+         fxR = f(xR);
+      %%
+      %  save the five nodes at the current level of the recursion
+         t = [a;xL;m;xR;b];
+      %%
+      %  get the h value for the whole interval
+         h = 0.5*(b-a);
+      %%
+      %  compute the Simpson's rule on the whole interval
+         S_coarse = h/3*(fa + 4*fm + fb);
+      %%
+      %  compute the Simpson's rule on the two subintervals
+         S_Left  = h/6*(fa + 4*fxL + fm);
+         S_Right = h/6*(fm + 4*fxR + fb);
+         S_fine  = S_Left + S_Right;
+      %%
+      %  error estimate
+         E = S_coarse - S_fine;
+      %%
+      %  check against the user set tolerance
+         if abs(E) < 10*tol
+      % exit when tolerence is met
+            In = S_fine;
+         else
+      %%
+      %  bisect again if the tolernce is note met
+            [IL,tL] = do_integral(a,fa,m,fm,xL,fxL,tol);
+            [IR,tR] = do_integral(m,fm,b,fb,xR,fxR,tol);
+            In = IL + IR;
+      % append the node edges together (this avoids duplicates)
+            t = [tL;tR(2:end)];
+         end
+      end
+end
diff --git a/forwardEuler.m b/forwardEuler.m
new file mode 100644 (file)
index 0000000..fa91267
--- /dev/null
@@ -0,0 +1,37 @@
+function [t,y] = forwardEuler(f,t0,T,y0,h)
+%%
+%  Forward Euler time integration technique to solve the
+%  initial value problem
+%
+%     y' = f(t,y),  y(t0) = y0
+%
+%  INPUT:
+%
+%     f  - right-hand-side function that can depend on
+%          the time t and the function y
+%     t0 - intial time
+%     T  - final time
+%     y0 - initial solution value
+%     h  - fixed time step size
+%
+%  OUTPUT:
+%
+%     t - vector of time values
+%     y - solution vector of approximate values
+%
+%%
+%  discretize the time variable and set the problem size
+   t = t0:h:T;
+   y = zeros(length(t),1);
+%%
+%  save the intial values
+   y(1) = y0;
+   t(1) = t0;
+%%
+%  use forward Euler with the given timestep h to solve the problem and
+%  save the information for later analysis
+   for j = 2:length(t)
+%  this is an explicit method so we always use values from the current solution
+      y(j) = y(j-1) + h*f(t(j-1),y(j-1));
+   end
+end
diff --git a/lab3_k.m b/lab3_k.m
new file mode 100644 (file)
index 0000000..372f02f
--- /dev/null
+++ b/lab3_k.m
@@ -0,0 +1,30 @@
+%lab3
+N = 1000;
+x = linspace(0,1,N+1);
+h = x(2) - x(1);
+f_prime = zeros(1, N+1);
+
+f_1 = @(x) (1);
+f_2 = @(x) (x);
+f_3 = @(x) (x.^2);
+f_4 = @(x) (exp(sin(4*x)));
+
+for i=1:N-1
+    f_prime(i+1) = (f_4(x(i+2)) - f_4(x(i))) / (2*h);
+end
+
+f_prime(1) = (f_4(x(2)) - (f_4(x(1)))) / (h);
+f_prime(end) = (f_4(x(end)) - f_4(x(end-1))) / h;
+
+%f_prime(1) = ((-3*f_4(x(1))) + (4*f_4(x(2))) - (f_4(x(3)))) / (2*h);
+%f_prime(end) = ((f_4(x(end-2))) - (4*f_4(x(end-1))) + (3*f_4(x(end)))) / (2*h);
+
+f_prime_exact = zeros(1, N+1);
+for i=1:N+1
+    f_prime_exact(i) = 4*exp(sin(4*x(i))) * cos(4*x(i));
+end
+
+f_prime
+f_prime_exact
+
+max(abs(f_prime_exact - f_prime))
\ No newline at end of file
diff --git a/plot_composite_quad.m b/plot_composite_quad.m
new file mode 100644 (file)
index 0000000..d854aaa
--- /dev/null
@@ -0,0 +1,24 @@
+function plot_composite_quad(f,t)
+%%
+%  Script that will plot a function and the interval boundaries for use
+%  with a compositie quadrature rule
+%
+%  INPUT:
+%
+%     f - integrand function
+%     t - set of nodes used for the composite quadrature rule
+%
+%  OUTPUT:
+%
+%     NONE - Produces plot to screen
+%
+%%
+   clf
+   p1 = fplot(f,[t(1) t(end)],'-k','LineWidth',2);
+   hold on
+   p2 = stem(t,f(t),'fill','color',[0.9100    0.4100    0.1700],'MarkerSize',9,'LineWidth',1.5);
+   set(gcf,'Position',[500, 60, 1250, 1250])
+   set(gca,'FontSize',16);
+   legend([p1 p2],{'$(x+1)^2\cos\left(\frac{2x+1}{x-3.3}\right)$','Quadrature nodes'}...
+          ,'interpreter','latex','fontsize',24,'Location','northwest')
+end
\ No newline at end of file