The Pig Problem
t = time (days)
w(t) = weight of the pig (lbs) at time t
p(t) = price for pigs ($/lbs)
Profit(t) = profit realized from sale of pig at time t ($)
> | restart: |
Define the model, and plot Profit as a function of time t ...
> | w := t -> 800 / ( 1 + 3 * exp(-t/30) ); |
> | p := t -> 65/100 - 1/100 * t; |
> | Profit := t -> w(t) * p(t) - 45/100 * t; |
> | plot( Profit(t), t=0..30, color=black ); |
It looks like the optimal time to sell lies around 12 days.
Attempt to solve analytically ...
> | dProfitdt := diff( Profit(t), t ); |
> | optTime := solve( dProfitdt=0, t ); |
Note that the derivative of the Profit function involves exponentials.
Setting it to zero yields a trancendental equation, which may be difficult to solve analytically if at all.
Indeed, the solution returned by Maple looks like a mess.
RootOf is a placeholder representing all roots of its argument.
We can ask Maple to give us all the roots by using the allvalues command ...
> | allvalues(%); |
> | evalf(%); |
Note that Maple returned a floating point number ,
indicating that Maple really solved for this equation numerically
(using Newton's method or a more sophisticated solver).
We could have asked Maple to solve for the optimal time to sell the pig directly by using the fsolve command ...
> | optTime := fsolve( dProfitdt=0, t ); |
Conclusion: the pig should be sold after about 12 days in order to maximize profit.
> | plot( dProfitdt, t=0..20, color=black ); |
Sensitivity analysis for the new model
(examining the effect of changing the 0.01 number on the time to sell).
There is no point asking Maple to solve for the optimal time analytically.
It could not do it earlier, and the structure of the equations has not changed.
Try for a numerical solution instead ...
> | p := t -> 0.65 - i*t; |
> | dProfitdt := diff( Profit(t), t ); |
> | optTime := fsolve( dProfitdt=0, t ); |
Error, (in fsolve) i is in the equation, and is not solved for
This error message indicates that Maple can't find a numerical solution, because it does not know the value of i.
That is, we won't be able to get an explicit expression for the optimal time to sell the pig.
Can we still make a graph of the optimal time to sell as a function of i?
One option is to do it by brute force, that is, by selecting some values of i, and calculating the corresponding optimal times to sell.
> | solve( subs( i=0.006, dProfitdt )=0, t ); |
> | solve( subs( i=0.007, dProfitdt )=0, t ); |
> | solve( subs( i=0.008, dProfitdt )=0, t ); |
> | solve( subs( i=0.009, dProfitdt )=0, t ); |
> | solve( subs( i=0.01, dProfitdt )=0, t ); |
> | solve( subs( i=0.011, dProfitdt)=0, t ); |
> | solve( subs( i=0.012, dProfitdt)=0, t ); |
> | solve( subs( i=0.013, dProfitdt)=0, t ); |
> | solve( subs( i=0.014, dProfitdt)=0, t ); |
Now use the pointplot command to plot the points. The argument for pointplot is a sequence, in curly brackets, of coordinates to be plotted. Each of the coordinate is of the form [x,y].
> | with(plots): |
> | pointplot({[0.006,34.97964195],[0.007,27.84522706],[0.008,21.82172741],[0.009,16.70625430],[0.01,12.33489144],[0.011,8.575274182],[0.012,5.320746522],[0.013,2.485414637],[0.014,0]}, labels=["i","opt time to sell"] ); |
CONCLUSION: the faster the price of pigs drops, the earlier the pig should be sold, as we concluded before.
How we got the graph obviously is rather ridiculous.
It's way too much work to execute the same type of command a number of times, then collect the results by hand.
Here's an automated way of accomplish the same thing ...
> | points := { seq( [ ivalue/1000, solve( subs( i=ivalue/1000, dProfitdt )=0, t ) ], ivalue=6..14 ) }; |
> | pointplot( points ); |
Last but not least, the quantitative sensitivity analysis.
We won't be able to compute an exact sensitivity value, but can approximate one, as shown in class.
We begin by computing an approximate value for the derivative of optTime with respect to i,
then computing the resulting approximation to the sensitivity value ...
> | i1 := 0.01; i2:= i1*1.01; |
> | optTime1 := solve( subs(i=i1,dProfitdt)=0, t ); |
> | optTime2 := solve( subs(i=i2,dProfitdt)=0, t ); |
> | dOptTimedi := ( optTime2-optTime1 ) / ( i2-i1 ); |
> | sens_optTime_i := dOptTimedi * i1 / optTime1; |
If we're into approximations anyways, why not obtain the approximation directly.
We need to know the percent change in the optimal time to sell the pig given a 1% change in the rate at which the price of the pig changes ...
> | (optTime2-optTime1)/optTime1*100; |