% PC profit example
%
% This code solves the PC profit optimization problem we covered in class,
% and calculates a few sensitivity coefficients.
%
% You can either run the code all at once, using the 'run' command, or in
% sections, using 'run section'. Sections are created using '%%'.
                                                
% Note: if you ever want to start the file over, type 'clear' to delete the
% workspace. 'clc' clears the command window.

%% Define variables and parameters/constants

% Note: the semicolon at the end of each line suppresses output.

n = 10000;  % Current sales [units/month]
m = 700;    % Manufactoring cost [$/unit]
w = 950;    % Wholesale unit price [$/unit]
d = 100;    % Amount of test discount [$/unit]
r = 0.5;    % Percent sales increase based on test discount [-]
a = 50000;  % Current advertizing budget [$/month]
b = 10000;  % Amount of test ad budget increase [$/month]
c = 200;    % Increase in sales due to ad budget test [units/month]
A = 100000; % Cap on ad budget [$/month]

%% Plot the model

D = linspace(-10,10,200);      % Create a test vector for D, since we haven't specified it. 
B = linspace(-1000,1000,200);
[Dmesh, Bmesh] = meshgrid(D,B);
N = n*(1 + r*Dmesh) +c*Bmesh;
q = w - m - d*Dmesh;
P = N.*q - a - b*Bmesh;       % The .* operation indicates vector multiplication. Matlab
                    % multiplies these vectors element-by-element.


mesh(D,B,P)

%% Symbolic differentiation

% Create symbolic variables and equations
syms sr sD          % This initializes sr and sD as symbolic variables
sN = n*(1 + sr.*sD);
sq = w - m - d*sD;
sP = sN.*sq - a 

% Differentiate: d(total profit)/d(discount)
dPdD = diff(sP,sD)  % Since sP and sD are symbolic, this calls Matlab's 
                    % symbolic differentiation. See
                    % www.mathworks.com/help/symbolic/diff.html

%% Find optimal discount

% Set dPdD=0 and solve for sD
% This gives us the optimal discount as a function of sr.
[sD] = solve(dPdD,sD)

%% Plot the optimal discount for some reasonable r-values
Dstar = subs(sD,sr,r)       % Substitute r = 0.5 into our equation for optimal discount
rvec = 0.41:0.01:0.6;
Dvec = subs(sD,sr,rvec);
plot(rvec,Dvec)

%%  Sensitivity analysis

%% Sensitivity of Dstar to r
dDdr = diff(sD,sr);             % d(optimal discount)/dr
dDdr_vec = subs(dDdr,sr,rvec);  % Evaluate this derivative at rvec
SDr = dDdr_vec.*rvec./Dstar;    % Sensitivity function S(Dstar,r)
plot(rvec,SDr)                  % Plot sensitivity as a function of r

%% Sensistivity of total profit P to r
 dPdr = diff(sP,sr)             % d(total profit)/dr
 dPdr_vec = subs(-10000*sD*(100*sD - 250),sr,rvec);  % evaluate dP/dr at rvec 
 Nstar= n.*(1+Dstar.*rvec);     % Number sold, evaluated at rvec 
 qstar = w-m-d*Dstar;           % Profit/unit, evaluated at rvec
 Pstar = Nstar.*qstar - a;      % Total profit, evaluated at rvec
 SPr = dPdr_vec.*rvec./Pstar;   % Sensitivity function S(total profit, r)
 plot(rvec,SPr)