I am new to this community as well as to scientific programming. I programmed a simple 4th order Runge-Kutta for the 1-D Cahn-Hilliard Equation for some first simple calculations on pattern forming systems. It turned out to be extremely sensitive to the change of the size of the spatial grid. With 100 grid elements the time step needs to be like 100 times smaller in order to converge ... This is the code i used:
import numpy as np
import matplotlib.pyplot as plt
# Parameters
k_c = 1.
eps = 0.5
L=12. # Domainsize
size = 50 # Gridsize
dx = L/size
T = 100 # total time
dt = .001 # time step
n = int(T/dt)
# Random Initial Conditions
U = np.random.normal(0,0.01,size) # Normalverteilte Zufallswerte um 0
# Laplace-Operator
def laplacian(Z):
Zleft = np.roll(Z,1)
Zright = np.roll(Z,-1)
Zcenter = Z
return (Zright + Zleft -2*Zcenter)/dx**2
def DoubleLaplacian(Z):
Z2left = np.roll(Z,2)
Z2right = np.roll(Z,-2)
Zleft = np.roll(Z,1)
Zright = np.roll(Z,-1)
Zcenter = Z
return (Z2left - 4* Zleft + 6*Zcenter - 4*Zright + Z2right)/dx**4
# Equation Body
def equation(U):
return laplacian(U**3)-2*eps/(k_c**2)*laplacian(U)-eps*(k_c**4)*DoubleLaplacian(U)
# Runge-Kutta-Step
def RK4Step(U, dt, equation):
k1 = dt*equation(U)
U_temp = U + k1/2; k2 = dt*equation(U_temp);
U_temp = U + k2/2; k3 = dt*equation(U_temp);
U_temp = U + k3; k4 = dt*equation(U_temp);
U += (k1 + 2*k2 + 2*k3 + k4)/6
# Integration
for i in range(n):
RK4Step(U, dt, equation)
Since I am new to this kind of stuff I wanted to ask whether there is a conceptual fault in my thinking or coding or it is simply due to the nonlinearity. I would be as well pretty happy to get some hints how to improve the performance of the code.