Tutorial A12

1 How does logging behave when you rerun the script without restarting the kernel? What happens if you delete the log-file before you rerun?

2 Play around with the datefmt keyword when configuring the logger to omit displaying year and months.

3 Below you find a script that propagates a particle in a harmonic potential numerically (using the euler integrator). You do not need to know about details of the scientific background yet. Just note so much: The script monitors the total energy of the particle at every propagation step. Whenever the change in energy is larger than 0.3, the velocity should be rescaled (multiplied) by a factor of 0.9. This part is not yet included in the script below. If this happens a message should be written to a log file, stating the old and the new speed. Your task is it, to setup a logger and insert a logging statement at the right position in the code.

Advanced: If you are already a bit familar with the background of task 3. Try to implement the euler integrator yourself (and ignore the definition of the propagate() function below.

4 Plot the result of the simulation as velocity vs. position.

# Example code for task 3

def force(x):
    """Compute the instantaneous force acting on the particle

    Args:
        x: Current position
        x0 (global): Position of the minimum of the potential
        k (global): Force constant
    """

    return -k * (x - x0)


def e_pot(x):
    """Compute the potential energy of the particle

    Args:
        x: Current position
        x0 (global): Position of the minimum of the potential
        k (global): Force constant
    """

    return 0.5 * k * (x - x0)**2


def e_kin(v):
    """Compute the kinetic energy of the particle

    Args:
        v: Particle velocity
        m (global): Particle mass
    """

    return 0.5 * m * v**2


def propagate(x, v):
    """Compute the position of the particle one timestep further

    Args:
        x: Current position
        v: Current velocity
        dt (global): Propagation timestep
        m (global): Particle mass

    Returns:
        New position and velocity
    """

    x_new = x + v*dt + force(x)/(2*m) * dt**2
    v_new = v + force(x)/m * dt

    return x_new, v_new


# Preparing lists to store the computed trajectories
x_t = []  # Position
v_t = []  # Velocity
e_t = []  # Energy

# Particle/Simulation parameters
x0 = 0  # Position of the potential minimum
k = 1.  # Force constant
m = 1.  # Particle mass

dt = 1e-1  # Simulation timestep
T = 10000  # Total number of steps

x = 1  # starting position
v = -2  # starting velocity

# Put them in the output lists
x_t.append(x)
v_t.append(v)
e_t.append(e_pot(x) + e_kin(v))

# Simulate
for i in range(T):
    # Update position and velocity
    x, v = propagate(x, v)
    x_t.append(x)
    v_t.append(v)

    # Update total energy
    e_t.append(e_pot(x) + e_kin(v))