First, download, test and install mc_sim. Follow the instructions on its GitHub. Now Monte Carlo simulations can be run multiple ways:

1. On a single CPU on a local computer: ./mc_sim -p ParameterFile -o OutputName
2. On multiple CPUs on a local computer: mpirun -n 4 mc_sim -p ParameterFile -o OutputName
3. On a cluster, see General Cluster Usage

For using a cluster, usually a job management system will be used. On a local computer, simulations can be started manually in the terminal, but for longer runs or many simulations it can be helpful to write a python wrapper script:

import subprocess
import multiprocessing
import numpy as np
 
rho_min = 0.1
rho_max = 0.8
rho_step = 0.1
 
N = 800
 
# formatted string
parameters = """
#== General_Parameters ==========================================
Initialization [type/filename]: fluid
BoxDimension [Lx,Ly,Lz]: X X X
BoxGeometry [type]: pbc
CellList [cutoff]: 3.0
MovesPerCycle [cycle]: 10000
NumParticleTypes [NumTypes]: 1
NumParticles [num_A,num_B,num_C,..]: Y
Temperature [temp]: Z
#== Potentials ==================================================
# LJ parameters  - sigma, cutoff, epsilon
Potential [type,type,name,parameters]: A A LJ 1.0 3.0 1.0
 
#== Output ======================================================
PrintSnapshot [f_cy]: 1
PrintEnergy [f_cy]: 1
PrintPressure [f_cy]: 1
 
#== Monte_Carlo_Move_Set ========================================
CanMove [f_A,f_B,f_C,..]: 1
"""
 
cmd_list=[]
#iterate over all parameters 
for T in [0.9,2.0]:
    for rho in np.arange(rho_min,rho_max,rho_step):
        V = N/rho
        L = (V)**(1/3.)
        # modify the formatted string to replace the parameters which change
        p = parameters.replace('BoxDimension [Lx,Ly,Lz]: X X X', 'BoxDimension [Lx,Ly,Lz]: %1.3f %1.3f %1.3f'%(L,L,L))
        p = p.replace('NumParticles [num_A,num_B,num_C,..]: Y','NumParticles [num_A,num_B,num_C,..]: %d'%(N))
        p = p.replace('Temperature [temp]: Z','Temperature [temp]: %1.3f'%(T))
        param_file = open('./Input_%1.3f_%1.3f.par'%(rho,T),'w')
        param_file.write(p)
        param_file.close() # important to not leave open file handles around
        cmd_list.append(["./mc_sim", "-p", "Input_%1.3f_%1.3f.par"%(rho,T), "-o", "T_%1.3f_L_%1.3f_N_%d"%(T,L,N),"-W", "00:05:00"])  #this just appends what to do to a list 
 
def worker(cmd):
    print("spawn simulation ",cmd)
    p = subprocess.Popen(cmd) # this executes the actual simulation
    p.wait() # wait for result 
 
pool = multiprocessing.Pool( processes = 4 ) # run 4 in parallel
results =[pool.apply_async(worker, [cmd]) for cmd in cmd_list]
ans = [res.get() for res in results]

This script runs a bunch of LJ NVT simulations each with one CPU at different state points. The temperature and density is changed for each, and four are run at the same time. Each simulation gets it's own parameter file and the output filename has the parameters in it, which is useful for plotting and sorting the data afterwards. Each simulation has a time limit of 5 minutes. If only one simulation is supposed to run at a time, subprocess.Popen() can be used in a simple for-loop with subprocess.wait(). MPI jobs can also be run by modifying the command cmd which is called by the subprocess package.

Useful packages are subprocess and multiprocessing.