diff --git a/HW6Problem4.m b/HW6Problem4.m
index fd9e0c7..e7712a9 100644
--- a/HW6Problem4.m
+++ b/HW6Problem4.m
@@ -1,27 +1,34 @@
 %part a.
-t = [0:.1:12];
 g = 9.81; %m/s^2
 cd = 0.25; %kg/m
 m = 60; %kg
 v_term = sqrt(m*g/cd);
 x_analytical = @(t) ((v_term^2)/g)*log(abs(cosh((g*t)/v_term))) + 100;
 
-
-plot(t, x_analytical(t))
-
 %part b.
-%dy = ode4(0, [1;0]);
-dt = [0.1 0.01];
-figure()
-hold on
-for k = 1:length(dt)
-    t = [0:dt(k):12]';
-    y_n = zeros(length(t), 2);
-    y_n(1, :) = [100 0]; %initial condition
-    for i = 2:length(t)
-        dy = ode4(t(i-1), y_n(i-1, :));
-        y_n(i, :) = y_n(i-1, :) + dt(k)*dy;
-    end
-    %[t23, y23] = ode23(@(t,y) phugoid_ode(t,y), [0,20], [10 0 0 2]);
-    plot(t, x_analytical(t), '.', t, y_n(:,1), 's');
-end
\ No newline at end of file
+dt = 0.01; %step size
+a = figure(1);
+set(0, 'defaultAxesFontsize', 16)
+set(0, 'defaultTextFontsize', 16)
+set(0, 'defaultLineLineWidth', 2)
+t = [0:dt:12]';
+y_n = zeros(length(t), 2);
+y_n(1, :) = [100 0]; %initial condition
+for i = 2:length(t)
+    dy = ode4(t, y_n(i-1, :));
+	y_n(i, :) = y_n(i-1, :) + dt(k)*dy;
+end
+plot(t, x_analytical(t), 'k-', t, y_n(:,1), 'o-');
+title('Analytical vs. Euler')
+legend('Analytical', 'Euler', 'Location', 'Northeast')
+saveas(a, 'euler4.png');
+
+%Part c.
+b = figure(2);
+dt = 0.01;
+t = [0:dt:12];
+y_heun = heun_sol_order2(@(t,y) ode4(t,y), dt, [100 0], [0 12]);
+plot(t, x_analytical(t), 'k-', t, y_heun(:,1), 'o-');
+title('Analytical vs. Heun')
+legend('Analytical', 'Heun', 'Location', 'Northeast')
+saveas(b, 'heun4.png');
\ No newline at end of file
diff --git a/README.md b/README.md
index 4dc6d88..168ac11 100644
--- a/README.md
+++ b/README.md
@@ -82,4 +82,50 @@ saveas(b, 'heun3.png');
 ### Part c.
 ![Analytical vs. Heun](./heun3.png)
 
-## Problem 4
\ No newline at end of file
+## Problem 4
+
+### Script:
+```matlab
+%part a.
+g = 9.81; %m/s^2
+cd = 0.25; %kg/m
+m = 60; %kg
+v_term = sqrt(m*g/cd);
+x_analytical = @(t) ((v_term^2)/g)*log(abs(cosh((g*t)/v_term))) + 100;
+
+%part b.
+dt = 0.01; %step size
+a = figure(1);
+set(0, 'defaultAxesFontsize', 16)
+set(0, 'defaultTextFontsize', 16)
+set(0, 'defaultLineLineWidth', 2)
+t = [0:dt:12]';
+y_n = zeros(length(t), 2);
+y_n(1, :) = [100 0]; %initial condition
+for i = 2:length(t)
+    dy = ode4(t, y_n(i-1, :));
+	y_n(i, :) = y_n(i-1, :) + dt(k)*dy;
+end
+plot(t, x_analytical(t), 'k-', t, y_n(:,1), 'o-');
+title('Analytical vs. Euler')
+legend('Analytical', 'Euler', 'Location', 'Northeast')
+saveas(a, 'euler4.png');
+
+%Part c.
+b = figure(2);
+dt = 0.01;
+t = [0:dt:12];
+y_heun = heun_sol_order2(@(t,y) ode4(t,y), dt, [100 0], [0 12]);
+plot(t, x_analytical(t), 'k-', t, y_heun(:,1), 'o-');
+title('Analytical vs. Heun')
+legend('Analytical', 'Heun', 'Location', 'Northeast')
+saveas(b, 'heun4.png');
+```
+
+### Part b.
+![Analytical vs. Euler](./euler4.png)
+
+### Part c.
+![Analytical vs. Heun](./heun4.png)
+
+
diff --git a/euler4.png b/euler4.png
new file mode 100644
index 0000000..38b90cb
Binary files /dev/null and b/euler4.png differ
diff --git a/heun3.png b/heun3.png
index f20321d..f355257 100644
Binary files a/heun3.png and b/heun3.png differ
diff --git a/heun4.png b/heun4.png
new file mode 100644
index 0000000..99d9e4f
Binary files /dev/null and b/heun4.png differ