# -*- coding: utf-8 -*-
"""
Created on Mon Jan  9 15:16:15 2017

@author: Mariusz
"""
import matplotlib.pyplot as plt
import numpy as np

def A(x): 
    return np.array([[1.,2.,3.,x],[4.,5.,x,6.],[7.,x,8.,9.],[x,10.,11.,12.]])

def f(y):
    return np.linalg.det(A(y)) - 1000.

# plot the function to "visually" bracket the root
X = np.linspace(9, 10, 100)
plt.plot(X, [f(z) for z in X])
plt.xlim(X[0], X[-1])
plt.grid(True)
plt.draw()

'''
Following computes approximate solution of f(x)=0, which is the bisection
method. Input: function f; a,b such that f(a)*f(b)<0, and tolerance tol
'''

tol = 1.0e-6
tol_fun = 1.0e-20
a = 9.
b = 10.

fa = f(a)
fb = f(b)
assert(fa*fb<0)
i = int(0)
while b-a > tol:
    i = i + 1
    c  = (b + a)/2.
    fc = f(c)
    print('Iter: %4d, a: %e, f(a): %e, c: %e, f(c): %e, b: %e, f(b): %e' % 
          (i, a, fa, c, fc, b, fb))
    if abs(fc) < tol_fun:
        break    
    if fa*fc < 0:
        b  = c
        fb = fc
    else:
        a  = c
        fa = fc

print('Backward error: %e' % fc)

plt.show()