Exam INF5620, fall 2012

• There is also a PDF version of this document.

We will throw a die and the number of eyes determines the problem to be presented. Thereafter, you will be given some questions, either about parts of your presention or facts from the other problems. After each presentation, the next candidate can throw the die and get about 10 min to collect the thoughts before presenting the assigned problem.

Problem 1: Falling body

The equations for the velocity \( v \) of a very small spherical body with density \( \varrho \) and radius \( r \) falling in a fluid with dynamic viscosity \( \mu \) and density \( \varrho_f \) is $$ \frac{4}{3}\pi r^3 \varrho v'(t) = -\frac{4}{3}\pi r^3 \varrho g - 6\pi r\mu v + \frac{4}{3}\pi r^3\varrho_f g ,$$ where \( g \) is the acceleration of gravity. The initial condition is \( v(0)=0 \).

1. Explain briefly how this equation arises from basic principles in physics and what the individual terms model.
2. Derive a Forward Euler, Backward Euler, and Crank-Nicolson scheme for the equation. Mention other possible schemes too.
3. Illustrate what kind of numerical artifacts that may appear when using the Forward Euler, Backward Euler, and a Crank-Nicolson schemes. Explain the reason for the artifacts (motivated by a mathematical analysis of the schemes).
4. Which one of the three schemes will you recommend for solving this equation with a) large time steps and b) small time steps?
5. Imagine that somebody claims that the 4th-order Runge-Kutta method is superior to the Forward Euler, Backward Euler, and a Crank-Nicolson schemes both for large and small time steps. Perform numerical experiments to assess if this claim is correct.
6. The equation above is not a good model if \( \varrho_f vr/\mu \) is much greater than 1, which is the case for a body that is not very small. How can the model be extended to cover this case? Suggest a numerical scheme for the modified equation.
7. Suppose the shape of the body is much more complicated than a sphere so that simple fluid resistance formulas are too inaccurate. Explain briefly how one can compute (in principle) an accurate drag force on the body.

Problem 2: Heat conduction with finite differences

The temperature distribution \( T(x,y,z,t) \) in a solid 3D body \( \Omega \) is governed by the heat equation $$ \varrho c T_t = \nabla\cdot (\kappa(x,y,z)\nabla T),\quad \pmb{x}\in\Omega, t\in (0,T] $$ Here, \( \varrho \) is the density of the body, \( c \) is a measure of material's heat capacity, \( \kappa \) is the heat conduction coefficient.

We consider heat conduction in a long cylinder where the sylindrical surface is isolated such that \( -\kappa\partial T/\partial n=0 \) here. The cross sections of the end of the cylinder are kept at constant temperatures. The left half of the cylinder is made of a material with constant heat capacity \( c_0 \) and constant heat conduction coefficient \( \kappa_0 \), while the right half has the corresponding constant values \( c_1 \) and \( \kappa_1 \). A time \( t=0 \) the two pieces, with different temperatures, are brought together such that the initial temperature field is \( T_0 \) in the left piece and \( T_1 \) in the right piece. Imagine that the left and right end points are kept at \( T_0 \) and \( T_1 \), respectively, at any time \( t>0 \).

1. Show that the simplification \( T=T(x,t) \) is possible in the described problem, where \( x \) is a coordinate along the cylinder (just insert \( T(x,t) \) in the original problem and see that it fulfills all equations). Set up the simplified PDE with proper boundary and initial conditions.
2. The 1D PDE problem is discretized by the Forward Euler, Backward Euler, or Crank-Nicolson schemes. Derive the discrete equations for one of these schemes.
3. Assume for simplicity that \( c_0=c_1 \) and that \( \kappa_0=\kappa_1 \). With a discontinuous initial conditions, numerical artifacts may appear in the solutions produced by the Backward Euler, Forward Euler, and Crank-Nicolson methods. Illustrate such artifacts. A possible program to play around with is demo_osc.py.
4. Present the ideas of an analysis that can explain the artifacts in the previous subproblem. Summarize the main findings from this analysis for the Backward Euler, Forward Euler, and Crank-Nicolson schemes. What are the stability restrictions of the various schemes?
5. To what extent are the artifacts and their explanations from the analysis from Problem 1 relevant here in Problem 2?

Problem 3: Heat conduction with finite elements

We consider the same PDE problem as in Problem 2, but this time solved by finite elements in space in 1D. For simplicity you can assume that \( c \) and \( \kappa \) are constants and that we are in 1D.

1. Use a finite difference scheme of your choice in time and derive a series of spatial problems. Derive variational formulations of these spatial problems.
2. Use P1 finite elements for the discretization in space. Show in detail how the element matrix and vector corresponding to the time derivative term or the diffusion term is computed. (Consider an arbitrary element in the interior of the mesh.)
3. Explain how the boundary conditions are incorporated in the finite element method. (You may want to compute one matrix for the time derivative term and one matrix for the diffusion term. At each time level, you form the coefficient matrix of the linear system and the right-hand side, and then you may incorporate Dirichlet conditions in this system. Alternatively, it is possible to prepare the element matrices and vectors so that they get the right form wrt Dirichlet conditions after assembly.)
4. Explain the method of manufactured solutions and how it can be used to verify an implementation of the numerical method in this problem (you may add a source term in the equation if desired).
5. Compare the difference equations arising from the finite element method with the corresponding equations arising from the finite difference method by expressing the finite element equations in terms of finite difference operators.

Problem 4: Wave equations with finite differences and elements

1. Set up a wave equation problem in 2D with zero normal derivative as boundary condition. Assume a variable wave velocity.
2. Sketch a physical problem where your mathematical model arises.
3. Present a finite difference discretization. Explain in particular how the boundary conditions and the initial conditions are incorporated in the scheme.
4. Explain (in princple) how the 2D discretization can be extended to 3D.
5. A pure Python implementation of a 2D or 3D scheme runs slowly compared to implementations in Fortran, C, or C++. Describe ways of speeding up the calculations of a Python implementation.
6. The scheme for the wave equation is perfect for parallel computing. Why? What are the principal ideas behind a parallel version of the scheme?
7. Numerical artifacts may occur when solving wave equations. Illustrate such artifacts. Sketch an analysis that can explain why the artifacts appear.
8. What is the stability restriction of the scheme?

Problem 5: Nonlinear diffusion problem

Consider the following problem: $$ - \frac{d}{dx}\left( \alpha(u)\frac{du}{dx}\right) = f(u),\quad x\in (0,1),\quad u(0)=0, u(1)=1 $$

1. Formulate a (standard) finite difference method for this problem using arithmetic averages of \( \alpha \) where needed.
2. Derive expressions for the linear system to be solved in each Newton iteration.
3. Derive a variational form for this problem.
4. Derive an expression for the Jacobian (in Newton's method) from the variational form.
5. Compute the discrete nonlinear equations using P1 elements and the Trapezoidal rule, and show that the equations are then identical to those coming from a finite difference method.
6. What is the group finite element method and how can it be applied here?

Problem 6: Advanced application

Derive a variational form for a system of differential equations of your choice. Explain how to apply the finite element method and what the overall computational algorithm for the problem is. You may find examples from anywhere to illustrate what type of computations that can be done with the mathematical model.

Some possibilities based on lectures are

• \( \mu\nabla^2 w = -\beta \), \( \kappa\nabla^2 T = -\mu \nabla w\cdot\nabla w \)
• \( w_t = \nabla\cdot (\mu ||\nabla w||^q \nabla w) + \beta \), \( \kappa\nabla^2 T = -\mu \nabla w\cdot\nabla w \)
• Time-dependent or stationary elasticity
• The Navier-Stokes equations