Approximations used in Crank-Nicolson method for solving PDEs numerically

In summary: I should have realized that. Anyway, if you wanted to do the central difference for the derivative, you would also need u_{i}^{n+2} if you didn't want to introduce any time legs. This means you would need twice as many equations. This would obviously slow down your algorithm but it might help...
  • #1
pafcu
10
0
Hi all,

A diffusion equation is of the form
[tex]\frac{\partial u}{\partial t}=k \frac{\partial^2u}{\partial x^2}[/tex]

Usually an equation like this seems to be solved numerically using the Crank-Nicolson method:
[tex]\frac{u_{i}^{n + 1} - u_{i}^{n}}{\Delta t} = \frac{k}{2 (\Delta x)^2}\left((u_{i + 1}^{n + 1} - 2 u_{i}^{n + 1} + u_{i - 1}^{n + 1}) + (u_{i + 1}^{n} - 2 u_{i}^{n} + u_{i - 1}^{n})\right)[/tex]

i.e the time-related part is approximated using a 2-point forward difference approximation and the space-part using an average of two 3-point central difference approximation. ([tex]u_{i}^{n}[/tex] means value of u at point i at timestep n, i.e u(i,n))

Is there a specific reason why a central difference is not used for the time? Since a central difference results in a smaller error, yet it's not usually used, what drawback does it have?

What about the central difference for the space-related part? Is there a reason for using the 3-point version instead of the 5-point version? Or is the simpler one usually used because it results in a tridiagonal matrix which is easier to solve?

edit: Clarified notation
 
Last edited:
Physics news on Phys.org
  • #2
pafcu said:
Is there a specific reason why a central difference is not used for the time? Since a central difference results in a smaller error, yet it's not usually used, what drawback does it have?

I don't know anything about this method but How would know what the value is in the future time step before the derivative was estimated? (I think though Rung Kutta Methods for ODE's try and do what you suggest).
 
  • #3
John Creighto said:
I don't know anything about this method but How would know what the value is in the future time step before the derivative was estimated? (I think though Rung Kutta Methods for ODE's try and do what you suggest).

I was unclear. The point is that I wish to know the value of [tex]u[/tex] at the next time step at every point, solving the mentioned equation for [tex]u_{i}^{n + 1}[/tex] for each [tex]i[/tex] (requires solving a set of linear equations) tells me just that.
 
  • #4
pafcu said:
I was unclear. The point is that I wish to know the value of [tex]u[/tex] at the next time step at every point, solving the mentioned equation for [tex]u_{i}^{n + 1}[/tex] for each [tex]i[/tex] (requires solving a set of linear equations) tells me just that.

In that case I'm surprised you have enough equations given it looks like there will be at least two unknowns per equation. However, I suppose you use an integration step to get that other unknown.
 
  • #5
John Creighto said:
In that case I'm surprised you have enough equations given it looks like there will be at least two unknowns per equation. However, I suppose you use an integration step to get that other unknown.

No, number of unknowns is not a problem.

Assuming [tex]u[/tex] is defined at m [tex]m[/tex] points (it has been discretized) we obviously end up with [tex]m[/tex] unknowns: [tex]u_i^{n+i}, 1 \leq i \leq m[/tex].
We get [tex]m-2[/tex] equations by plugging in [tex]i=2,3,4,...,m-1[/tex] into the equation in the original post. Plugging in [tex]i=1[/tex] doesn't work because [tex]u_{0}^{n+1}[/tex] required by the approximation isn't defined. Similarly you can't plug in the value [tex]i=m[/tex] because [tex]u_{m+1}^{n+1}[/tex] isn't defined.

Is all hope lost? No, because if we additionally specify two boundary conditions (say [tex]u_0^{n+1}=0[/tex] and [tex]u_m^{n+1}=0[/tex]) we end up with the required [tex]m[/tex] equations and the system has a unique solution.

There is no question this method doesn't work. My specific problem was with some of the details, why certain approximations are used instead of others that should result in more accurate results.
 
  • #6
pafcu said:
No, number of unknowns is not a problem.

Assuming [tex]u[/tex] is defined at m [tex]m[/tex] points (it has been discretized) we obviously end up with [tex]m[/tex] unknowns: [tex]u_i^{n+i}, 1 \leq i \leq m[/tex].
We get [tex]m-2[/tex] equations by plugging in [tex]i=2,3,4,...,m-1[/tex] into the equation in the original post. Plugging in [tex]i=1[/tex] doesn't work because [tex]u_{0}^{n+1}[/tex] required by the approximation isn't defined. Similarly you can't plug in the value [tex]i=m[/tex] because [tex]u_{m+1}^{n+1}[/tex] isn't defined.

Is all hope lost? No, because if we additionally specify two boundary conditions (say [tex]u_0^{n+1}=0[/tex] and [tex]u_m^{n+1}=0[/tex]) we end up with the required [tex]m[/tex] equations and the system has a unique solution.

There is no question this method doesn't work. My specific problem was with some of the details, why certain approximations are used instead of others that should result in more accurate results.

I should have realized that. Anyway, if you wanted to do the central difference for the derivative, you would also need [tex]u_{i}^{n+2}[/tex] if you didn't want to introduce any time legs. This means you would need twice as many equations. This would obviously slow down your algorithm but it might help accuracy.
 
  • #7
Yes, there is a reason. It would be unstable. A quick way to check. The instability mode is usually an alternating wave, +-+-+-+... It will be modified by a factor h say in each step. Substitute that in the algorithm - ie in successive timesteps, (+1,-1,+1...), (+h,-h,+h..), (+h^2,-h^2,+h^2..).
For a two-step algorithm like C-N, you'll get a linear equation for h, which will be <1. The wave decays, which is what you want and expect.
For a three-step algorithm like the one you propose, you'll get a quadratic in h. You'll quickly see that the product of the roots is 1. One is less than 1, and is also the one you expect. In fact, it's more accurate than C-N. But there is a spurious second root that is greater than 1. That causes such waves to grow exponentially, and that is why the solutions are unstable. It's a common problem with multi-step explicit algorithms.
 
  • #8
pafcu said:
Hi all,

A diffusion equation is of the form
[tex]\frac{\partial u}{\partial t}=k \frac{\partial^2u}{\partial x^2}[/tex]

Usually an equation like this seems to be solved numerically using the Crank-Nicolson method:
[tex]\frac{u_{i}^{n + 1} - u_{i}^{n}}{\Delta t} = \frac{k}{2 (\Delta x)^2}\left((u_{i + 1}^{n + 1} - 2 u_{i}^{n + 1} + u_{i - 1}^{n + 1}) + (u_{i + 1}^{n} - 2 u_{i}^{n} + u_{i - 1}^{n})\right)[/tex]

i.e the time-related part is approximated using a 2-point forward difference approximation and the space-part using an average of two 3-point central difference approximation. ([tex]u_{i}^{n}[/tex] means value of u at point i at timestep n, i.e u(i,n))

Is there a specific reason why a central difference is not used for the time? Since a central difference results in a smaller error, yet it's not usually used, what drawback does it have?

What about the central difference for the space-related part? Is there a reason for using the 3-point version instead of the 5-point version? Or is the simpler one usually used because it results in a tridiagonal matrix which is easier to solve?

edit: Clarified notation

The reason you go backward in time yet a central difference in space is precisely due to the ease involved in solving. If you went central in time as well, you would have an extra u evaluated at n-1, which would disturb the tri-diagonal structure of the RHS.

The reason the 3 point scheme is used for space is also for the same reason that you mentioned, for ease in computation. There are problems where a 5-point difference scheme is used, but then the physics is different there.

As a good exercise, you could try using a 3 point difference scheme for time and compare your results.

Another reason I came up with is that in a diffusion problem, the state at any time t would only depend on time t-1. Think about it, the temperature at a point now does not care what the temperature at that same point was say ten delta time back. All in all, its just physics!
 

1. What is the Crank-Nicolson method?

The Crank-Nicolson method is a numerical method for solving partial differential equations (PDEs) that uses a combination of the forward and backward Euler methods. It is a second-order accurate method that is unconditionally stable.

2. How does the Crank-Nicolson method work?

The Crank-Nicolson method works by discretizing the PDE into a system of linear equations and then solving for the unknown values at each time step. It uses a time-centered approach, taking into account both the current and previous time steps, to improve the accuracy of the solution.

3. What are the approximations used in the Crank-Nicolson method?

The Crank-Nicolson method uses a central difference approximation for the spatial derivative and a time-centered average for the time derivative. This leads to a tridiagonal system of equations that can be solved efficiently using various numerical methods, such as Gaussian elimination or the Thomas algorithm.

4. What are the advantages of using the Crank-Nicolson method?

The Crank-Nicolson method has several advantages, including its second-order accuracy, unconditional stability, and ability to handle a wide range of PDEs. It also produces more accurate results compared to other numerical methods, such as the forward Euler method.

5. Are there any limitations of the Crank-Nicolson method?

While the Crank-Nicolson method is a versatile and accurate numerical method, it does have some limitations. It may not be suitable for solving PDEs with highly non-linear terms, and it can be computationally expensive for large systems due to the need to solve a tridiagonal system of equations at each time step.

Similar threads

  • Differential Equations
Replies
8
Views
3K
  • Differential Equations
Replies
1
Views
816
  • Differential Equations
Replies
3
Views
1K
  • Differential Equations
Replies
23
Views
2K
Replies
1
Views
1K
  • Differential Equations
Replies
3
Views
2K
Replies
1
Views
1K
  • Differential Equations
Replies
4
Views
1K
  • Differential Equations
Replies
2
Views
1K
  • Differential Equations
Replies
1
Views
2K
Back
Top