creator = particlesim.utils.config_parser.ProblemCreator("config/200_particle_nacl_rand.cfg")

system_config = creator.generate_problem()

plot_systemconfig(system_config)

sampler = particlesim.api.Sampler(system_config)

"""
Explanation: Config File:
[general]
box-size = 12
[particle_class_1]
label = Natrium
type = Na
charge = 1
distribution = uniform
number = 100
[particle_class_2]
label = Chlor
type = Cl
charge = -1
distribution = uniform
number = 100
Generate system configuration from cfg-file with "problem creator" and a related sampler
End of explanation
"""


iteration_number = 20000
step = 0.002
beta = [0.09,1]

"""
Explanation: Global parameter for simulated annealing
End of explanation
"""


%%time
traj,pot = sampler.metropolis_sa(iteration_number=iteration_number,step=step,beta=beta)

plot_nacl(traj=traj,left=0,right=12,num_na=100,traj_sample=-1)

beta_values = 1.0 / np.linspace(1.0 / beta[1], 1.0 / beta[0], iteration_number)[::-1]

fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(pot)
ax.set_xlabel(r"Iteration", fontsize=15)
ax.set_ylabel(r"Potential [kcal/mol]", fontsize=15)
ax.set_ylim([pot[-1]-1000,-10000])
ax2 = ax.twinx()
ax2.plot(beta_values,c='r')
ax2.set_ylabel(r"$\beta$ [mol/kcal]",fontsize=15)
ax.set_title('Potential curve over iterationsteps')
print(pot[-1])

plot_beta(beta_values,T=True)

"""
Explanation: Call metropolis simulated annealing function from sampler
End of explanation
"""


beta_list=[[0.1,0.9]]+[[0.3,1.]]
iter_list=[iteration_number//2,iteration_number//2]
#beta_list=[beta]
#iter_list= [iteration_number]
beta_many = multiple_beta(beta_list,iter_list, linspace=True)
plot_beta(beta_many,T=True)

%%time
traj_many, pot_many = sampler.metropolis_sa(iteration_number=iteration_number,step=step,beta=beta_many)

fig = plt.figure()
ax = fig.add_subplot(111)
ax.plot(pot_many)
ax.set_xlabel(r"Iteration", fontsize=15)
ax.set_ylabel(r"Potential [$kcal$/$mol$]", fontsize=15)
ax.set_ylim([pot_many[-1]-1000,20000])
ax2 = ax.twinx()
ax2.plot(beta_many,c='r')
ax2.set_ylabel(r"$\beta$ [mol/kcal]",fontsize=15)
print(pot_many[-1])

plot_nacl(traj=traj_many,left=0,right=12,num_na=100,traj_sample=-1)

system_config_mc = creator.generate_problem()

system_config_mc.xyz = traj[-1]

sampler_mc = particlesim.api.Sampler(system_config_mc)

"""
Explanation: create beta values for double simulated annealing of system configuration
End of explanation
"""


%%time
beta_mc = particlesim.utils.conversion.kelvin_to_beta(200)
iteration_mc = 2500
step_mc = 0.001
traj_mc_200, pot_mc_200 = sampler_mc.metropolis(iteration_number=iteration_mc,step = step_mc, beta= beta_mc) 

beta_mc = particlesim.utils.conversion.kelvin_to_beta(800)
traj_mc_800, pot_mc_800 = sampler_mc.metropolis(iteration_number=iteration_mc,step = step_mc, beta= beta_mc) 

beta_mc = particlesim.utils.conversion.kelvin_to_beta(1200)
traj_mc_1200, pot_mc_1200 = sampler_mc.metropolis(iteration_number=iteration_mc,step = step_mc, beta= beta_mc)

fig, axes = plt.subplots(1, 2, figsize=(10, 4))
#fig = plt.figure()
#ax = fig.add_subplot(111)
axes[0].plot(pot_mc_200, c='r', label = r"200 $K$")
axes[0].plot(pot_mc_800, c='b', label = r"800 $K$")
axes[0].plot(pot_mc_1200,c='g', label = r"1200 $K$")
axes[0].set_xlabel(r"Iteration", fontsize=15)
axes[0].set_ylabel(r"Potential [$kcal$/$mol$]", fontsize=15)
axes[0].legend()

#ax.set_ylim([pot_mc[-1]-200,0])

#ax2 = fig.add_subplot(222)
def create_hist(pot_mc):
    hist, edges = np.histogram(pot_mc, bins=50)
    rho = hist.astype(np.float64) / np.sum(hist)
    return edges,rho
edges200, rho200 = create_hist(pot_mc_200)
edges800, rho800 = create_hist(pot_mc_800)
edges1200, rho1200 = create_hist(pot_mc_1200)

axes[1].plot(edges200[1:], rho200, c='r', label = r"200 $K$")
axes[1].plot(edges800[1:], rho800, c='b', label = r"800 $K$")
axes[1].plot(edges1200[1:], rho1200, c='g', label = r"1200 $K$")

axes[1].set_xlabel(r"$E$ [$kcal$/$mol$]", fontsize=15)
axes[1].set_ylabel(r"$\rho(E)$", fontsize=15)
axes[1].legend()
fig.tight_layout()

"""
Explanation: Do plain metropolis monte carlo with last configuration from simulated annealing with three different temperatures

T_1 = 200 Kelvin
T_2 = 800 Kelvin
T_3 = 1200 Kelvin
End of explanation
"""


