apps130704_negativeion.app04_3_runmulti

import sys
import time
import datetime
import logging
logging.basicConfig()
logging.getLogger().setLevel(logging.DEBUG)


import numpy as np

from . import app04_1_initializer
from . import app04_2_simulator

def execmain(bsw, vsw, m, q, dt, n=1000, output=None):
    '''Main script'''

    if output == None:
        output = sys.stdout
    
    poslist = app04_1_initializer.position_creator2(maxsize=n)
    vellist = app04_1_initializer.direction_creator_rnd2(maxsize=n)
    n = min([poslist.shape[1], vellist.shape[1]])
    print('#Bsw=', bsw, file=output)
    print('#Vsw=', vsw, file=output)
    print('#M=', m, file=output)
    print('#Q=', q, file=output)
    print('#dt=', dt, file=output)
    
    poslist = poslist[:, :n] * 1e3
    vellist = vellist[:, :n] * np.random.uniform(100, 300, size=n) * 1e3

    tracer = app04_2_simulator.Tracer(bsw, vsw, m, q, dt)

    finalpos = []
    finalvel = []

    ta = 0
    for i, (pos, vel) in enumerate(zip(poslist.T, vellist.T)):
        if np.dot(pos, vel) < 0:   # Only upflowing velocity considered.
            vel = -vel

        t0 = time.time()
        print(i, ('%.2f ' * 3) % tuple(pos / 1000), ('%.2f ' * 3) % tuple(vel / 1000), end=' ', file=output)

        tlist, plist, thlist, vhlist = tracer.run(pos, vel)
        print('%.3f' % tlist[-1], ('%.2f '*3) % tuple(plist[-1, :] / 1000), ('%.2f ' * 3) % tuple(vhlist[-1, :] / 1000), len(tlist), file=output)

        ta += (time.time() - t0)
        if i % 100 == 99:
            avg = ta / (i + 1)
            left = datetime.timedelta(seconds=(n - i) * avg)
            
            logging.debug('time exec = %f (avg (%d). %f)' % (time.time() - t0, i + 1, ta / (i + 1)))
            logging.debug('.... Expected end %s' % (datetime.datetime.now() + left))


def main():
    # Parameter roughly on 2008-02-14 (DOY 35) from ACE
    bsw = np.array([3e-9, -3e-9, 0])    
    vsw = np.array([-600e3, 0, -40e3])

    m = 1.67e-27
    q = 1.60e-19

    dt = 0.005  # seconds.
    n = 2000000  # To examine (max)

    execmain(bsw, vsw, m, -q, dt, n, output=open('dat080214_n.txt', 'w'))
    execmain(bsw, vsw, m, q, dt, n, output=open('dat080214_p.txt', 'w'))

if __name__ == "__main__":
    main()