Thermodynamics Assignment Help Professional

Calculating Derivatives

clc, clear all, close all

%% Part 1

%% Question 1

syms q p
C(q) = 610 - 0.2*q + 0.15*q^2; % Cost function
D(p) = 260 - 20*p; % Demand function

% Revenue is equal to the amount of units multiplied by the demand (price)
% function

R(p) = p*D(p)

% The profit is equal to the Revenue minus the Cost

P(p) = R(p) - subs(C(q), q, D(p))

% The marginal cost is the first derivative of the cost function
MC(q) = diff(C, q)

% The average cost is the Cost function divided by the total units
AC(q) = C(q)/q

% The marginal revenue is the first derivative of the Revenue function
MR(p) = diff(R,p);

%% Question 2
figure
ps = 0:0.1:13;
plot(ps, D(ps), 'linewidth', 2);
xlabel('Price ($)');
ylabel('Demand (q)');
title('Demand vs. Price');
grid on

%% Question 3
figure
subplot(1,2,1)
qs = 1:1:D(0);
plot(qs, C(qs), 'linewidth', 2), grid on
xlabel('Quantity (q)');
ylabel('Cost ($)');
title('Cost, Marginal Cost, Marginal Revenue and Average Cost vs. Quantity');

hold on
plot(qs, MC(qs), 'linewidth', 2);
plot(qs, AC(qs), 'linewidth', 2);
plot(qs, MR(qs), 'linewidth', 2);
legend('Cost', 'Marginal Cost', 'Average Cost', 'Marginal Revenue');

subplot(1,2,2)
plot(qs, MC(qs), 'linewidth', 2), hold on;
plot(qs, AC(qs), 'linewidth', 2);
plot(qs, MR(qs), 'linewidth', 2);
legend('Marginal Cost', 'Average Cost', 'Marginal Revenue');
xlabel('Quantity (q)');
ylabel('Cost ($)');
title('Marginal Cost, Average Cost and Marginal Revenue vs. Quantity');
grid on
%% Question 4
% To calculate the point of elasticity of demand, we select 2 points of the
% curve

% For example, we selec the points (8, 100) and (9, 80)

% Change in quantity
q_change_perc = (9-8)/(9+8)/2 *100;
p_change_perc = (80-100)/(80+100)/2 *100;

p_elast_D = abs(q_change_perc/p_change_perc)

% Because the value is < 1, the demand is inelastic

% Now, we will calculate the intervals from P = 1 to 13, where demand is
% elastic and inelastic

pvec = 0:0.01:13;
N = length(pvec);
elast_vec = [];
for i = 1:N-1
    P1 = [pvec(i) D(pvec(i))];
    P2 = [pvec(i+1), D(pvec(i+1))];
    q_change_perc = (P2(1)-P1(1))/(P2(1)+P1(1))/2 *100;
    p_change_perc = (P2(2)-P1(2))/(P2(2)+P1(2))/2 *100;

    p_elast_D = abs(q_change_perc/p_change_perc);
    elast_vec = [elast_vec, p_elast_D];
end

temp = pvec(find(elast_vec > 1));
elastic_interval = [temp(1), temp(end)]

temp = pvec(find(elast_vec < 1));
inelastic_interval = [temp(1), temp(end)]

temp = pvec(find(elast_vec == 1));
if length(temp)>0
    unitary_interval = [temp(1), temp(end)]
else
    unitary_interval = []
end
%% Question 6

% The revenue is a negative parabola, so we find its vertex
a = -20;
b = 260;
c = 0;

pmax = -b/(2*a);
qmax = subs(D,p,pmax);

%% Question 7

% The Average Cost is a function that is inversely proportional to the
% quantity

% To minimize this function, we will use the MATLAB's function fmincon

% We set the bounds for the quantity to: 0 < q < Inf
type ObjectiveAverageCost
[qmin, fval, exiflag] = fmincon('ObjectiveAverageCost', 10, [], [], [], [], 0, Inf);

qmin_int = round(qmin)

pmin = double(solve(D==qmin))
pmin_int = double(solve(D==qmin_int))

%% Question 8

% The profit function is a negative parabola, so the price that maximizes
% it is its vertex

a = -80;
b = 1816;

p_profit_max = -b/(2*a)
q_profit_max = D(p_profit_max)


%% Part 3

%% Question 10

% This time, the supply function is a price function, where the independent
% variable is the quantity.

% we can express the quantity in function of price as:
S(p) = p-8;

% The equilibrium point is where the demand and supply are the same

p_equil = solve(D == S)
q_equil = D(p_equil)

figure
ps = 0:0.1:20;
plot(ps, D(ps), 'linewidth',2), hold on
plot(ps, S(ps), 'linewidth', 2)
plot(p_equil, D(p_equil), 'r*', 'linewidth', 2);
xlabel('Price ($)');
ylabel('Quantity (q)');
title('Supply and Demand vs. Price');
legend('Demand', 'Supply');
grid on

%% Question 11

% The consumer surplus is the integral between the equilibrium point and
% the Demand function

CS = int((D(p)-q_equil), p, 0, p_equil)

%% Question 12

% The producer surplus is the integral between the Supply function and the
% equilibrium point
PS = int((q_equil - S(p)), p, 0, p_equil)

function f = ObjectiveAverageCost(x)

    f = 3*x/20 + 610/x - 1/5;
    
end

Solving Equations Using Matlab

clc, clear all, close all

%% Question 1
syms x y z

%% Question 2
fprintf("\n\nQuestion 2:\n");
% let's solve using symbolic solver
xsol = double(solve(x^4 + 2*x^3 - 8*x^2 - 9*x+18==0))

%% Question 3
fprintf("\n\nQuestion 3:\n");
xsol = double(solve(x^3 - 5*x^2 == -x-15))

%% Question 4
fprintf("\n\nQuestion 4:\n");
% The intersection points are the solutions of the equation
sol = solve([4*x^2+5*(y-1.5)^2 == 40;y-3*x+1==0], [x,y]);
xsol = double(sol.x)
ysol = double(sol.y)

%% Question 5
fprintf("\n\nQuestion 5:\n");
% Define the system of equations
equations = [2.1*x+6.2*y-3.1*z==205;
             -3.7*x+10.8*y+1.1*z == -107;
             x+2*y-3*z==23];
sol = solve(equations, [x,y,z]);
xsol = double(sol.x)
ysol = double(sol.y)
zsol = double(sol.z)

%% Question 6
fprintf("\n\nQuestion 6:\n");
% We will use x = F1 and y = F2
equations = [(1 + cosd(30))*x + 10.3*cosd(25)*y == 0;
             5*sind(30)*x - 5.6*sind(25)*y == 30];
         
sol = solve(equations, [x,y]);
F1 = double(sol.x)
F2 = double(sol.y)

%% Question 7
fprintf("\n\nQuestion 7:\n");
y(x) = tan(4-3*x)-sqrt(3);
% Now solve
xsol = double(solve(y(x)+2*x==7, x))

%% Question 8
fprintf("\n\nQuestion 8:\n");
clear y
syms y
z1(x,y) = 4*x^2+5*(y-1.5)^2 -40;
z2(x,y) = y-3*x+1;

% The intersection points are where z1== z2
sol = solve(z1==z2, [x,y]);
xsol = double(sol.x)
ysol = double(sol.y)

%% Question 9
fprintf("\n\nQuestion 9:\n");
% We will use x and y as our variables

% Define first equation
eq1 = x+y == 11; % the sum of both numbers is 11
eq2 = 10*y+x + 27 == 10*x+y; 
sol = solve([eq1;eq2], [x,y]);
number = double(10*sol.x + sol.y)

%% Question 10
fprintf("\n\nQuestion 10:\n");
d = 250; %from A to B
va = 60;
vb = 0;
ab = -20;
aa = 0;

% Equation for position of a:

eq1 = va*x + aa*x^2 /2;
eq2 = d + vb*x + ab*x^2 /2;

sol = solve(eq1 == eq2)
t = double(sol);
t = t(t > 0); % take only the positive value
Xa = double(va*t + aa*t^2 /2)
Xb = double(vb*t + ab*t^2 /2)

fprintf("The trains met when they are at at distance of %.3f of station A\n", Xa)

Thermodynamics Professional

Moffat