  The equilibrium phase behavior of microphase-forming systems is notoriously
difficult to obtain because of the extended metastability of their modulated
phases. In this paper we present a systematic simulation methodology for
studying layered microphases and apply the approach to two prototypical
lattice-based systems: the three-dimensional axial next-nearest-neighbor Ising
(ANNNI) and Ising-Coulomb (IC) models. The method involves thermodynamically
integrating along a reversible path established between a reference system of
free spins under an ordering field and the system of interest. The resulting
free energy calculations unambiguously locate the phase boundaries. The simple
phases are not observed to play a particularly significant role in the devil's
flowers. With the help of generalized order parameters, the
paramagnetic-modulated critical transition of the ANNNI model is also studied.
We confirm the XY universality of the paramagnetic-modulated transition and its
isotropic nature. Interfacial roughening is found to play at most a small role
in the ANNNI layered regime.
