#!/usr/bin/env python3

import os
import sys
import numpy as np
import phoebe

from astropy import units as u

deg = np.pi/180.0
atol = 4.0e-16

def test_xyz_2body(verbose=False):
    """
    Test of b.default_binary().
    Non-equal masses (q = 0.1).
    Integration spanning 2 d.
    Using xyz.py integrator.

    """

    # to re-compute default_binary()
    phoebe.conf.devel_on()

    dir_ = os.path.dirname(os.path.realpath(__file__))

    b2 = phoebe.default_binary()

    b2.set_value('q@binary', value=1.0)
    b2.set_value('requiv@primary', value=2.0*u.solRad)
    b2.set_value('requiv@secondary', value=2.0*u.solRad)

    b2.add_dataset('mesh', compute_times=[0.25], columns=['visibilities'], dataset='mesh01')

    b2.add_compute(dynamics_method='xyz', integrator='ias15', stepsize=1.0e-4, epsilon=1.0e-12, geometry='hierarchical')
    b2.run_compute(compute='phoebe02')

    fig, plt = b2.plot(fc='visibilities', ec='None', y='vs', show=True)
    plt.savefig("rebound.png", dpi=300)
    sys.exit(1)

if __name__ == "__main__":
    logger = phoebe.logger(clevel='INFO')
#    logger = phoebe.logger(clevel='DEBUG')

    test_xyz_2body(verbose=True)