import numpy as np
import matplotlib.pyplot as plt

plt.style.use('ggplot')  

thermostats = {
    'be': {'color': '#2c7bb6', 'name': 'Berendsen Thermostat'},
    'nh': {'color': '#d7191c', 'name': 'Nose-Hoover Thermostat'},
    'sd': {'color': '#fdae61', 'name': 'Stochastic Dynamics'}, 
    'vr': {'color': '#1a9641', 'name': 'Velocity Rescale'}
}

for thermo_code, params in thermostats.items():
    t = np.loadtxt(f'et_{thermo_code}_potential_en.xvg')[:,0]
    U = np.loadtxt(f'et_{thermo_code}_potential_en.xvg')[:,1]
    K = np.loadtxt(f'et_{thermo_code}_kinetic_en.xvg')[:,1]
    
    fig, ax = plt.subplots(figsize=(10, 5))
    
    ax.plot(t, U, color=params['color'], alpha=0.8, 
            linewidth=2, label='Potential Energy')
    ax.plot(t, K, color='#333333', linestyle='--', 
            alpha=0.7, linewidth=1.5, label='Kinetic Energy')
    
    ax.set_title(f"{params['name']}", pad=15, fontsize=13)
    ax.set_xlabel('Simulation Time (ps)', labelpad=10)
    ax.set_ylabel('Energy (kJ/mol)', labelpad=10)
    ax.legend(framealpha=0.9, loc='upper right')
    
    ax.axhline(np.mean(U), color=params['color'], linestyle=':', linewidth=1)
    ax.axhline(np.mean(K), color='#333333', linestyle=':', linewidth=1)
    
    plt.tight_layout()
    plt.show()

import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import gaussian_kde

plt.style.use('seaborn-v0_8-whitegrid')

colors = ['#3498db', '#e74c3c', '#2ecc71', '#9b59b6']

for idx, (code, name) in enumerate([('be', 'Berendsen'), 
                                  ('nh', 'Nose-Hoover'),
                                  ('sd', 'Stochastic Dynamics'),
                                  ('vr', 'Velocity Rescale')]):
    
    bond_lengths = np.loadtxt(f'bond_{code}.xvg')[:,1] * 10  
    
    fig, ax = plt.subplots(figsize=(8, 5))
    
    ax.hist(bond_lengths, bins=15, color=colors[idx], 
           alpha=0.7, edgecolor='white', density=True)
    
    kde = gaussian_kde(bond_lengths)
    x_vals = np.linspace(min(bond_lengths), max(bond_lengths), 200)
    ax.plot(x_vals, kde(x_vals), color=colors[idx], linewidth=2)
    
    mean_val = np.mean(bond_lengths)
    ax.axvline(mean_val, color='#34495e', linestyle='--', linewidth=1.5,
              label=f'Mean: {mean_val:.2f} Å')
    
    ax.set_title(f'C-C Bond Length Distribution\n({name} Thermostat)', pad=15)
    ax.set_xlabel('Bond Length (Å)', labelpad=10)
    ax.set_ylabel('Probability Density', labelpad=10)
    ax.legend()
    
    plt.tight_layout()
    plt.show()

Изучение работы методов контроля температуры в молекулярной динамике¶

Файл топологии этана¶

Кинетическая и потенциальная энергиии¶

Распределение длин С-С связей¶