|Using EXCEL Spreadsheets to Evaluate the Fourier Series Solutions to a 1D Heat Equation|
The goal of this tutorial is to create an EXCEL spreadsheet that calculates the Fourier series solution to the following initial-boundary value problem for the one-dimensional heat equation:
The basic idea of finding a series solution is to expand the unknown function u(x, t) in a series of eigenfunctions that satisfy the same boundary conditions as the original problem. For the boundary-value problem described above, the general expansion of u over a spatial domain [0, L] is
Note that the coefficients bn are the coefficients of f(x) in a Fourier series. Thus, for a given initial condition f(x), the function u(x, t) is completely determined.An Example of Initial Condition
Suppose the initial condition has the following form
Using the above expression for the Fourier coefficients bn with L = 1, one finds
While substituting the expression for bn into the expansion of u, it is worth noticing that
sin(nπ/2) = 0 for even values of n
sin(nπ/2) = (-1)n for odd values of n
Hence, half of the terms in the expansion of u are null and the expression reduces to
Basic Steps to Set Up a Spreadsheet
These are the basic steps needed in setting up an EXCEL spreadsheet that allows you to calculate the first few terms in the series solution derived above. The terms will be calculated at specific points in the spatial domain and at a given time.
The first step is to tabulate the sample points in the spatial interval. These will be the points at which the expansion of u will be evaluated.
Subdivide the interval [0, 1] into N+1 equally-spaced sample points a distance Δx apart. In a blank EXCEL spreadsheet, enter the sample points along a row, as shown below. For this tutorial, let Δx = 0.05 and N = 19. Note that the value of Δx is shown in cell C2.
Each sample point is calculated as
for n = 1, ... N and x0 = 0.
For example, in the spreadsheet below the value of the first sample point x1 = x0 + 0.05 is shown in cell E5: it is computed by addying the fixed space increment (cell C2) to the value contained in the cell to the immediate left of E5 (cell D5). The expression displayed in the formula box shows the specific EXCEL formula used to compute the value in cell E5. The same formula can be copied and pasted in the cells along row 5 to calculate the desired sample points. Note that the value of the space increment Δx should always be read from cell C2, therefore the address of cell C2 should be displayed with dollar signs in the formula (see absolute referencing).
The second step is to tabulate the mode number corresponding to the index k in the expansion of u.
The first mode number is k=0. In the EXCEL spreadsheet enter 0 in a column different from those used to store the spatial points, for example choose cell B7, as shown in the image below. The sequence of increasingly higher mode numbers will be entered in the rows below cell B7. Note that you do not need to enter each mode number by hand. Simply add 1 to the previous mode number. For example, the value in cell B8 is obtained by adding 1 to the value in cell B7, as displayed in the formula box. The same formula used in cell B8 can be copied and pasted into as many cells in column B as the modes you want to compute.
The kth term of the expansion of u has the form
Its value depends on four quantities: the mode number k; the point x and the time t at which the expansion is calculated; the value of the constant α2 characteristic of the physical properties of the system under consideration.
In this tutorial, take α2 = 1 and t = 0.05; these values will be stored in some cells in row 2 so that they are easily accessible during the computation. The other quantities, x and k, will be determined by the sample points stored in row 5 and the mode numbers coded in column B.
Each term of the expansion will be stored in a cell. For example, the first term calculated at x0 = 0 will be stored in cell D7, that is in the same row as the corresponding mode number and in the same column as the corresponding sample point, as shown in the image below where the spatial sample points are highlighted in blue for clarity.
The EXCEL formula used to compute the value in cell D7 is displayed in the formula box. The formula is obtained from the expression for the kth term shown above in which the value of k is read from cell B7 and the value of x is read from cell D5. Specifically, the values used in coding the formula in cell D7 are:
Note how some cell addresses in the formula contain dollar signs "$". Even though it is not necessary to use dollar signs in order to code a formula that computes just one mode in the expansion of u, it will become extremely useful as you attempt to fill in the rest of the spreadsheet. In fact, in order to effortlessly compute a large number of terms in the expansion it is convenient to simply paste and copy this formula into other cells, which makes cell referencing a critical aspect of how the formula should be coded, as explained next.
For an efficient use of the spreadsheet, it is essential to store information in an orderly manner. The sample points are stored in adjacent columns along one row, the mode numbers are entered in subsequent rows in one column: the same structure should be used to compute the terms of the expansion. For example, if the first term of the expansion at x0 = 0 is stored in cell D7, then it is convenient to store the second term computed at the same point x0 in the cell immediately below D7, or cell D8, consistently with the way mode numbers are tabulated. Similarly, the value of the first mode of the expansion computed at the next sample point x1 should be stored in the cell to the right of cell D7, or cell E7, following the tabulation of the sample points. This way of storing information in columns and rows will allow you to organize your computation in a structured manner.
As mentioned earlier, instead of typing a new formula for each term of the expansion, the same formula used to compute the first term at x0 can be copied and pasted into other cells of the spreadsheet. When you do so, EXCEL will automatically update the cell addresses in the formula based on the position of the new cell in which the formula is pasted relative to the original cell. So, if the formula in cell D7 is copied and pasted into, say, cell D8, that is, into a cell in the same column of but one row below the original cell, all cell addresses in the formula will be updated to reflect the change of row, with the exception of those containg dollar signs. Any row or column label preceded by a dollar sign will not be changed when the formula is pasted into a new cell.
For example, in the formula used in cell D7 there are three quantities whose value you do not want to change as a result of copying and pasting the formula into a different cell. These are the number π, the value of the constant α2, and the fixed value of t. In order to keep the addresses of the cells containing these values unchanged when the formula is pasted into a different cell, a dollar sign is added in the cell address in front of both its column and row label. Thus, e.g. the cell containing the number π is called by $J$2 instead of J2. This is called absolute referencing of cells. An EXCEL formula, however, can also be coded using mixed cell referencing. This happens when a cell address has only one dollar sign, either in front of the column or the row label. In this case, when the formula is pasted into a new cell only the label that is not preceded by the dollar sign is updated. For example, the formula in cell D7 contains references to cells B7 and D5, which store the mode number and the sample point, respectively. If no dollar signs were used when referencing these cells and the formula were copied and pasted, say, into cell E8 to compute the second term of the expansion at the next sample point x1, EXCEL would increase all row and column labels not preceded by $ by one, that is, any reference to B7 would change to C8, and any reference to D5 would change to E6. In this case the computation would be incorrect. To perform the correct computation, instead, the formula should read the new mode number from cell B8 and the new spatial point from cell E5. Thus, any reference to B7 should change to B8, and any reference to D5 should change to E5. In other words, the column label B in the addresses referencing the cell containing k should be fixed, and so should be the row label 5 in the addresses referencing the cell containing x. Therefore, in order to fill in the spreadsheet by copying and pasting the formula used to compute the first term of the expansion at x0, the formula should be coded using $B7 instead of B7, and D$5 instead of D5.
Keep in mind that the Edit|Fill|Fill Down and Edit|Fill|Fill Right features will allow you to fill in cells quickly and thus calculate as many modes as you like effortlessly.
By filling in the cells in columns D through X in the spreadsheet described in the previous steps, one finds the individual terms of the expansion of u at different points in the domain at a fixed time. The number of modes computed will depend on the number of rows filled in the spreadsheet. Eventually, all modes considered should be added to find the value of the (truncated) expansion at each sample point. It is easy to observe how the modes of the expansion descrease quickly to zero. Thus it is not necessary to calculate a large number of modes to find an accurate solution u(x, t) to the problem.
To find the value of u(x, t) at each sample point, the numbers in each column should be added together. EXCEL allows you to do that quickly by using the built-in function SUM, as shown in the image below where the values of u(x, t) at each sample point are calculated for k = 30 in cells D39 through X39 and highlighted in red for clarity. The correct syntax for the function SUM used in cell D39 is shown in the formula box. Recall that in this tutorial t = 0.05.
An approximated solution to the initial-boundary value problem described in this tutorial can be obtained numerically by calculating a finite difference approximation to the one-dimensional heat equation, as described in Tutorial 1. The two solutions can then be plotted on the same graph for comparison.
Begin by plotting the series solution calculated in Step 4:
Then, create a new worksheet in your spreadsheet by clicking on the "Sheet 2" tab. Use this worksheet to compute the finite difference solution up to a time level t = 0.05. To fill in this new worksheet you can simply copy the spreadsheet described in Tutorial 1 and paste it into your new worksheet. Make sure that the finite difference calculation solves the same exact initial-value boundary problem as the one described in this tutorial: check that the initial condition, the boundary conditions, and the value of the constant α2 are the same as those used to find the Fourier series solution. For increased numerical accuracy, choose Δt = 0.0005.
Now add the plot of the solution calculated using the finite difference scheme to the chart showing the series solution that was created in Sheet 1:
The resulting plot will look like the image below. Note the excellent agreement between the two solutions.