The Perfect Brownie Pan

Riley Dulin Jun Seok Yeo Nick Fuerstenberg

March 2015

1 Abstract

This paper examines different shapes for a brownie pan such that it maximizes the area covered in an oven, and the temperature distribution across the brown- ies. We also examine the relationship between the perimeter of the pan and the temperature distribution. We used a combination of analytical and numerical solutions to find these values, and wrote MATLAB code to simulate the tem- perature distribution. Through these methods, we hope to discover the optimal shape for a brownie pan.

2 Introduction

When baking in a rectangular pan, heat is concentrated in the 4 corners and the product gets overcooked at the corners (and to a lesser extent at the edges). In a round pan the heat is distributed evenly over the entire outer edge and the product is not overcooked at the edges, and there are no corners at which the product will be extremely overcooked. However, since most ovens are rectangu- lar in shape using round pans is not efficient with respect to using the space in an oven. Develop a model to show the distribution of heat across the outer edge of a pan for pans of different shapes – rectangular to circular and other shapes in between.

Assume:

1. A width to length ratio of W/L for the oven which is rectangular in shape

2. Each pan must have an area of A

3. Pans in the oven must have a space d between each other and the sides of the oven

4. Initially two racks in the oven, evenly spaced

Develop a model that can be used to select the best type of pan (shape) under the following conditions:

1. Maximize number of pans that can fit in the oven (N)

2. Maximize even distribution of heat (H) for the pan

3. Optimize a combination of conditions (1) and (2) where weights p and (1−p) are assigned to illustrate how the results vary with different values of W/L and p

[1]

3 Assumptions

First The ideal shape of the brownie pan is a member of the family of Superel- lipses (also known as Lamé Curves)∣∣∣x

∣∣∣n + ∣∣∣y b

∣∣∣n = 1 (1) n ∈ [2,∞) (2)

Where a and b are the positive half-width and half-length respectively of the shape. Below is a picture of the ”squircle”, a superellipse with n = 4,a = 1,b = 1:

Figure 1: A ”Squircle”: x4 + y4 = 1

This is justified because the family of superellipses go from perfect ellipse to perfect rectangle (n = 10 is sufficiently close to infinity to make a near- perfect rectangle), and in between they form rounded rectangles, which would be the perfect compromise between the two extremes.

3.1 Area Assumptions

Second The optimal layout of the pans in the oven is in a rectangular spacing pattern

Third a and b have the same ratio to each other as the oven’s width to length ratio. This is justified because in order to maximize the rectangular spac- ing pattern, an equal amount of rows and columns of brownie pans are necessary because a square has a higher number of pans than rectangles, and in order to get a square pattern each pan must match the dimensions of the oven

Fourth The second rack is exactly the same in alignment and heat distribution as the first rack, so we will only analyze the situation for one rack. This is justified because an oven heats the top and bottom rack in approxi- mately the same way, and we assume the brownies on the top rack will not interfere with the radiation received by the bottom rack

3.2 Heat Assumptions

Fifth Heat flux only flows from high temperature to low temperature and does not flow back

Sixth The pan and the brownies are perfect black-bodies

Seventh The temperature and the radiation in the oven are equally distributed

Eigth The temperature of the pan is equally distributed

Ninth All energy that the pan and the brownies receive from radiation converts to heat with 100% transfer (ideal transfer rate)

Tenth Thermal conductivity is constant even when there is a change in temper- ature

Eleventh Oven remains the same temperature and provides unlimited heat until equilibrium

Twelfth The oven has a vacuum inside, and thus there are no convection currents, and only radiative heat transfer. We assume this because it makes calcu- lations significantly easier to make

3.3 Strategies

We had several strategies of how to model this problem. The problem was split into two halves, the maximization of the ”Area Score”, and the maximization of the ”Heat Score”. These are termed as scores because in order to compare things as different as area and temperature distributions, we needed to normalize the values. So we give each field a score between 0 and 1 by subtracting the minimum and dividing by the range of the area and temperature values, so that the minimum scoring shape is assigned 0, the maximum shape is assigned 1,

and the rest are scaled to be between those two numbers while still retaining the shape of their plot. We then use those values to compare against each other and decide what is the optimal shape.

During our research we found that an alternative exists for calculating the Heat Score, namely, the Perimeter Score. More details are found in the perime- ter section.

3.4 Constants and Variables

3.4.1 Constants

W Width of the oven (m) L Length of the oven (m) r Ratio of W/L A Area of the pan (m 2) d Distance between each pan (m) p Weight of Area Score (Heat Score is 1 −p) h Height of the brownie pan (m) pt Thickness of the brownie pan (m) Kb Thermal conductivity of brownies (W/(m*K)) cb Specific heat of brownies (J/(kg*K)) cp Specific heat of the pan (J/(kg*K)) σ Stefan-Boltzmann’s Constant (W/(m2 ∗K4)) ρb Density of the brownies (kg/m

3) ρp Density of the pan (kg/m

3) TR Room temperature (K) TO Oven temperature (K) ∆t Time step (s) t0 Start time (s) tf End time (s) ∆x Spatial step (m) ccube Specific heat for one cube (J/(kg*K)) Qtop Heat from the oven radiation (J/kg)

Qbottom Heat from the pan (J/kg) Qside,j Heat from sides of a cube (J/kg)

3.4.2 Variables

n Parameter of the superellipse a Half-width of the superellipse (m) b Half-length of the superellipse (m) x Spatial x coordinate (m) y Spatial y coordinate (m) t Time value (s) N Number of pans P Perimeter of the superellipse Tp Temperature of the pan (K) Tm Temperature for one cube on the mth second (K)

Tside,j Temperature of the jth side of a brownie (K) T(x,y,t) Temperature function (K)

SA(W,L,a,b,d) Area Score SH (a,b,n) Heat Score SP (a,b,n) Perimeter Score

4 Model

The formal solution to this problem is the value of n that maximizes the follow- ing expression

p∗SA(W,L,a,b,d) + (1 −p) ∗SH (a,b,n) (3)

where n is the parameter of a superellipse given by (1) and (2) and each S is normalized in the following way SA:

S = V −min(V)

max(V) −min(V) (4)

SH,SP :

S−1 = min(V) −V

max(V) −min(V) + 1 (5)

where V is a vector of values that are determined by each respective score and S is the normalized score vector. SH and SP use the ”inverse” normalization function because their values are opposite to the desired measurement. See Sections 4.1.3, 4.2.5, and 4.3.1 for individual explanations of why the scores have their formulae.

4.1 Area Score

To find the Area Score we followed this outline:

1. Find a and b

2. Calculate N, the number of pans

3. Normalize the values

4.1.1 Find a and b

To find a and b, we use a formula for finding the area of a superellipse:

A = 4ab∗ ( Γ ( 1 + 1

))2 Γ ( 1 + 2

) (6) a = br (7)

We then solve these equations, and thus we have found the a and b for this value of n. We will use this for all of the scores we compute.

4.1.2 Find N

To find N, we must form equations for how many pans can fit in the oven. We know a and b, along with the constants that represent the oven size, and we form these two linear equations:

(Ncol + 1) ∗d + Ncol ∗ 2a = W (8)

(Nrow + 1) ∗d + Nrow ∗ 2b = L (9) N = Nrow ∗Ncol (10)

4.1.3 Normalize

Once we finish calculating N for each n, we normalize them using (4), because the number of pans that can fit in the oven is proportional to our desired be- havior. In other words, as the number of pans increases, so should the score, so we use the first normalization function.

4.2 Heat Score

The Heat Score is the most complicated of all the scores to calculate, and requires the creation of a simulation of the cooking brownies. To simulate the event, we break the brownie mixture down into very small cubes and find the energy transfer through each cube. To accomplish this we wrote some code in MATLAB that simulated the transfer of heat across the brownie mixture. See Appendix A for more details on the simulation, and see Appendix B for our consideration of an analytical solution

4.2.1 Heat Transfer for the pan

As shown in the diagram, the pan accepts heat from the oven in the form of radiation according to Stefan–Boltzmann’s Law:

∂Qpan ∂t

= AσT 4O (11)

Qpan = cp ∗ ∆Tp (12)

∆Tp = AσT 4O cp

∆t (13)

Figure 2: Diagram of Heat Transfer

4.2.2 Heat Transfer for a Single Brownie Cube

Figure 3: Segment of brownie into cubes

The total energy delivered to a single brownie cube of width and length ∆x is given by the sum of the heat transferred by all of its neighboring brownies, or the pan if it is touching the pan, plus the heat transfer from the top of the oven and the heat transfer from the bottom of the pan. To determine if the cube is touching the pan, we use (1) with the x and y of the brownie and see if that is > 1. If so, then that side will be set to the Tp. The heat transfer equations from radiation (from the oven on the top of the cube) and conduction (from the pan on the bottom of the cube) are as follows:

4.2.3 Heat Transfer for a Single Brownie Cube

Figure 4: Diagram of Heat Transfer

Figure 4 shows how the heat is transferred onto one cube from multiple sources

Qtop = ∆t∆x 2σT 4O (14)

Qbottom = ∆t∆x 22kb

Tp −Tm h

(15)

Qside,j = ∆thkb (Tside,j −Tm) (16)

Since ccube = cbρb∆x

2h (17)

Q = ccube ∗ ∆T (18)

the equation for the change in temperature for one cube from the top and bottom is:

∆T1 = ∆t∆x2

ccube

( σT 4O + 2kb

Tp −Tm h

) (19)

and the equation for the change in temperature for one cube from its neighboring cubes is:

∆T2 = ∆tkbh

ccube

  4∑

j=1

(Tside,j) − 4Tm

  (20)

Combine them together to get the net change in temperature for one brownie cube:

Tm+1 = Tm + ∆T1 + ∆T2 (21)

T0 = TR (22)

We split the mixture into a 200×200 matrix of these cubes, and calculate this change in temperature for each cube for m = 1200, or 20 minutes (the average baking time for a brownie recipe). [2] [3]

4.2.4 Average Variance

In order to assign a particular number for a heat distribution to be described by, we decide to use the average standard deviation of the temperature over the whole time period. The average standard deviation is defined as the square root of the average variance:

σT =

√√√√√ tf∑m=0 σ2Tm tf

(23)

This number increases approximately linearly as the n value increases, then plateaus at around n > 8.

4.2.5 Normalize

Similarly to the other scores, we then normalize this score using (5), because an increase in average variance means that the heat is not as well distributed, which goes against our desired characteristic. Thus we use the ”inverse” normalization function to make low variances higher valued.

4.3 Perimeter Score

We developed a third score, the Perimeter score, as a response to the length of time it took to run the heat simulation. We assume that the strength of the Heat Score is dependent only on one factor, the perimeter of the superellipse. This is a logical assumption because brownies typically only burn at the edges and corners, and so if the area is fixed, then the minimal perimeter would have the best Heat Score. We will verify if this is true by comparing it to the Heat Score and seeing if it is a decent approximation. Finding the Perimeter Score is fairly straightforward, after finding a and b during the Area Score part, we simply need to find the arc length of the superellipse to get the perimeter:

4 ∗ ∫ a

√√√√ 1 +

( xn

an −xn

)2(n−1) n

dx (24)

We can solve this using numerical integration software.

4.3.1 Normalize

This perimeter value is then normalized using (5), because as perimeter in- creases, the heat is not as well distributed, which goes against our desired characteristic. Thus we use the ”inverse” normalization function to make low variances higher valued.

5 Verification

After running the simulation and obtaining our results, we find that the perfect shape of the brownie is given by n = 3.3. The shape of this pan and the plot of the scores are shown in Figures 6 and 5, respectively. Figure 7 shows the simulation of the heat distribution over time for n = 3.3 in a 3D plot and a 2D contour plot.

Figure 5: Area Scores(Blue), Heat Scores(Orange), and Perimeter Scores(Red)

Figure 6: The Shape of the Pan

(a) t = 8 (b) t = 8

(e) t = 56 (f) t = 56

(g) t = 80 (h) t = 80

Figure 7: Temperature Graphs (t in minutes)

5.1 Perimeter Perspective

As you can see in Figure 5, the perimeter and the heat score are somewhat similar in their distribution, and the heat score moves closer to the perimeter score as ∆x is decreased. And if the perfect n is recalculated using perimeter instead, we find n = 3.3 again. Therefore we conclude that the perimeter score is a suitable approximation for the heat score. For future tests then, we would use the perimeter instead of our simulation because finding the perimeter is orders of magnitude faster than finding the heat distributions.

6 Discussion

The true goal of modelling some situation is that in the future, the model could be used instead of doing any actual physical experiment. In order to assess whether our model meets these goals, we look at the predictive power of the model, and the model’s behavior under changes in the initial conditions or constants. After coming to a conclusion about these characteristics, we will then assess what are the weaknesses and strengths of our model, and what we would change if we were to do this again.

6.1 Predictive Power

The model’s predictive power is relatively strong as long as the assumptions made under section 3 are used. For any given value of n, we can see how many brownie pans will fit in the oven, and what the temperature distribution across the pan will be; however, beyond the realm of superellipses, we would need to adjust our model to allow other shapes to be analyzed in the same way. The changes would be minimal though, because regardless of the shape, the Area Score and the Heat Score or even the Perimeter Score would still yield similar results, and could still be used to find the maximum point of a shape. A possible extension to this model would be to change the assumed perfect shape of the pan to instead have the form of the superformula:

r(ϕ) =

(∣∣∣∣∣cos ( mϕ 4

) a

∣∣∣∣∣ n2

∣∣∣∣∣sin ( mϕ 4

) b

∣∣∣∣∣ n3)− 1n1

(25)

This equation allows for many more various shapes, and would be a much more general solution of this problem. If we used this formula, it would become almost infeasible to test all of the shapes possible, but we would instead rely wholly on an analytical solution, by finding the maximum of the perimeter and the minimum of the dimensions. Due to the amount of time it takes to produce a Heat Score and the difficulty of solving a partial differential equation for any given boundary conditions, we would focus on perimeter instead.

6.2 Changes in Initial Conditions

Our model does not handle changes in the initial conditions well. The Area Score is highly sensitive to changes in W , L, A, and d. By changing those values we saw behavior that varied from a constant function (every shape took up the same number of pans in the oven, this happens when the size of the oven is too small compared to the area of the pan, so that minute changes in the dimensions of the pan don’t help as much) to a function with a large gap (when d was large enough that a slightly smaller pan size allowed it to fit one more pan in, and cause a large jump in the graph). To make the Area Score more noticeable, we increased the size of a standard oven by a factor of 7, or in other words turned it into an industrial sized oven, and this produced a much more dynamic Area Score. The Heat Score is also sensitive to changes, the most dramatic of which is changing ∆x, the size of the cubes. If this value gets too large (> 0.001) the average variance actually decreases as n increases, which produces a curve in the exact opposite shape. In order for our model to be the most correct, the limit should be taken as lim∆x→0 SH . The same thing also applies to ∆t, and it should also be as close to 0 as possible. In this way, our model is sensitive to changes in the initial conditions, and thus will not adapt well to different parameters.

6.3 Simulation Error

Another note to add is that our Heat Score simulation stops at t = 1200 sec- onds, or equivalently 20 minutes; however, in section 5 it is noticeable that the equilibrium occurs after 80 minutes. We had tested this before, and the length of time observed changed the variances (in relation to each other) very little, so for the sake of the time it would take to simulate until equilibrium we short- ened it to be until a normal brownie baking time. One possible reason that the oven took so long was because under the Twelfth Assumption, we considered convection currents caused in the oven to be negligible, but they could account for a much larger percentage of the heat transfer than we realized. We would probably include convection currents if we were to do this model again.

A Simulation Analysis

In our numerical simulation in MATLAB, we have to consider any possible sources of error. One of the largest possible sources is numerical instability. Our simulation creates 2 200×200 matrices (one for current and one for the last temperature) and runs operations on it m = 1200 times. This means that the simulation is O(ijm) (i is the number of rows in one matrix, j is the number of columns in one matrix, and m is the number of times to iterate over). Since this is a large number of calculations (due to ∆x and ∆t being close to 0), the chance for error is non-negligible. In the calculation of the temperature change from the side cubes (20), we are subtracting two very similar numbers, Tside,j − Tm and summing those 4 differences. Subtraction of very similar numbers creates numerical instability; however, we could not see a way to mitigate these subtractions. Another possible source is when we compute the average variance (23), we subtract the sum of squares and the square of sums over m, which are two large numbers being subtracted from each other that are relatively close in magnitude, so there could be some instability there as well.

B Analytical Solution

An analytical solution exists for the brownies if we simplify the model more. If we assume that the pan reaches the oven temperature almost instantly (plausible given the thermal diffusivity, α2, of aluminum pans is much greater than the thermal diffusivity for chocolate brownies), and that the change in temperature from the oven emitting radiation on the top of the brownies directly is negligible (because it hits all brownies evenly, it does not affect the distribution). Under these assumptions, this is the partial differential equation we came up with to describe the analytic solution:

α2b (Txx + Tyy) = Tt (26)

Subject to the following boundary and initial conditions:

( x,b

√ 1 −

an , t

) = TO (27)

T(x,y, 0) = TR (28)

Unfortunately, we do not know of a way to solve this equation, and the under- lying assumptions are questionable, so we decided not to pursue this route.

References

[1] The ultimate brownie pan. Mathematical Contest in Modeling, 2013.

[2] R. Byron Bird, Warren E. Stewart, and Edwin N. Lightfoot. Transport Phenomena. John Wiley Sons, Inc, 1960.

[3] H. S. Carslaw and J. C. Jaeger. Conduction of Heat in Solids. Oxford University Press, 1959.