Computational Physics - Department of Physics

310 10 Partial Differential Equations

y(i) = u(i) = func(delta_x*i);

}

// Boundary conditions (zero here)

y(n) = u(n) = u(0) = y(0);

// Matrix A,onlyconstants

a = c = - alpha;

b = 1 + 2*alpha;

// Time iteration

for(int t = 1; t <= tsteps; t++){

// here we solve the tridiagonal linear set of equations,

// see chapter 6

tridag(a, b, c, y, u, n+1);

// boundary conditions

u(0) = 0;

u(n) = 0;

// replace previous time solution with new

for(int i = 0; i <= n; i++){

y(i) = u(i);

}

// You may consider printing the solution at regular time intervals

.... //printstatements

} //endtime iteration

}

10.2.3Crank-Nicolson scheme

It is possible to combine the implicit and explicit methods in a slightly more general approach.
Introducing a parameterθ(the so-calledθ-rule) we can set up an equation

a(t)

t

g(x)

b(t)

x

ui− 1 ,j+ 1 ui,j+ 1 ui+ 1 ,j+ 1

ui,j

✲

✻

Fig. 10.2Calculational molecule for the implicit scheme.