Gauss Legendre Method for Implicit Integration

Question

Methods that are usually adopted for time integration in transport phenomena problems are either:

Euler (explicit, first-order accurate)

$\frac{dY}{dt}=f(t,Y)$

$Y^{n+1}=Y^n+\Delta t f(t,Y^n)$

Backward (implicit, first-order accurate)

$Y^{n+1}=Y^n+\Delta t f(t+\Delta t,Y^{n+1})$

or Crank-Nicholson (implicit, second-order)

$Y^{n+1}=Y^n+\frac{\Delta t}{2}(f(t,Y^{n})+f(t+\Delta t,Y^{n+1})$

I am trying to use a higher order scheme so that larger time steps can be used at the same amount of error. I have recently been looking at the Gauss-Legendre Method that has a butcher tableau of this from the article here.

1. Given the Butcher tableau for Gauss-Legendre, is it possible to formulate an expression similar to that above for Euler, Backward Euler, and CN? Right now I have

$Y^{n+1} = Y^n+\frac{1}{2}k_1+\frac{1}{2}k_2$

where

$k_1=\Delta tf(t+(\frac{1}{2}-\frac{1}{6}\sqrt(3)\Delta t,Y^n+\frac{1}{4}k_1+(\frac{1}{4}-\frac{1}{6}\sqrt{3})k_2)$ $k_2=\Delta tf(t+(\frac{1}{2}+\frac{1}{6}\sqrt(3)\Delta t,Y^n+(\frac{1}{4}+\frac{1}{6}\sqrt{3})k_1+\frac{1}{4}k_2)$

but I'm not sure where the $Y^{n+1}$ appears in this formulation on the RHS of $k_1$, $k_2$ or the first equation.

2. Are their any limitations to the transport equations? Ideally it is something like $f(t,C) = -\nabla\cdot(\vec{U}C)+D\nabla^2C$

3. Are there any other higher order implicit time integration schemes that are a) a lower number of stages and b) higher order of accuracy than 2?

Help is much appreciated and or advice on higher-order time integration schemes would be nice

Answer 1

1. If I understand your first question, the answer is "no". $ k_1 $ and $ k_2 $ are defined implicitly in terms of each other, so it is not easy to write a simpler one-line expression. You could write \begin{align*} k_1 &= \Delta t f( t + (1/2 - \sqrt{3}/6) \Delta t, (1/2) Y^{n+1} + (1/2)Y^n - (\sqrt{3}/6) k_2) \\\\ k_2 &= \Delta t f( t + (1/2 + \sqrt{3}/6) \Delta t, (1/2) Y^{n+1} + (1/2)Y^n + (\sqrt{3}/6) k_1) \end{align*} but I'm not sure that helps you much.

2. Yes, the structure of your equation (in particular the eigenvalues of whatever discrete operator you use to represent the transport equation) will interact with the time stepping scheme you use to determine what size time steps you can use for stability. One search term you need is "stability domain".

3. There are a whole zoo of different time integration schemes, but I don't think there are any that fit your requirements. You might want to consider linear multistep schemes (the major category besides Runge-Kutta), but it really depends on what you want to do. You could look at Chapter 6 of Atkinson's Numerical Analysis book, but really any numerical analysis textbook would get you started.

Answer 2

It's been a long time since I looked at this sort of thing, but some of the references here might be useful:

ftp://ftp.math.ucla.edu/pub/camreport/cam04-26.pdf

The report and its references describe how to construct RK schemes with large regions of absolute stability.