  We propose an algorithm to obtain numerically approximate solutions of the
direct Ising problem, that is, to compute the free energy and the equilibrium
observables of spin systems with arbitrary two-spin interactions. To this
purpose we use the Adaptive Cluster Expansion method, originally developed to
solve the inverse Ising problem, that is, to infer the interactions from the
equilibrium correlations. The method consists in iteratively constructing and
selecting clusters of spins, computing their contributions to the free energy
and discarding clusters whose contribution is lower than a fixed threshold. The
properties of the cluster expansion and its performance are studied in detail
on one dimensional, two dimensional, random and fully connected graphs with
homogeneous or heterogeneous fields and couplings. We discuss the differences
between different representations (Boolean and Ising) of the spin variables.
