apps170608_collision.binary2_timeseriesΒΆ

Plot the time series

""" Plot the time series"""

from irfpy.util import collision as coll
import numpy as np

import matplotlib.pyplot as plt
import matplotlib.patches

def main():
    d = 0.2e-10    # impact parameter in meter
    r = 1e-10    #  To define the maximum impact parameter, or radius of particles
    m1 = 2.5
    m2 = 4.8
    v = 100    # in m/s

    do_check(m1, v, m2, d, r, 'binary2_nominal.png')
    do_check(m1, v, m2, r * 0.99, r, 'binary2_far.png')
    do_check(m1, v, m2, r * 0.01, r, 'binary2_head.png')
    do_check(m1, v, m2 * 10000, d, r, 'binary2_toheavy.png')
    do_check(m1 * 10000, v, m2, d, r, 'binary2_tolight.png')


def do_check(m1, v, m2, d, r, figname):
    vx1, vy1, vx2, vy2 = coll.binary_collision_stationary(m1, v, m2, d, r)

    m1q = m1 ** (1/3.)
    m2q = m2 ** (1/3.)
    r1 = m1q * r / (m1q + m2q)
    r2 = m2q * r / (m1q + m2q)

    ### Let's see time step
    dt = r / v * 0.25

    l = np.sqrt(r ** 2 - d ** 2)

    body1 = matplotlib.patches.Circle((-l, d), radius=r1, color='r')
    body2 = matplotlib.patches.Circle((0, 0), radius=r2, color='b')

    fig = plt.figure()
    ax = fig.add_subplot(111)
    ax.add_patch(body1)
    ax.add_patch(body2)
    plt.plot(-l - np.arange(10) * dt * v, d + np.zeros([10]), 'k.')    # Original velocity vector
    plt.plot(-l + np.arange(10) * dt * vx1, d + np.arange(10) * dt * vy1, 'ko')    # New vlocity vector
    plt.plot(np.arange(10) * dt * vx2, np.arange(10) * dt * vy2, 'ko')

    ax.grid()

    ax.set_xlim(-r*2, r*2)
    ax.set_ylim(-r*2, r*2)
    ax.set_aspect(1)

    plt.savefig(figname)

if __name__ == '__main__':
    main()