diff --git a/integration_tests/benchmark/configurationA/README.md b/integration_tests/benchmark/configurationA/README.md index ba4db84..8d521b3 100644 --- a/integration_tests/benchmark/configurationA/README.md +++ b/integration_tests/benchmark/configurationA/README.md @@ -1,11 +1,16 @@ -# Configuration A +# Configuration A - Particles in the NVT ensemble -Particle in the NVT ensemble. +## Parameters -nmb_1= 50 # Define atom number -sig_1 = 3 * ureg.angstrom # Define LJ parameters (sigma) -eps_1 = 0.1 * ureg.kcal/ureg.mol # Define LJ parameters (epsilon) -mss_1 = 10 * ureg.gram/ureg.mol # Define atom mass -L = 20 * ureg.angstrom # Define box size -rc = 2.5 * sig_1 # Define cut_off -T = 300 * ureg.kelvin # Pick the desired temperature +nmb_1= 50 # Define atom number +sig_1 = 3 * ureg.angstrom # Define LJ parameters (sigma) +eps_1 = 0.1 * ureg.kcal/ureg.mol # Define LJ parameters (epsilon) +mss_1 = 10 * ureg.gram/ureg.mol # Define atom mass +L = 20 * ureg.angstrom # Define box size +rc = 2.5 * sig_1 # Define cut_off +T = 300 * ureg.kelvin # Pick the desired temperature + +## Measurements + +Epot = -3.532 ± 0.021 kcal/mol +press = 298.069 ± 0.832 atm diff --git a/integration_tests/benchmark/configurationA/nvt.lmp b/integration_tests/benchmark/configurationA/nvt.lmp index 2613810..f29f27f 100644 --- a/integration_tests/benchmark/configurationA/nvt.lmp +++ b/integration_tests/benchmark/configurationA/nvt.lmp @@ -1,6 +1,8 @@ -variable dump equal 5000 -variable thermo equal 5000 -variable steps equal 100000 +variable dump equal 50000 +variable thermo equal 50000 +variable steps equal 1000000 +variable eqs equal 100000 +variable av equal ${dump}/10 variable nmb_1 equal 50 # Define atom number variable sig_1 equal 3 # Define LJ parameters (sigma) @@ -34,16 +36,17 @@ timestep 0.25 thermo ${thermo} dump mydmp all custom ${dump} dump.lammpstrj id type x y z vx vy vz -run ${steps} # equilibration +run ${eqs} # equilibration variable Epot equal pe variable Ekin equal ke variable Etot equal v_Epot+v_Ekin variable pressure equal press variable temperature equal temp -fix myat1 all ave/time ${dump} 1 ${dump} v_Epot file Epot.dat -fix myat2 all ave/time ${dump} 1 ${dump} v_Ekin file Ekin.dat -fix myat3 all ave/time ${dump} 1 ${dump} v_Etot file Etot.dat -fix myat4 all ave/time ${dump} 1 ${dump} v_pressure file pressure.dat -fix myat5 all ave/time ${dump} 1 ${dump} v_temperature file temperature.dat +fix myat1 all ave/time 10 ${av} ${dump} v_Epot file Epot.dat +fix myat2 all ave/time 10 ${av} ${dump} v_Ekin file Ekin.dat +fix myat3 all ave/time 10 ${av} ${dump} v_Etot file Etot.dat +fix myat4 all ave/time 10 ${av} ${dump} v_pressure file pressure.dat +fix myat5 all ave/time 10 ${av} ${dump} v_temperature file temperature.dat + run ${steps} diff --git a/integration_tests/benchmark/configurationA/update_README.py b/integration_tests/benchmark/configurationA/update_README.py new file mode 100644 index 0000000..6b02d07 --- /dev/null +++ b/integration_tests/benchmark/configurationA/update_README.py @@ -0,0 +1,44 @@ +import numpy as np +import re + +# Function to calculate the average and standard error of the second column in a data file +def calculate_statistics(filename): + try: + data = np.loadtxt(filename, usecols=1) # Load the second column + mean = np.mean(data) + std_error = np.std(data, ddof=1) / np.sqrt(len(data)) # Standard error of the mean + return mean, std_error + except Exception as e: + print(f"Error reading {filename}: {e}") + return None, None + +# Filenames +config_file = "README.md" # Your configuration file +energy_file = "Epot.dat" +pressure_file = "pressure.dat" + +# Calculate statistics +average_epot, error_epot = calculate_statistics(energy_file) +average_press, error_press = calculate_statistics(pressure_file) + +if average_epot is None or average_press is None: + print("Failed to calculate averages or errors. Exiting.") + exit(1) + +# Update the configuration file +try: + with open(config_file, 'r') as file: + config_lines = file.readlines() + + with open(config_file, 'w') as file: + for line in config_lines: + if re.match(r"^Epot\s*=", line): + file.write(f"Epot = {average_epot:.3f} ± {error_epot:.3f} kcal/mol\n") + elif re.match(r"^press\s*=", line): + file.write(f"press = {average_press:.3f} ± {error_press:.3f} atm\n") + else: + file.write(line) + + print(f"Updated {config_file} with Epot = {average_epot:.3f} ± {error_epot:.3f} kcal/mol and press = {average_press:.3f} ± {error_press:.3f} atm") +except Exception as e: + print(f"Error updating {config_file}: {e}") diff --git a/integration_tests/test_monte_carlo_move.py b/integration_tests/test_monte_carlo_move.py index 0548b37..74c58f8 100644 --- a/integration_tests/test_monte_carlo_move.py +++ b/integration_tests/test_monte_carlo_move.py @@ -25,12 +25,12 @@ def setUp(self): L = 20 * ureg.angstrom # Define box size rc = 2.5 * sig_1 # Define cut_off T = 300 * ureg.kelvin # Pick the desired temperature - displace_mc = sig_1/2 # choose the displace_mc + displace_mc = sig_1/4 # choose the displace_mc # Initialize the MonteCarlo object self.mc = MonteCarlo( ureg = ureg, - maximum_steps = 10000, + maximum_steps = 50000, thermo_period = 1000, dumping_period = 1000, number_atoms = [nmb_1],