advertisement

18.086 - Spring 2016, Problem Set 1 February 9, 2016 Out: Feb. 9 Due: Feb. 23 Problem 1 (10 points) Consider the one-dimensional heat equation ut = c2 uxx , (1) where −∞ < x < ∞ and u(x, t) is the temperature field. Suppose the initial heat profile at t = 0 is u(x, 0) = f (x) for some function f (x). In the following we denote byˆthe Fourier transform of a function. 1. Apply the Fourier transform method to Eq. (1) to obtain an ODE in 2 2 time. Show that û(k, t) = fˆ(k)e−c k t is its solution satisfying the initial conditions (fˆ(k) is the Fourier transform of f (x)). 2. Transforming back, show that the solution to Eq. (1) is given by Z ∞ 2 2 1 u(x, t) = fˆ(k)e−c k t+ikx dk 2π −∞ Problem 2 (20 points) Consider the nonlinear PDE ut = −uxxxx − 2uxx + ( − 1)u − cu3 , (2) with a fixed control parameter. For the moment, set c = 0. 1. Use the ansatz u(x, t) = Aeikx+σt to show that the trivial solution u(x, t) = 0 is stable for < 0, i.e. no modes have σ(k) > 0 for < 0. 2. Using the dispersion relation σ(k) from part 1, show that for > 0 modes with kc2 = 1 are the first to become unstable. Note: These solutions blow up, which is why one needs the cu3 term to limit the growth for large amplitudes. For c = 1, it can be shown that the above equation has a stable, stationary solution with a wavelength λ = 2π/kc . 3. Confirm these findings numerically: Let u be defined in the interval [−20, 20] with periodic boundary conditions and choose a spatial resolution of at least N = 60 grid points. a) Derive a centered five-point stencil for the space derivative uxxxx . Hint: Taylor series! b) Use a) to implement a Forward-Euler time integration to solve Eq. 1 (2) with c = 1. Start with an initial solution u0 (x, t) = R(x), where for all x ∈ [−20, 20], R(x) is a small random number, R(x) ∈ [−0.01, 0.01]; a time step of ∆t ≈ 0.02, and show that, indeed, for < 0, the solution u(x, t) = 0 is obtained, while > 0 leads to a periodic pattern with wavelength λ = 2π/kc . Problem 3 (20 points) Consider the classical, planar 2-body problem of two point-masses m1 and m2 obeying Newton’s second law of motion, ~ ¨1 ~x = − mf1 ¨2 ~x = f~ = (3) f~ m2 1 m2 G |~xm x1 x2 |3 (~ 1 −~ (4) − ~x2 ) . (5) f~ is the gravitational force between the two bodies and G the gravitational constant. ~x1 (t) and ~x2 (t) are the positions of the two bodies and G is the gravitational constant. In this problem, we are interested in the trajectories for the initial conditions ~x1 (0) = (−1, 0), ~x2 (0) = (1, 0), ~x˙ 1 (0) = (0, 1), ~x˙ 2 (0) = (0, −1), and m1 = m2 = 0.1, G = 100. 1. Reformulate the problem as an ODE u̇ = F (u) where u is now a vector and F a vector-valued function. Integrate this system in time using the forward Euler method. Perform at least N = 10000 iterations with time step ∆t = 0.005 (depending on your implementation you might need to choose a smaller timestep!) and compare the trajectories to Runge-Kutta 4. 2. The leapfrog method seen in lecture can be applied separately to mass point 1 and 2. For mass point 1, it reads ~xj,n+1 − 2~xj,n + ~xj,n−1 f~ =− 2 ∆t m1 (6) and similar for mass point 2. A more intuitive form often found in literature is treating position and velocities separately. It reads: ~ x ~ ~vj,n+1/2 = vj,n−1/2 − ∆t mf1 ~xj,n+1 = ~xj,n + ∆t~vj,n+1/2 (7) −~ x j,n where ~vj,n+1/2 = j,n+1 . As you can see, velocities are calculated ∆t at half-time steps, and position updates mix values from integer and halfinteger time steps. This has (somehow) inspired people to call this method the leapfrog algorithm. Show that both formulations are equivalent! 3. Implement the position-velocity leapfrog scheme.PUse forward Euler for 2 m2 the starting step. Compare the total energy, E = i=1 21 mi~vi2 − G |~xm11−~ x2 | , with the RK4 solution for ∆t = 0.03 after at least 30000 time steps. Which method would you use a) to calculate the trajectory of a comet to determine how close to earth it will pass, and b) to study if our solar system remains stable if the sun had twice its mass? 2