32chsCk.html

Partial Differential Equations PowerTool

by Dr. Jim Herod

Section 7.4: Recipe for a Cheese Cake

Maple Packages for Section 7.4

> restart:

> with(plots):

Warning, the name changecoords has been redefined

>

Here is a recipe for a cheese cake. If you have to worry about cholesterol, just stay away from this one. If you don't, this one is really great. Fix it a day before you want to impress guests and serve it with the best coffee you have. Wow! Guaranteed a success if you follow the directions.

I cook mine in a ten inch pan. It's about 1.5 inches thick. The color of the real cake has a range from brown, to orange, to yellow.

> with(plots):

> J:=plot3d([r,theta,1.5],r=0..5,theta=-Pi..Pi,coords=cylindrical,
scaling=constrained,color=yellow,orientation=[45,65],style=hidden):

> K:=plot3d([5,theta,z],theta=-Pi..Pi,z=0..1.5,coords=cylindrical,
scaling=constrained,color=orange,style=hidden):

> display3d({J,K});

>

Here is the recipe.

Ingredients: Herod's Cheese Cake

2 1/12 lb. cream cheese ( five 8 oz. pkg. )

1 3/4 cup sugar

3 t. flour

1 1/2 t. grated orange rind

1 1/2 t. grated lemon rind

1/4 t. vanilla

5 eggs

2 egg yolks

1/4 c. heavy cream.

Put cheese in mixer bowl. Beat at low speed. Add sugar gradually, then the remainder of ingredients in order . (Eggs should be added one at a time.) When blended and smooth, pour in a lined pan and place in a pre-heated 550 degree oven. Bake 12 - 15 minutes. Reduce heat to 200 degrees. ( Cool oven quickly by leaving the door open and fanning as necessary.) Continue baking for 1 hour. Cool before cutting. The cake is better after refrigeration. Lasts indefinitely ... if you can resist it ... in the refrigerator.

The Mathematical Question

Now, as an applied mathematician, here is your question: Eggs congeal at about 140 degrees. All the ingredients of the cheese cake before cooking are about 46 degrees. When the ingredients are mixed and placed in the hot oven, how long does it take to get the center of the cooking cake to 140 degrees?

The Mathematical Model.

This is a heat diffusion problem. We model it as

PDE: = c u

side boundary: u(t, 5, , z) = 550, t > 0, , 0 < z < 1.5,

top and bottom boundary: u(t, r, , 0) = 550 = u(t, r, , 1.5),

initial condition: u(0, r, , z) = 46.

That is, we assume the heat diffuses into the cheese cake as modeled by the diffusion equation for a cylinder, that the cylinder is initially at 46 degrees, and the surrounding temperature is 550 degrees.

We ask: at what time is the temperature in the middle 140 degrees? That is, compute t so that

u(t, 0, 0, 3/4 ) = 140.

Solution:

Recall that when the Laplacian operator is written out for a cylinder the equation is

+ = .

The number c above is the rate of diffusion for heat through the cream cheese/egg mixture. For this problem, since the initial condition and boundary condition are independent of theta, the solution will be independent of theta. Thus, we can rewrite the partial differential equation as

= + .

Agree that the steady state solution for the equation will be 550 degrees. (As a cook, you must not let the cake achieve this state!) Thus, we solve the problem with homogeneous boundary conditions and add 550. For ease of typing, take a = 5, the radius of the cake, and b = 1.5, the height of the cake.

Perform a separation of variables on u. That is, assume that u(t, r, z) = T(t) R(r) Z(z). We get that

(r R ') ' Z T + Z '' R T = T ' R Z.

Dividing by R Z T gives

(r R ') ' / R + Z '' / Z = T ' / T.

As usual, we make the argument that each of the terms of the sum on the left side of the above equation must be constant. We have that

(1) (r R ') ' = R , R(a) = 0.

(2) Z '' = Z, Z(0) = Z(b) = 0.

and

(3) T ' = - c ( ) T.

We handle these one at a time.

Analysis of equation (1):

We rewrite equation (1) as

(r R ') ' = r R , R(a) = 0.

R '' + r R ' + R = 0, with R(a) = 0.

We verify that the solution for this last equation is among the Bessel functions.

> R:=r->BesselJ(0,mu*r);
r^2*diff(R(r),r,r)+r*diff(R(r),r)+mu^2*r^2*R(r);
simplify(%);

Analysis of equation (2):

Equation (2) is recognized as an equation, with boundary conditions, which defines the sine functions.

> Z:=z->sin(lambda*z);
diff(Z(z),z,z)+lambda^2*Z(z);

Analysis of equation (3):

Equation (3) is recognized as an equation for the exponential function.

> T:=t->exp(-c*(mu^2+lambda^2)*t);
diff(T(t),t)+c*(mu^2+lambda^2)*T(t);

Products of solutions of (1), (2), and (3).

Products of solutions for equations (1), (2), and (3) should form a solution for the original equation.

> u:=(t,r,z)->BesselJ(0,mu*r)*sin(lambda*z)*exp(-c*(mu^2+lambda^2)*t);

> 1/r*diff(r*diff(u(t,r,z),r),r)+diff(u(t,r,z),z,z)-1/c*diff(u(t,r,z),t);

> simplify(%);

>

General solution

We add constants times products of solutions for (1), (2), and (3) to make the general solution of the original equation:

u(t, r, z) =

We verify that this sum is a solution. For simplicity, take only 9 terms.

> u:=(t,r,z)->sum(sum(A[n,m]*BesselJ(0,mu[m]*r)*sin(lambda[n]*z)*
exp(-c*(mu[m]^2+lambda[n]^2)*t),m = 1 .. 3),n = 1 .. 3);

> 1/r*diff(r*diff(u(t,r,z),r),r)+diff(u(t,r,z),z,z)-1/c*diff(u(t,r,z),t):

> simplify(%);

>

Computing 's and 's.

We identify the eigenvalues and ,

For , we have that sin( b) = 0, so that b = n .

> for n from 1 to 40 do
lambda[n]:=n*Pi/(3/2):
end do:
n:='n':

>

For , we have that BesselJ(0, a) = 0, so that is the zero of BesselJ(0,x) divided by a.

> for m from 1 to 40 do
mu[m]:=evalf(BesselJZeros(0,m))/5:
end do:
m:='m':

Computing coefficients

We choose the coefficients so that when t = 0,

46 - 550 = .

These coefficients will come from the Fourier coefficients:

= .

To speed the calculations, we do all the computations for n and then all the computations for m and multiply these. Here are the computations for the sine terms, for the terms involving n.

> a:=5; b:=1.5;

> for n from 1 to 40 do
T[n]:=int(sin(lambda[n]*z),z = 0 .. 3/2)/
int(sin(lambda[n]*z)^2,z = 0 .. 3/2):
end do:
n:='n':

>

Here are the computation for the Bessel terms, for the terms involving m.

> for m from 1 to 30 do
B[m]:=int(BesselJ(0,mu[m]*r)*r,r = 0 .. 5)/
int(BesselJ(0,mu[m]*r)^2*r,r = 0 .. 5):
end do:
m:='m':

>

Here is multiplying these together to get .

>    for n from 1 to 40 do
  for m from 1 to 30 do
     A[n,m]:=evalf((46-550)*T[n]*B[m]):
  end do:
end do:
n:='n':m:='m':

>

From here on, the computations take some time and space. Some how, this is not a surprise: if you eat much of this cheese cake, you will take a lot of space! On my little home computer, the total time to execute the entire worksheet with Maple 8 was about 100 seconds. The space used was 15.6 Bytes.

A graph for an approximation of the initial value

> u0:=(r,z)->sum(sum(A[n,m]*BesselJ(0,mu[m]*r)*sin(lambda[n]*z),n=1..40),m=1..30);

> plot(u0(r,3/4)+550,r=-5..5, numpoints=200);

>

Here is the definition for the solution. This is needed in what follows.

> u:=(t,r,z)->sum(sum(A[n,m]*BesselJ(0,mu[m]*r)*sin(lambda[n]*z)* exp(-c*(mu[m]^2+lambda[n]^2)*t),n=1..40),m=1..30);

>

We are ready to find when the middle of the cake achieves 140. To do this, you need to know c for this cheese cake. It is

> c:=0.0077;

A graphical method for finding when the center is 140 degrees.

To find the time for the center -- u(t, 0, 3/4) -- to be 140, we could draw make a graphical approximation.

> plot([140,u(t,0,.75)+550],t=12..13);

It seems clear that 140 degrees is achieved between 12.9 and 13 minutes.

> plot([140,u(t,0,.75)+550],t=12.9..13);

Now, I surmise that the desired t solution is between 12.96 and 12.98.

> plot([140,u(t,0,.75)+550],t=12.96..12.98);

>

That defines the value of t pretty closely. The value is about 13 minutes.

A numerical procedure for finding t.

To find the time for the center to be 140, we could ask Maple to solve the following equation for t.

> Digits:=4;
fsolve(u(t,0,3/4)+550=140,t,12.9..13);

> Digits:=10;

>

EMAIL: herod@math.gatech.edu or jherod@tds.net

URL: http://www.math.gatech.edu/~herod

(The receipe for the cheese cake is not included in the copyright. I have no idea where my grandmother got it.)