% Electric Motor Modeling Example
% Carnegie Mellon University, 24-370: Engineering Design I
% Professor Steve Collins
% Spring, 2011

% initial conditions (state is [current, angle, angular velocity]')
initial_state = [0; 0; 0];
% time range of simulation in seconds
timespan = 0:2e-5:1;  
% set some options for the forward-integration
ode_options = odeset('RelTol',1e-7,'AbsTol',1e-6);
% forward-integrate the equations of motion
[t,x] = ode45(@motorodefun,timespan,initial_state); % t is time, x is state
% rename for bookkeeping
I = x(:,1);
theta = x(:,2);
dtheta = x(:,3);

% Reconstruct state-dependent results
T = zeros(size(t));  V = zeros(size(t));
for i = 1:length(t)
    V(i) = voltagefun(t(i),x(i,:)');
    T(i) = torquefun(t(i),x(i,:)');
end

% Calculate power measures
% electrical power
Pe =    % electrical power is the product of current and voltage
% mechanical power
Pm =    % mechanical power is the product of 

% Plot the results
% motor current and voltage
figure(1); plot(t,I); xlabel('time'); ylabel('current (Amperes)');
figure(2); plot(t,V); xlabel('time'); ylabel('voltage (Volts)');
% motor kinematics and torque
figure(3); plot(t,theta); xlabel('time'); ylabel('motor angle (radians)');
figure(4); plot(t,dtheta); xlabel('time'); ylabel('motor angular velocity (radians per second)');
figure(5); plot(t,T); xlabel('time'); ylabel('motor torque (Newton meters)');
% power
figure(6); plot(t,[Pe Pm]); xlabel('time'); ylabel('power (Watts)'); legend('Electrical Power','Mechanical Power')
