Contents

function matlab_intro_part_2_assignment

MATLAB: Introduction - Part 2: Assignment

clear all
close all
clc

Problem 1

y = randn(2,3,40,12);

% Part (a)
ya = y(:,:,1,1);
size(ya)

% Part (b)
yb = y(:,:,1:end/2,5);
size(yb)

% Part (c)
y(1,1:2,1:10,3:5) = exp(1);

% Part (d)
try
    y(2,:,:,12) = ones(3,40);
    disp('Problem 1 (d): No exception caught');
catch EX
    disp(['Problem 1 (d): Exception caught: ' EX.message]);
end

% Part (e)
try
    y(:,:,15,1) = eye(3);
    disp('Problem 1 (e): No exception caught');
catch EX
    disp(['Problem 1 (e): Exception caught: ' EX.message]);
end
ans =

     2     3


ans =

     2     3    20

Problem 1 (d): No exception caught
Problem 1 (e): Exception caught: Subscripted assignment dimension mismatch.

Problem 2

N = 50;
s = 0;
for i = 1:N
    s = s + 1/i + 1/((i + 2)*(i + 3));
end
sprintf('Sum = %g',s)
ans =

Sum = 4.81367

Problem 3

n = 8;
H = zeros(8); % Preallocate H

for i = 1:n
    for j = 1:n
        H(i,j) = 1/(i + j - 1);
    end
end
H
H =

  Columns 1 through 7

    1.0000    0.5000    0.3333    0.2500    0.2000    0.1667    0.1429
    0.5000    0.3333    0.2500    0.2000    0.1667    0.1429    0.1250
    0.3333    0.2500    0.2000    0.1667    0.1429    0.1250    0.1111
    0.2500    0.2000    0.1667    0.1429    0.1250    0.1111    0.1000
    0.2000    0.1667    0.1429    0.1250    0.1111    0.1000    0.0909
    0.1667    0.1429    0.1250    0.1111    0.1000    0.0909    0.0833
    0.1429    0.1250    0.1111    0.1000    0.0909    0.0833    0.0769
    0.1250    0.1111    0.1000    0.0909    0.0833    0.0769    0.0714

  Column 8

    0.1250
    0.1111
    0.1000
    0.0909
    0.0833
    0.0769
    0.0714
    0.0667

Problem 4

x = linspace(-3,3,200);
y = zeros(size(x)); % Preallocate y

for i = 1:length(y)
    y(i) = prob4f(x(i));
end

figure('Name','Problem 4','Color','white');
plot(x,y);
xlabel('x');
ylabel('f(x)');

Problem 5

u = cmu.units;

% Methane component
methane.Pc = 46.4*u.bar;
methane.Vc = 99E-6*u.m^3/u.mol;
methane.a = 3*methane.Pc*methane.Vc^2;
methane.b = methane.Vc/3;

% Ranges for T and V
Tmin = 50*u.K;
Tmax = 600*u.K;
Vmin = 10E-4*u.m^3/u.mol;
Vmax = 200E-4*u.m^3/u.mol;

% Calculate P based on a meshgrid of T and V
T = linspace(Tmin,Tmax,300).';
V = linspace(Vmin,Vmax,300).';
[TT VV] = meshgrid(T,V);
P = vdWP(methane,VV,TT);

% Plot P with contours
figure('Name','Problem 5: Surface Plot','Color','white');
meshc(TT,VV,P/u.bar);
title('PVT Diagram for Methane');
xlabel('T [K]');
ylabel('V [m^3/mol]');
zlabel('P [bar]');

% Isotherms
Tiso = [100; 190; 250; 400]*u.K;
Piso = zeros(length(V),length(Tiso))*u.bar;

% Calculate P for each Tiso and plot the isotherm
figure('Name','Problem 5: Isotherms Plot','Color','white');
legendcell = cell(size(Tiso));
colors = lines(length(Tiso));
for i = 1:size(Piso,2)
    Piso(:,i) = vdWP(methane,V,Tiso(i));

    plot(V,Piso(:,i)/u.bar,'Color',colors(i,:));
    hold on;
    legendcell{i} = sprintf('T = %g',Tiso(i));
end
hold off;

title('Isotherms for Methane');
xlabel('V [m^3/mol]');
ylabel('P [bar]');
legend(legendcell{:});


%%%%%%%%%%%%%%%%%%%%%%%%%%%%% Subfunctions %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function fx = prob4f(x)
if (x < 0)
    fx = 0;
elseif (x >= 0 && x <= 1)
    fx = x;
elseif (x >= 1 && x <= 2)
    fx = 2 - x;
else
    fx = 0;
end


function P = vdWP(component,V,T)
u = cmu.units;
R = 8.314472E-5*u.m^3*u.bar/(u.K*u.mol);

P = R*T./(V - component.b) - component.a./V.^2;