Skip to content

Commit

Permalink
Browse files Browse the repository at this point in the history
editing
  • Loading branch information
nin13001 committed Dec 6, 2017
1 parent ad5e7ed commit f97cd0a
Show file tree
Hide file tree
Showing 15 changed files with 102 additions and 8 deletions.
36 changes: 31 additions & 5 deletions README.md
Expand Up @@ -16,9 +16,10 @@ for i=2:length(x)
dy=prob2_ode(x,y_n(i-1));
y_n(i)=y_n(i-1)+ dt(2)*dy;
end
plot(x,exp(-x),'k-',x,y_n,'o');
plot(x,exp(-x),x,y_n)
legend('Analytical','Euler','Location','Northeast')
```
![](assets/README-66741291.png)
![](assets/README-be88b588.png)

## Part C

Expand Down Expand Up @@ -46,7 +47,7 @@ y(1)=y0;
end
end
```
![](assets/README-23f54e97.png)
![](assets/README-8b7eb1c9.png)

# Problem 3
## Part A
Expand All @@ -64,6 +65,31 @@ for i=2:length(x)
dy=prob3_ode(x,y_n(i-1,:));
y_n(i,:)=y_n(i-1,:)+ dt(2)*dy;
end
plot(x,cos(3*x),'k-',x,y_n(:,1),'o');
plot(x,cos(3*x),x,y_n(:,1));
```
![](assets/README-4e4bb548.png)
![](assets/README-bec1ab77.png)

## Part C
```
function y = prob3ptC(prob3_ode, dt, y0,tlimits)
t = [tlimits(1):dt:tlimits(2)];
y = zeros(length(t),length(y0));
y(1,:) = y0;
for i = 2:length(t)
dy = prob3_ode(t(i-1),y(i-1,:));
y(i,:) = y(i-1,:)+dt*dy;
dyc = (dy + prob3_ode(t(i),y(i,:)))/2;
y(i,:) = y(i-1,:)+dt*dyc;
n = 1;
while (1)
n = n+1;
yold = y(i,:);
dyc = (dy+prob3_ode(t(i),y(i,:)))/2;
y(i,:) = y(i-1,:)+dt*dyc;
ea = abs(y(i,:)- yold)./y(i,:)*100;
if max(ea) < 0.0001 || n>100; break, end
end
end
end
```
![](assets/README-fcbec9a1.png)
Binary file removed assets/README-23f54e97.png
Binary file not shown.
Binary file removed assets/README-4e4bb548.png
Binary file not shown.
Binary file removed assets/README-66741291.png
Binary file not shown.
Binary file added assets/README-8b7eb1c9.png
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file added assets/README-be88b588.png
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file added assets/README-bec1ab77.png
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Binary file added assets/README-fcbec9a1.png
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
2 changes: 1 addition & 1 deletion prob2plot.m
Expand Up @@ -11,7 +11,7 @@ for i = 2:length(t)
end
a = figure(1);
y_heun = Prob2ptC (@(t,y) prob2_ode(t,y), 0.1, 1, [0 3]);
plot(t,exp(-t),'k-',t,y,'Xg',t,y_heun,'Or');
plot(t,exp(-t),t,y,t,y_heun);
legend('Analytical','Euler','Heun','Location','Northeast')
xlabel('Time')
ylabel('Position')
12 changes: 12 additions & 0 deletions prob2ptB.asv
@@ -0,0 +1,12 @@
dt=[0.5 .1 .001];
figure();
hold on
x=[0:dt(2):3]';
y_n=zeros(size(x));
y_n(1)=1;
for i=2:length(x)
dy=prob2_ode(x,y_n(i-1));
y_n(i)=y_n(i-1)+ dt(2)*dy;
end
plot(x,exp(-x),x,y_n)
legend('Analytical','Euler','Location','Northeast')
3 changes: 2 additions & 1 deletion prob2ptB.m
Expand Up @@ -8,4 +8,5 @@ for i=2:length(x)
dy=prob2_ode(x,y_n(i-1));
y_n(i)=y_n(i-1)+ dt(2)*dy;
end
plot(x,exp(-x),'k-',x,y_n,'o');
plot(x,exp(-x),x,y_n)
legend('Analytical','Euler','Location','Northeast')
17 changes: 17 additions & 0 deletions prob3plot.asv
@@ -0,0 +1,17 @@
set(0, 'defaultAxesFontSize', 16)
set(0,'defaultTextFontSize',14)
set(0,'defaultLineLineWidth',1)
dt = 0.01;
t = [0:dt:3];
y = zeros(length(t),2);
y(1,:) = [1,0];
for i = 2:length(t)
dy = prob3_ode(t,y(i-1,:));
y(i,:) = y(i-1,:)+dt*dy;
end
b = figure(2);
y_heun = prob3ptC (@(t,y) prob3_ode(t,y),0.01,[1,0],[0 3]);
plot(t,cos(3*t),t, y(:,1),t,y_heun(:,1))
legend('Analytical','Euler','Heun','Location','Northeast')
xlabel('Time')
ylabel('Position')
17 changes: 17 additions & 0 deletions prob3plot.m
@@ -0,0 +1,17 @@
set(0, 'defaultAxesFontSize', 16)
set(0,'defaultTextFontSize',14)
set(0,'defaultLineLineWidth',1)
dt = 0.01;
t = [0:dt:3];
y = zeros(length(t),2);
y(1,:) = [1,0];
for i = 2:length(t)
dy = prob3_ode(t,y(i-1,:));
y(i,:) = y(i-1,:)+dt*dy;
end
b = figure(2);
y_heun = prob3ptC (@(t,y) prob3_ode(t,y),0.01,[1,0],[0 3]);
plot(t,cos(3*t),t, y(:,1),t,y_heun(:,1))
legend('Analytical','Euler','Heun','Location','Northeast')
xlabel('Time')
ylabel('Position')
3 changes: 2 additions & 1 deletion prob3ptB.m
Expand Up @@ -8,4 +8,5 @@ for i=2:length(x)
dy=prob3_ode(x,y_n(i-1,:));
y_n(i,:)=y_n(i-1,:)+ dt(2)*dy;
end
plot(x,cos(3*x),'k-',x,y_n(:,1),'o');
plot(x,cos(3*x),x,y_n(:,1))
legend('Analytical','Euler','Location','Northeast')
20 changes: 20 additions & 0 deletions prob3ptC.m
@@ -0,0 +1,20 @@
function y = prob3ptC(prob3_ode, dt, y0,tlimits)
t = [tlimits(1):dt:tlimits(2)];
y = zeros(length(t),length(y0));
y(1,:) = y0;
for i = 2:length(t)
dy = prob3_ode(t(i-1),y(i-1,:));
y(i,:) = y(i-1,:)+dt*dy;
dyc = (dy + prob3_ode(t(i),y(i,:)))/2;
y(i,:) = y(i-1,:)+dt*dyc;
n = 1;
while (1)
n = n+1;
yold = y(i,:);
dyc = (dy+prob3_ode(t(i),y(i,:)))/2;
y(i,:) = y(i-1,:)+dt*dyc;
ea = abs(y(i,:)- yold)./y(i,:)*100;
if max(ea) < 0.0001 || n>100; break, end
end
end
end

0 comments on commit f97cd0a

Please sign in to comment.