# Load the configuration
configuration = nlread("equilibrated_pma.hdf5")[0]

# Create and set the calculator
builder = OPLSPotentialBuilder()
calculator = SoftMatterCalculator(
    builder, use_lennard_jones_pme=True, assign_atom_types=True
)
configuration.setCalculator(calculator)

# Set the simulation
atomic_constraints = [FixCenterOfMass()]
method = NPTNoseHooverMonteCarlo(initial_velocity=ConfigurationVelocities())
simulation = SoftMatterDynamicsSimulation(
    configuration, method, atomic_constraints=atomic_constraints
)
simulation.setLogging(1000)
simulation.setTrajectory(1000, "pma_trajectory.hdf5")

# Measure the volume during the simulation
measurement = SoftMatterDynamicsMeasurements(volume=True, call_interval=10)

# Set the pressure profile
pressure_profile = SimulationPressureProfile(
    initial_value=0*MPa,
    initial_gradient=200*MPa
)
pressure_profile.addQuantityValue(0.1, 20*MPa, gradient=0*MPa)
pressure_profile.addQuantityValue(0.2, 20*MPa, gradient=200*MPa)
pressure_profile.addQuantityValue(0.3, 40*MPa, gradient=0*MPa)
pressure_profile.addQuantityValue(0.4, 40*MPa, gradient=200*MPa)
pressure_profile.addQuantityValue(0.5, 60*MPa, gradient=0*MPa)
pressure_profile.addQuantityValue(0.6, 60*MPa, gradient=200*MPa)
pressure_profile.addQuantityValue(0.7, 80*MPa, gradient=0*MPa)
pressure_profile.addQuantityValue(0.8, 80*MPa, gradient=200*MPa)
pressure_profile.addQuantityValue(0.9, 100*MPa, gradient=0*MPa)

# Run the simulation
simulation.simulate(100000, hook_functions=[measurement], profiles=[pressure_profile])