work in progress solution for exercise 1 of coursework

6 files changed, 138 insertions, 0 deletions

diff --git a/coursework17/error_script.m b/coursework17/error_script.m
new file mode 100644
index 0000000..29988bf
--- /dev/null
+++ b/coursework17/error_script.m
@@ -0,0 +1,13 @@
+% This script will carry out error analysis
+T = 150e-6;
+f = 1/T;
+
+A=6;
+
+data_points=10e3;
+
+% TODO find exact solution
+exact = @(t) ();
+
+for j=1:data_points
+	my_error = exact(time_array(j)) - Vout_array(j);
diff --git a/coursework17/heun.m b/coursework17/heun.m
new file mode 100644
index 0000000..ab7edfd
--- /dev/null
+++ b/coursework17/heun.m
@@ -0,0 +1,25 @@
+function [time_array, Vout_array] = heun(R, L, Vin, current_initial, step, t_final)
+	% Determine how many times to step forward
+	N = round(t_final/step);
+
+	% Initialise some arrays for the output
+	current_array=zeros(1,N);
+	time_array=zeros(1,N);
+	current_array(1) = current_initial;
+	Vout_array(1) = Vin(0) - R*current_initial;
+
+	% Given the input voltage as a function of time we can express 
+	% dI/dt (dcurrent_dt) as a function of time as well (given R and L)
+	dcurrent_dt = @(t, I) (Vin(t) - I*R)/L;
+
+	% Loop as many times as step allows and:
+	% calculate the gradient at point
+	% Use it to emtimate next two points
+
+	for j=1:N-1 % loop N-1 times
+		grad = feval(dcurrent_dt, time_array(j), current_array(j));
+		current_array(j+1)=current_array(j) + step/2*(grad + feval(dcurrent_dt, time_array(j)+step, current_array(j)+step*grad)) ; 
+		time_array(j+1) = time_array(j) + step;
+		Vout_array(j+1) = Vin(time_array(j+1)) - R*current_array(j+1);
+	end
+
diff --git a/coursework17/midpoint.m b/coursework17/midpoint.m
new file mode 100644
index 0000000..44ccdae
--- /dev/null
+++ b/coursework17/midpoint.m
@@ -0,0 +1,25 @@
+function [time_array, Vout_array] = midpoint(R, L, Vin, current_initial, step, t_final)
+	% Determine how many times to step forward
+	N = round(t_final/step);
+
+	% Initialise some arrays for the output
+	current_array=zeros(1,N);
+	time_array=zeros(1,N);
+	current_array(1) = current_initial;
+	Vout_array(1) = Vin(0) - R*current_initial;
+
+	% Given the input voltage as a function of time we can express 
+	% dI/dt (dcurrent_dt) as a function of time as well (given R and L)
+	dcurrent_dt = @(t, I) (Vin(t) - I*R)/L;
+
+	% Loop as many times as step allows and:
+	% calculate the gradient at point
+	% Use it to emtimate next two points
+
+	for j=1:N-1 % loop N-1 times
+		grad = feval(dcurrent_dt, time_array(j), current_array(j));
+		current_array(j+1)=current_array(j) + step*feval(dcurrent_dt, time_array(j)+step/2, current_array(j)+step/2*grad); 
+		time_array(j+1) = time_array(j) + step;
+		Vout_array(j+1) = Vin(time_array(j+1)) - R*current_array(j+1);
+	end
+
diff --git a/coursework17/octave-workspace b/coursework17/octave-workspace
new file mode 100644
index 0000000..e69de29
--- /dev/null
+++ b/coursework17/octave-workspace
diff --git a/coursework17/ralston.m b/coursework17/ralston.m
new file mode 100644
index 0000000..d910442
--- /dev/null
+++ b/coursework17/ralston.m
@@ -0,0 +1,25 @@
+function [time_array, Vout_array] = ralston(R, L, Vin, current_initial, step, t_final)
+	% Determine how many times to step forward
+	N = round(t_final/step);
+
+	% Initialise some arrays for the output
+	current_array=zeros(1,N);
+	time_array=zeros(1,N);
+	current_array(1) = current_initial;
+	Vout_array(1) = Vin(0) - R*current_initial;
+
+	% Given the input voltage as a function of time we can express 
+	% dI/dt (dcurrent_dt) as a function of time as well (given R and L)
+	dcurrent_dt = @(t, I) (Vin(t) - I*R)/L;
+
+	% Loop as many times as step allows and:
+	% calculate the gradient at point
+	% Use it to emtimate next two points
+
+	for j=1:N-1 % loop N-1 times
+		grad = feval(dcurrent_dt, time_array(j), current_array(j));
+		current_array(j+1)=current_array(j) + step/4*(grad + 3*feval(dcurrent_dt, time_array(j)+2/3*step, current_array(j)+2/3*step*grad)) ; 
+		time_array(j+1) = time_array(j) + step;
+		Vout_array(j+1) = Vin(time_array(j+1)) - R*current_array(j+1);
+	end
+
diff --git a/coursework17/ralston_script.m b/coursework17/ralston_script.m
new file mode 100644
index 0000000..9db929f
--- /dev/null
+++ b/coursework17/ralston_script.m
@@ -0,0 +1,50 @@
+% Vin = aI' + bI
+
+R = 0.5;    % 0.5Ohm
+L = 0.0015; % 1.5mH
+data_points = 10000;
+
+% Go on for 8 time constants
+time_constant = L/R;
+step = time_constant*8/data_points;
+
+% Signal 1 Vin = 5.5V
+Vin = @(t) 5.5 + 0*t;  % 5.5V
+current_initial=0;
+
+[t, Vout] = ralston(R, L, Vin, current_initial, step, data_points*step);
+plot(t, Vout);
+
+% Signal 2a Vin = 5.5*exp(-t^2/r) V
+Vin = @(t) 5.5*exp(-(t*160*10^(-6))^2)  % 5.5V
+current_initial=0;
+
+[t, Vout] = ralston(R, L, Vin, current_initial, step, data_points*step);
+plot(t, Vout);
+
+% Signal 2b Vin = 5.5*exp(-t/r)
+Vin = @(t) 5.5*exp(-t*160*10^(-6))  % 5.5V
+current_initial=0;
+
+[t, Vout] = ralston(R, L, Vin, current_initial, step, data_points*step);
+plot(t, Vout);
+
+% Signal 3
+T(1) =  20e-6;  % 20us
+T(2) =  450e-6; % 450us
+T(3) =  1000e-6;% 1000us
+
+for j=1:3
+	f = 1/T(j);
+	Vin = @(t) sin(2*pi*f*t);
+	[t, Vout] = ralston(R, L, Vin, current_initial, step, data_points*step);
+	plot(t, Vout);
+
+	Vin = @(t) sawtooth(2*pi*f*t);
+	[t, Vout] = ralston(R, L, Vin, current_initial, step, data_points*step);
+	plot(t, Vout);
+
+	Vin = @(t) square(2*pi*f*t);
+	[t, Vout] = ralston(R, L, Vin, current_initial, step, data_points*step);
+	plot(t, Vout);
+end