MAT 402 - EBM 607 Mini Project In this project we will go over all steps of scientific computing. We will first derive the differential equation of a physical phenomenon, namely, unsteady heat conduction. Then, with some simplifying assumptions, we will derive a one-dimensional (1D) form from the three-dimensional (3D) differential equation. Next, we will devise a numerical solution method, solve the problem in Matlab and analyze the results. 1. 3D Heat Equation. Derive the heat (or diffusion) equation in 3D, which is given by T 1 t c kT f (1) c where, T T x , y , z , t is the temperature field, c is the specific heat capacity of the medium (J/(kg/K)), is the density of the medium, is the thermal conductivity of the medium (W/(mK)), and x , y , z , t is the rate of heat generated in the medium. To derive the heat equation, pick an arbitrary volume V in 3D medium and let S be the surface of V. If the mass contained in V is M, then heat in volume V is M c T , where T is the absolute temperature in Kelvin (K). Then, rate of accumulation of heat in V is k f f d d cM T dt The density of heat flux (J/m2) in a direction dt n c T dV V c V T t dV . (2) is given by Fourier's law: (3) q kT n Find the total heat entering the volume V in the form of a surface integral through S. Finally, express the total heat generated within V using f and a volume integral. To ensure energy balance, rate of increase of heat within V given by (2) must equal the rate of heat entering S (given by surface integral that you will find) plus total heat generated within V. Then you must convert your surface integral through S to a volume integral in V using the divergence theorem. This will give you an integral equation where all quantities are expressed inside integrals. Finally, since the volume is chosen arbitrarily, the equation must hold for the integrands. That equality will give you the heat equation given by (1). 2. 1D simplification. Assume that T does not change in x and y dimensions. Also, assume that the heat conductivity is constant (i.e., we have a homogenous medium) and there is no heat generation within the volume of interest. Then, the partial differential equation given by (1) becomes a 1D heat equation: T t k T 2 c x 2 f c . (4) However, we will have some further simplification and assume that no heat is generated in our 1D domain. Hence, (4) becomes T t T 2 2 x 2 (5) Note that if steady state is achieved, (5) turns out to be the Laplace equation 2 T 0 . With appropriate boundary conditions, finite discretization of Laplace equation results a tridiagonal matrix system in 1D whose solution provides the steady-state temperature values. The story is quite different for the heat equation as we will see. 3. Numerical Method. For notational simplicity we will denote temperature with u and use subscript notation for the derivative, e.g., u x denotes u / x . Then, our complete mathematical problem definition with boundary conditions and an initial condition will be 0 u xx u t , 2 x L, 0 t u 0, t p t 0 t u L,t q t 0 t u x, 0 f x 0 (6) x L To make use of the finite difference method, divide L into N equal parts of length x L / N and define x j j x , j 0 ,1, , N . Further, choose a time increment t and discretize time using t k k t , k 0 ,1, 2 schemes. Denote . Next, discretize u x j , tk to Euler's method (i.e., as U U j ,k 1 j,k U u xx using central difference and ut using forward difference for notational simplicity. Finally, obtain a time-update scheme similar j,k ) that will allow you to evolve the solution in time. Implement your numerical method in Matlab. You need to write a function that will input parameters and functions in (6) and return the values of u. You can pick small numbers for N and choose a simple problem to test your code. Explain this step in detail and show that your code works. 4. Test case. When you implement a numerical method or an algorithm, you need a test case to debug your code. The test case is a specific problem for which you can find the analytical solution of the heat equation. If 2 1, L 1, p t q t 0 , a n d f ( x ) 1 0 0 then the solution to heat 1D heat equation is u x,t 400 n 1 , 3 ,... 1 n s in n x e n 2 t . Choose N 5 0 , t 0 .0 0 0 1 8 and compare your numerical solution with this analytical solution. Find U for k 1 0 , 5 0 ,1 0 0 , a n d 5 0 0 and plot your numerical solutions on top of the analytical solution. Comment on your results. The contribution of this project to your overall grade will be 30 %. Send your detailed Word report showing your work in each step defined here with codes until the deadline of the assignment, which is 3rd January 2016, Sunday, 23:59. MAT 402 - EBM 607 Mini Project In this project we will go over all steps of scientific computing. We will first derive the differential equation of a physical phenomenon, namely, unsteady heat conduction. Then, with some simplifying assumptions, we will derive a one-dimensional (1D) form from the three-dimensional (3D) differential equation. Next, we will devise a numerical solution method, solve the problem in Matlab and analyze the results. 1. 3D Heat Equation. Derive the heat (or diffusion) equation in 3D, which is given by T 1 t c kT f (1) c where, T T x , y , z , t is the temperature field, c is the specific heat capacity of the medium (J/(kg/K)), is the density of the medium, is the thermal conductivity of the medium (W/(mK)), and x , y , z , t is the rate of heat generated in the medium. To derive the heat equation, pick an arbitrary volume V in 3D medium and let S be the surface of V. If the mass contained in V is M, then heat in volume V is M c T , where T is the absolute temperature in Kelvin (K). Then, rate of accumulation of heat in V is k f f d d cM T dt The density of heat flux (J/m2) in a direction dt n c T dV V c V T t dV . (2) is given by Fourier's law: (3) q kT n Find the total heat entering the volume V in the form of a surface integral through S. Finally, express the total heat generated within V using f and a volume integral. To ensure energy balance, rate of increase of heat within V given by (2) must equal the rate of heat entering S (given by surface integral that you will find) plus total heat generated within V. Then you must convert your surface integral through S to a volume integral in V using the divergence theorem. This will give you an integral equation where all quantities are expressed inside integrals. Finally, since the volume is chosen arbitrarily, the equation must hold for the integrands. That equality will give you the heat equation given by (1). 2. 1D simplification. Assume that T does not change in x and y dimensions. Also, assume that the heat conductivity is constant (i.e., we have a homogenous medium) and there is no heat generation within the volume of interest. Then, the partial differential equation given by (1) becomes a 1D heat equation: T t k T 2 c x 2 f c . (4) However, we will have some further simplification and assume that no heat is generated in our 1D domain. Hence, (4) becomes T t T 2 2 x 2 (5) Note that if steady state is achieved, (5) turns out to be the Laplace equation 2 T 0 . With appropriate boundary conditions, finite discretization of Laplace equation results a tridiagonal matrix system in 1D whose solution provides the steady-state temperature values. The story is quite different for the heat equation as we will see. 3. Numerical Method. For notational simplicity we will denote temperature with u and use subscript notation for the derivative, e.g., u x denotes u / x . Then, our complete mathematical problem definition with boundary conditions and an initial condition will be 0 u xx u t , 2 x L, 0 t u 0, t p t 0 t u L,t q t 0 t u x, 0 f x 0 (6) x L To make use of the finite difference method, divide L into N equal parts of length x L / N and define x j j x , j 0 ,1, , N . Further, choose a time increment t and discretize time using t k k t , k 0 ,1, 2 schemes. Denote . Next, discretize u x j , tk to Euler's method (i.e., as U U j ,k 1 j,k U u xx using central difference and ut using forward difference for notational simplicity. Finally, obtain a time-update scheme similar j,k ) that will allow you to evolve the solution in time. Implement your numerical method in Matlab. You need to write a function that will input parameters and functions in (6) and return the values of u. You can pick small numbers for N and choose a simple problem to test your code. Explain this step in detail and show that your code works. 4. Test case. When you implement a numerical method or an algorithm, you need a test case to debug your code. The test case is a specific problem for which you can find the analytical solution of the heat equation. If 2 1, L 1, p t q t 0 , a n d f ( x ) 1 0 0 then the solution to heat 1D heat equation is u x,t 400 n 1 , 3 ,... 1 n s in n x e n 2 t . Choose N 5 0 , t 0 .0 0 0 1 8 and compare your numerical solution with this analytical solution. Find U for k 1 0 , 5 0 ,1 0 0 , a n d 5 0 0 and plot your numerical solutions on top of the analytical solution. Comment on your results. The contribution of this project to your overall grade will be 30 %. Send your detailed Word report showing your work in each step defined here with codes until the deadline of the assignment, which is 3rd January 2016, Sunday, 23:59