#!/usr/bin/python
# -*- coding: utf-8 -*-
"""Find an intersection of two circles by hill climbing with adaptive step sizes

This kind of hill-climbing isn’t gradient descent; it’s gradient-free.
However, for the notion of a “step size” to be meaningful, it does
need to be possible to multiply the step by a scalar, so the step type
needs to be the vector type of a vector space (over ℝ).  In theory the
points where we evaluate the algorithm could be from the torsor of
that vector space, but this implementation uses the same step function
to generate the starting point, so they have to be vectors too.

"""
from __future__ import print_function, division
import sys

import numpy


def distance(p1, p2):
    return ((p1 - p2)**2).sum()**0.5

def circles(c1, r1, c2, r2):
    def badness(guess):
        return ((distance(guess, c1) - r1)**2 +
                (distance(guess, c2) - r2)**2)
    return badness

def climb(f, step):
    size, growth, here = 1, 2, step()
    cost = f(here)
    
    while True:
        yield here, cost  # This may be the same as last time

        where = size * step()
        new = here + where
        new_cost = f(new)
        if new_cost > cost:
            where = -1 * where
            new = here + where
            new_cost = f(new)

        if new_cost > cost:
            size /= 2
            continue

        where *= growth
        grown_new = here + where
        grown_new_cost = f(grown_new)
        if grown_new_cost < new_cost:
            new, new_cost = grown_new, grown_new_cost
            size *= growth
        else:
            growth = 1/growth

        here, cost = new, new_cost

if __name__ == '__main__':
    if sys.argv[1] == '--flagellate':
        sys.argv[1:2] = []
        import flagellate
        climb = flagellate.climb

    x1, y1, r1, x2, y2, r2 = map(float, sys.argv[1:7])
    n = int(sys.argv[7]) if len(sys.argv) > 7 else 128
    for i, (here, cost) in enumerate(climb(circles((x1, y1), r1, (x2, y2), r2),
                                           lambda: numpy.random.rand(2) - 0.5)):
        print(i, here, cost)
        if i >= n:
            break