Rationale behind using Monte-Carlo simulations to obtain excat p-values
The problem is to obtain an exact p-value for the one-sample chi-square for small sample situations where the Normal approximation does not hold.
This can be done by generating random samples (groups) of particular sizes conditional on the sum of these samples having the same total sample size as the observed data (ie generating observations from a multinomial distribution). If this is done a large number of times, say 5000, the p-value is equal to the proportion of times that the observed chi-square value is greater than or equal to the chi-square value observed from the data. A consequence of this method is that a different p-value will be generated with each time this procedure is run although the differences should be negligible provided each run comprises of a large enough number of simulated random samples e.g. 5000.
An alternative strategy is to evaluate all possible combinations of groups whose sizes sum to the total sample size of the observed data and express the p-value as the proportion of times the computed chi-square values are greater than or equal to the one observed from the data.
-_Reference Edgington E., Onghena P. (2007) Randomization tests. 4th edition. Blackwell Publishing. (gives examples of computing p-values using exact tests for a variety of commonly occurring statistical tests which are useful for small sample sizes)