This directory contains the code that was used to trace the coexistence curves as seen in figure 2. This is an older version of the Monte Carlo code, so it is a little less user friendly than the free energy version.  Here is the complete work-flow required to obtain the density histograms that can be used to trace the coexistence curve.

1. There is no makefile. Rather, the code is compiled with
	mpicc gc_mc_pt.c -o gc_mc_pt.out -lm

2. For each activity of interest, create a folder called activity0, activity1... (see the example input decks).

3. The code performs replica exchange in temperature, so prepare each activity folder with a deck of options files (see example), named options0.txt, options1.txt... These should be identical except for the value of 'epsilon'. The options files are fully annotated, you should be able to simply edit the examples for your own needs.

4. Copy the executable into each activity folder and run it with one process for each each options file.

6. In the auxiliary code folder, compile histogramsPT.c (for the sphere) or histogramsPT_FLT.c (for the plane with
	cc histogramsPT.c -o histogramPT.out
Copy this, hist.sh and sort.sh into the outer folder (which contains activity0, activity1...).

7. Run sort.sh. This simply sorts the data from each simulation rank so that all the data for each temperature is together. This creates folders containing the sorted data called data0, data1... You may have to edit sort.sh so it loops over the correct number of activities (outer loop) and temperatures (inner loop).

8. Run hist.sh. This runs histogramsPT.out on each data-set and collects the results in a folder called 2Dhist. If you are running the simulation on a plane, edit hist.sh to reference histogramsPT_FLT.out. Please also ensure that it loops over the correct number of activity folders.

	OUTPUT: in 2Dhist there will be a file called density_average.csv which contains the average density of particles (third column) at each epsilon (first column) and activity (second 	column).

9. Do multiple histogram reweighting. Again in the auxiliary code folder, compile histS and weight2D_sphere.c with
		cc hist_sphere2.c nrutil.c -o histS -lm
		cc weight2Dsphere.c -o weight2Dsphere.out -lm
Despite the names, these work on both the plane and the sphere. Copy both into the 2Dhist folder. Run
		./histS histograms.2D

10. Generate density histograms. To obtain the interpolated density histogram for any epsilon and activity within the range of the simulation, run
	 		./weight2Dsphere.out <new epsilon> <new activity> log_omegNM.dat

	OUTPUT: newhist.dat is the estimated histogram at this point in parameter space. Column 1 is the density and column 2 is the frequency.

Note: If the histogram is very noisy or has a handful of very high, 1 bin wide peaks, try removing lines from histograms.2D so that it only refers to the data points closest to the target epsilon and activity. Then rerun histS and weight2Dsphere. We find it often works best with only ~10 input histograms.
