Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

[Questions] How to report to stdout while running ReplicaExchangeSampler #737

Open
4doom4 opened this issue Jul 24, 2024 · 2 comments
Open

Comments

@4doom4
Copy link

4doom4 commented Jul 24, 2024

I am using MultiStateReporter ("file.nc", checkpoint_interval=interval) as the storage in a simulation from ReplicaExchangeSampler. This works to have the states as netcdf once the simulation is done. Is there a way to attach multiple reports to the simulation like in native openmm to report to stdout information on the progress of the simulation? I saw https://openmmtools.readthedocs.io/en/stable/devtutorial.html#id18 but I am not sure how to use this.

Thanks.

@ijpulidos
Copy link
Contributor

Some progress of the simulation should already be given in stdout by default. Can you share the code that you are using that's not reporting anything on stdout?

@4doom4
Copy link
Author

4doom4 commented Aug 22, 2024

I basically followed the REST2 protocol from here: https://openmm.github.io/openmm-cookbook/latest/notebooks/tutorials/Running_a_REST_simulation.html

The relevant part of the script:

print("Configuring REST2 protocol...")
# Set temperatures for each thermodynamic state
n_replicas = 16  # Number of temperature replicas
T_min = 300 * kelvin  # Minimum temperature (i.e., temperature of desired distribution)
T_max = 600 * kelvin  # Maximum temperature
temperatures = [T_min + (T_max - T_min) * (math.exp(float(i) / float(n_replicas-1)) - 1.0) / (math.e - 1.0)
                for i in range(n_replicas)]


# Create reference thermodynamic state
rest_state = RESTState.from_system(rest_system)
thermostate = ThermodynamicState(rest_system, temperature=T_min)
compound_thermodynamic_state = CompoundThermodynamicState(thermostate, composable_states=[rest_state])


# Create thermodynamic states
sampler_state =  SamplerState(coord.positions, box_vectors=rest_system.getDefaultPeriodicBoxVectors())
beta_0 = 1/(kB*T_min)
thermodynamic_state_list = []
sampler_state_list = []
for temperature in temperatures:
    # Create a thermodynamic state with REST interactions scaled to the given temperature
    beta_m = 1/(kB*temperature)
    compound_thermodynamic_state_copy = copy.deepcopy(compound_thermodynamic_state)
    compound_thermodynamic_state_copy.set_rest_parameters(beta_m, beta_0)
    thermodynamic_state_list.append(compound_thermodynamic_state_copy)

    # Generate a sampler_state with minimized positions
    context, integrator = cache.global_context_cache.get_context(compound_thermodynamic_state_copy)
    sampler_state.apply_to_context(context, ignore_velocities=True)
    openmm.LocalEnergyMinimizer.minimize(context)
    sampler_state.update_from_context(context)
    sampler_state_list.append(copy.deepcopy(sampler_state))


# Set up sampler
swap_attempts = 500
move = mcmc.LangevinDynamicsMove(timestep=0.002*picoseconds, n_steps=1000) # each move is 2 ps
simulation = ReplicaExchangeSampler(mcmc_moves=move, number_of_iterations=swap_attempts)

# Run repex
start, start_time = timer()
print(f"Start REST2 protocol at {start_time}...")
reporter = multistate.MultiStateReporter(f"{prefix}_rest2.nc", checkpoint_interval=5)

simulation.create(thermodynamic_states=thermodynamic_state_list,
                  sampler_states=sampler_state_list,
                  storage=reporter)

simulation.run()

end, end_time = timer()
print(f"REST2 protocol completed at {end_time}. Elapsed time: {end - start:.2f} seconds")

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants