2 difference from LCDM results and the bias between CDM and baryons; § III D briefly describes the essential mod- ifications one must bear in mind when identifying viri- alized halos from the simulation outputs and shows the mass functions for our models; § III E we pick out two halos from our simulation box and analyzes their total internal profiles, as well as their baryonic/CDM density profiles. We finally summarize in § IV. In this section we first describe the method to simulate structure formation with two differently coupled matter species and the appropriate equations to be used. Those equations for a single matter species have been discussed in details previously in [23], but the inclusion of different matter species requires further modifications and we list all these for completeness. instead simplifies the simulation by assuming that the scalar field’s effect is to change the value of the gravi- tational constant, and presenting an justifying argument for such an approximation. As discussed in [23, 24], this approximation is only good in certain parameter spaces and for certain choices of the scalar field potential, and therefore full simulations tare needed to study the scalar field behaviour more rigorously. Recently there have also appeared N -body simulations of the f ( R ) gravity model [25, 26], which do solve the scalar degree of freedom explicitly. However, embedded in the f ( R ) framework there are some limitations in the generality of these works. As a first thing, f ( R ) gravity model (no matter what the form f is) only corresponds to the couple scalar field models for a specific value of coupling strength [27]. Second, in f ( R ) models the cor- rection to standard general relativity is through the mod- ification to the Poisson equation and thus to the gravi- tational potential as a whole [25], while in the coupled scalar field models we could clearly separate the scalar fifth force from gravity and analyze the former directly [23]. Also, in f ( R ) models, as well as the scalar-tensor theories, the coupling between matter and the scalar field is universal (the same to dark matter and baryons), while in the couple scalar field models it is straightforward to switch on/off the coupling to baryons and study the ef- fects on baryonic and dark matter clusterings respectively (as we will do in this paper). Correspondingly, the gen- eral framework of N -body simulations in coupled scalar field models could also handle the situation where the chameleon effect is absent and/or scalar field only couples to dark matter, and thus provide a testbed for possible violations of the weak equivalence principle. The Lagrangian for our coupled scalar field model is where R is the Ricci scalar, κ = 8 πG with G Newton’s constant, ϕ is the scalar field, V ( ϕ ) is its potential energy and C ( ϕ ) its coupling to dark matter, which is assumed to be cold and described by the Lagrangian L $_{CDM}$. L S in- cludes all other matter species, in particular our baryons . The contribution from photons and neutrinos in the N - body simulations (for late times, i.e. , z ∼ O (1)) is negli- gible, but should be included when generating the matter power spectrum from which the initial conditions for our N -body simulations are obtained (see below). The dark matter Lagrangian for a point-like particle with bare mass m$_{0}$ is In this paper we shall go beyond [23] and consider the case where the chameleon scalar field couples differently to different species of matter. To be explicit, we consider two matter species, and let one of them have no coupling to the scalar field. Because it is commonly believed that normal baryons, being observable in a variety of experi- ments, should have extremely weak (if any) coupling to scalar fields, we call the uncoupled matter species in our simulation ”baryons”. It is however reminded here that this matter species is not really baryonic in the sense that it does not experience normal baryonic interactions. The inclusion of true baryons will make the investigation more complicated and is thus beyond the scope of the present work. The paper is organized as follows: in § II we list the essential equations to be implemented in the N -body sim- ulations and describe briefly the difference from normal LCDM simulations. § III is the main body of the paper, in which § III A gives the details about our simulations, such as code description and parameter set-up, § III B displays some preliminary results for visualization, such as baryon/CDM distribution, potential/scalar field con- figuration and the correlation between the fifth force (for CDM particles) and gravity; § III C quantifies the nonlin- ear matter power spectrum of our model, especially the Also, because g$_{ab}$ ˙ x a 0 ˙ x b 0 ≡ g$_{ab}$u $^{a}$ub = 1 where u a is the four velocity of the dark matter particle, the Lagrangian could be rewritten as where y is the coordinate and x 0 is the coordinate of the centre of the particle. From this equation it can be easily derived that L $_{CDM}$( y ) = − m$_{0}$ √ − g δ ( y − x $_{0}$) √ g$_{ab}$ ˙ x a 0 ˙ x b 0 (2) T ab CDM = m$_{0}$ √ − g δ ( y − x $_{0}$) ˙ x a 0 ˙ x b $_{0}$. (3) L $_{CDM}$( y ) = − m$_{0}$ √ − g δ ( y − x $_{0}$) , (4) which will be used below. L = 1 2 [ R κ − ∇ $^{a}$ϕ ∇ $_{a}$ϕ ] + V ( ϕ ) − C ( ϕ ) L CDM + L $_{S}$(1) A. The Fundamental Equations II. THE EQUATIONS