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 );

w := proc (t) options operator, arrow; 800/(1+3*exp(-1/30*t)) end proc

p := proc (t) options operator, arrow; 13/20-1/100*t end proc

Profit := proc (t) options operator, arrow; w(t)*p(t)-9/20*t end proc

[Maple Plot]

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 );

dProfitdt := 80/(1+3*exp(-1/30*t))^2*(13/20-1/100*t)*exp(-1/30*t)-8/(1+3*exp(-1/30*t))-9/20

optTime := -30*RootOf(480*_Z*exp(_Z)+506*exp(_Z)-169-81*exp(_Z)^2)

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(%);

-30*RootOf(480*_Z*exp(_Z)+506*exp(_Z)-169-81*exp(_Z)^2,-.4111630481)

12.33489144

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 );

optTime := 12.33489144

Conclusion:  the pig should be sold after about 12 days in order to maximize profit.

>    plot( dProfitdt, t=0..20, color=black );

[Maple Plot]

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 );

p := proc (t) options operator, arrow; .65-i*t end proc

dProfitdt := 80/(1+3*exp(-1/30*t))^2*(.65-i*t)*exp(-1/30*t)-800/(1+3*exp(-1/30*t))*i-9/20

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 );

34.97964105

27.84522706

21.82172741

16.70625430

12.33489144

8.575274182

5.320746522

2.485414637

-0.

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"] );

[Maple Plot]

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 );

points := {[1/125, 21.82172741], [11/1000, 8.575274182], [13/1000, 2.485414637], [7/500, -0.], [3/250, 5.320746522], [3/500, 34.97964105], [7/1000, 27.84522706], [9/1000, 16.70625430], [1/100, 12.33489...
points := {[1/125, 21.82172741], [11/1000, 8.575274182], [13/1000, 2.485414637], [7/500, -0.], [3/250, 5.320746522], [3/500, 34.97964105], [7/1000, 27.84522706], [9/1000, 16.70625430], [1/100, 12.33489...

[Maple Plot]

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;

i1 := .1e-1

i2 := .101e-1

optTime1 := 12.33489144

optTime2 := 11.93334490

dOptTimedi := -4015.465400

sens_optTime_i := -3.255371496

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;

-3.255371496