#!/usr/bin/python
import matplotlib.pyplot as plt
import numpy as np

def copt(N):
    n=np.ceil(np.sqrt(N/2.))
    k=2*n-1
    return 2*N*k + k*(N*1./n)**2


fig = plt.figure()
ax = fig.add_subplot(1,1,1)
N=np.arange(1,10**5,1)
ax.plot(N, N**2,label="C_N")
ax.plot(N,N,label="X_N")
ax.plot(N,1./N,label="X_N/C_N")
ax.set_xscale('log')
ax.set_yscale('log')
#ax.axis((1,10**5,1/100.,10**10))
ax.legend()
fig.show()

raw_input('press return to continue')

fig = plt.figure()
ax = fig.add_subplot(1,1,1)
N=np.arange(1,10**5,1)
ax.plot(N, N**2,label="C_N")
ax.plot(N,copt(N),label="C_N,opt")
ax.plot(N,N**2*1./copt(N),label="C_N/C_N,opt")
ax.set_xscale('log')
ax.set_yscale('log')
#ax.axis((1,10**5,1/100.,10**10))
ax.legend()
fig.show()

raw_input('press return to continue')
