Hello, dear friend, you can consult us at any time if you have any questions, add WeChat: daixieit

EMAT30008:  SCIENTIFIC COMPUTING

Part 3:  Finite Dierence Methods and PDEs

 

1    Overall aim

Remember that the aim of the coursework is to integrate the work you do from each week into a single piece of software.  The software should be a general numerical continuation code that can track, under variations of a parameter in the system,

● limit cycle oscillations of arbitrary ordinary differential equations (of any number of dimensions), and

● steady-states of second-order diffusive PDEs.

Moreover,  you  should provide time  simulation  codes  for both  ODEs  and  PDEs to enable the  solu- tions to be verified. Full details are on the course website: https://engmaths.github.io/emat30008/ assessment/.

 

2    Steps to get there

It’s best to start with the time simulation code, and work up from there. Start slowly, change one thing at a time, and keep your code as modular as possible.  Don’t worry about the numerical continuation part till you’re comfortable with the time-domain solver.

Some suggestions for dividing the overall task into manageable pieces are given below. You could aim to complete items 1-3 by the end of week 20, item 4 by the end of week 21, and then focus on items 5 and 6 after that. You should be able to work on item 7 throughout!

If it helps, please feel free to team up and share ideas (see, e.g., https://en.wikipedia.org/wiki/ Pair_programming), for items 1–3.

1. First, start off by playing with the Python code for the forward Euler scheme, as available on Blackboard.   Experiment with different discretisations  (mesh parameters mx and mt) as well as different problem parameters (L, T, and s), so you’re confident you know how it works.  Use the code to verify the forward Euler stability criterion for yourself.

2. What is the steady state of the PDE problem that is solved by the code on Blackboard.  Does it make (physical) sense to you? What happens if we change the initial condition — for example, to t(u, 0) = sinp (mu) for some integer 夕?

3.  Modify the Python code to replace the component-wise forward Euler method with the matrix/vec- tor approach (see page 15 of the lecture notes).  Make sure the solutions you obtain are the same as from the component-wise method.

4.  Modify the code so that it uses backward Euler and/or Crank-Nicholson (nothing complicated; just ordinary linear solve will do). You should get the right solution but with no stability problems, and the Crank-Nicholson solver should converge faster than either of the Euler methods.

5. Adapt the code so that it can deal with more general parabolic PDE problems.  Some extensions that you might wish to consider (there’s no need to do them all!):

● non-homogenous Dirichlet boundary conditions, so (for example) t(0, r) = t0 (r)atau

● more general boundary conditions, e.g. a Robin (mixed) boundary condition αt + β      = g

at u = 0, L, or periodic boundary conditions t(0, r) = t(L, r) for all r,at        a2t

at            a2 t

  the cable equation Cm ar  = gmau2  + gL (tL − t), where Cm  > 0, gM  > 0, gL  > 0, and tL are constants

 Fisher’s equation       = s       + pt(1 − t) where s > 0 and p > 0 are constants.

6. Adapt your numerical continuation code (from the ODEs part of the course) so that it can also track, under variations of a parameter in the system, steady-states of second-order diffusive PDEs. What parameter a user could ask to vary is one of your design choice; e.g., the value of t on a boundary, a parameter in a heat source function, or any other appropriate part of the PDE problem.

7.  Turn the code into good  Python code, that you can call to solve a parabolic PDE problem.  You should be  applying  sound software engineering principles,  as  for the  ODE part of the  course, throughout.

You could, for example,

● vectorize the code, rather than using loops,

● break it down into subroutines/internal function calls,

● create a function library,

● improve the computational efficiency of the linear algebra; e.g., linear solve for a tridiagonal matrix, sparse matrix operations, etc.

● add visualisation functionality; e.g., plotting, snapshots, animation, diagnostics, etc.,

● make the design easily adaptable for future modifications,

● make sure the code is fully documented.

It would be a good idea to check your code is working by testing against known (analytic) solutions wherever possible. Some information on explicit solutions of the 1D heat equation (and other PDEs) can be found here: https://people.maths.bris.ac.uk/~mazvs/supplementary.pdf.