Small scale problems of the $\Lambda$CDM model: a short review

The $\Lambda$CDM model, or concordance cosmology, as it is often called, is a paradigm at its maturity. It is clearly able to describe the universe at large scale, even if some issues remain open, such as the cosmological constant problem , the small-scale problems in galaxy formation, or the unexplained anomalies in the CMB. $\Lambda$CDM clearly shows difficulty at small scales, which could be related to our scant understanding, from the nature of dark matter to that of gravity; or to the role of baryon physics, which is not well understood and implemented in simulation codes or in semi-analytic models. At this stage, it is of fundamental importance to understand whether the problems encountered by the $\Lambda$DCM model are a sign of its limits or a sign of our failures in getting the finer details right. In the present paper, we will review the small-scale problems of the $\Lambda$CDM model, and we will discuss the proposed solutions and to what extent they are able to give us a theory accurately describing the phenomena in the complete range of scale of the observed universe.


Introduction
Despite the ΛCDM model being successful, according to the largest part of the cosmology community, in describing the formation and evolution of the large scale structure in the Universe, the state of the early Universe and the abundance of different forms of matter and energy [1][2][3][4][5], its predictive power -already checked against new discoveries (e.g., lensing of the CMB [6,7], B-mode polarisation [8] the kinetic SZ effect) -it presents several difficulties. Among the most famous, we recall the "cosmological constant fine tuning problem", and the "cosmic coincidence problem" [9,10].
The first problem is connected to the fact that most quantum field theories predict a huge cosmological constant from the energy of the quantum vacuum at present, more than 100 orders of magnitude too large, see Refs. [9][10][11]. More precisely the theoretical expectations give ρ Λ 10 71 GeV 4 , in contradiction with the cosmological upper bounds giving ρ Λ 10 −47 GeV 4 which gives rise to an extreme fine-tuning problem. It also entails fine tuning at Planck scale era, thus in the initial conditions of dark energy. The second is connected to the reason why dark energy and dark matter energy densities are approximately equal nowadays (see Ref. [12]).
Tensions of unknown origin are also present between the 2013 Planck parameters [13] and σ 8 obtained from cluster number counts and weak lensing, the actual value of the Hubble parameter, The present review will be restricted to the small-scale problems of the ΛCDM model (hereafter SSPΛCDM) connected to the formation of cusps, and to satellites. As we will see hereafter, these issues are strictly connected.
From one side, unified models have shown [59,73,74] that a mechanism which can transform cusps into cores conversely helps the solution of the MSP. Several authors (e.g., [75][76][77]) noticed that the effects of a parent halo's tidal forces on a satellite depend fundamentally on the shape of the latter. A cuspy profile allows the satellite to retain most of its structure, when entering the main halo. Inversely, for a cored profile, the tidal field of the main halo can easily strip the satellite from its gas and even destroy it in some cases [77]. As a result, such satellite will not end up visible, either because it was destroyed or because it lacks the gas to make stars.
From the other side, the satellites most puzzling issue is now recognised as having shifted: rather than the number of satellites, it is related to their inner mass density, specifically to their density profiles being flatter than those of N-body simulations.
In Sections 2, and 3, we will summarise the discussions around, respectively, the satellite plane problem and the Baryonic Tully-Fisher Relation. In Section 4, we will review the CC problem. We will shortly review the early discussions on that problem, then the solution related to the role of baryonic physics, as well as solutions dealing with changes of the dark matter properties or related to modified theories of gravity. In Sections 5, and 6, we review the MSP, and the Too-Big-To-Fail (TBTF) problem. In Section 7, we will discuss a unified baryonic solution to the quoted problems, and Section 7 is devoted to conclusions.

The Baryonic Tully-Fisher Relation
The Tully-Fisher relation was discovered to empirically correlate spiral galaxies luminosities with their HI line-width [60], that is their stellar mass with their rotation velocities [61,120,121]. This correlation is well fitted by a power law with small scatter for late-types of high mass and velocities [61,120], while it no longer follows a single power law at low mass and velocities because cold gas starts to contribute significantly. It has therefore been generalised by using baryonic mass instead of just stellar mass, and the resulting Baryonic Tully-Fisher Relation (BTFR) fits well with a single power law over several orders of magnitude [67] and with small scatter for selected high-quality data [70,122].
The concordance cosmology framework understands the BTFR as an imprint of the finite age of the Universe on the halo mass-rotation velocity relation (see, e.g., [123,124]), that translates into a constant density contrast, and thus implies linear scaling of virial radius with velocity, i.e., from geometry, M ∝ V 3 . The problems for the ΛCDM model coming from the scatter and slope of the BTFR highlighted by [122] (and references therein) hinge on that relation. The link with the BTFR is made by assuming proportionality of galaxies baryonic mass and rotation velocity with their host haloes total virial mass and velocity, respectively. The latter assumption is not trivially satisfied since the observations used for the BTFR clearly provide a very limited baryonic mass and velocity probe compared with equivalent haloes in simulations [125]. In fact, so few halo baryons assemble in actual galaxies that the correlation with virial mass is unclear, while, as the disk sizes remain so small compared to the estimated virial sizes of galaxy haloes, scaling between their two characteristic velocities appears unrealistic. In addition, baryonic feedback mechanisms determining the galaxies' baryonic masses seems too variable to ensure similar galactic baryon fractions for all haloes. The tightness of the observed BTFR is therefore difficult to explain in the ΛCDM model, contrary to the fundamental acceleration scale included in the MOdified Newtonian Dynamic (MOND) proposition [70]. Indeed, [122] claim that observation of a select sample of disk galaxies, assuming constant stellar mass to light ratio, displays a significantly lower scatter than in the ΛCDM simulations of Ref. [126], and that the residuals correlations with the radius or surface brightness of galaxies are not following the ΛCDM semi-analytic predictions of Ref. [127], a conclusion which puts the ΛCDM model at odds with observations. For many years, numerical galaxy formation simulations were unable to produce morphologically realistic galaxies, let alone reproduce the BTFR (see, e.g., [128,129], and references therein), and even semi-analytic models, with empirical inputs, had a hard time predicting it correctly (e.g., [130]). Recently, however, this state of affairs has evolved with improvements in modeling of the baryonic feedback mechanisms allowing production of realistic rotation disk galaxies [131][132][133][134][135]. The small statistics of these models did not allow for predictions on the BTFR, while attempts at reproductions led to controversy on the impact of those baryonic processes on the dark halo, some claiming drastic feedback was needed to obtain the BTFR [136,137] while others reported no need for such feedback (see, e.g., [138,139]).
Recently, the combination of the large EAGLE simulation programme [140,141] -calibrated on small scales on observed galaxy stellar mass function and present radius but not on the BTFR, with multiple realisations of Local Group-like galaxies in smaller volume, i.e. the APOSTLE project [119,142] -claimed to have successfully reproduced the BTFR over four decades [143], even reproducing its break-down at the faint end, as was indeed observed [68,144].
This, after the other claim of successful model from Ref. [145] (a semi-empirical model coupling observed Halo Abundance Matching baryon mass fractions with ΛCDM haloes and claiming to generate a realistic BTFR), seems to have found a path to solve the problem.

The Cusp/Core Problem
Flores & Primack [37] and Moore [31] ruled out cuspy profiles from DDO galaxies' rotation curves, and showed them to be well approximated by cored (or pseudo-) isothermal density profiles. The problem then lies in the cuspy profiles produced in dissipationless simulations of the CDM model (see Figure 1).
Gentile et al. [184][185][186] decomposed the total rotation curves of some spiral galaxies in stellar, gaseous, and dark matter components. Fitting the density with various models they found that constant density core models are preferred over cuspy profiles. Similar results were obtained by Oh et al. [36] using 7 dwarf galaxies from THINGS (The HI Nearby Galaxy Survey) galaxies. The comparison of the RCs (and the density profiles) with the NFW, and pseudo-isothermal (ISO) profiles is plotted in Figure 3. The plot clearly shows that the ISO profile is a better fit to the RCs.
Similar results were also presented using the LITTLE THINGS galaxies [187]. These studies were mainly related to dwarf or low surface brightness galaxies. In the case of high-surface brightness objects, or large galaxies, determining their inner density structure is more complicated. Therefore, stating the nature of the inner density profile for all galaxies, cored or cuspy, is not so obvious nowadays. While, according to [188], high-surface brightness galaxies are cored, other authors (e.g., [156,170,[189][190][191]) conclude differently. The THINGS sample tends to be better described by isothermal profiles (ISO) for low luminosity galaxies, M B > −19, while for M B < −19, cuspy or cored profiles describe them equally well.
The situation is even more confusing as dwarfs do not always have flat slopes, as seen in Simon et al. [189]. They studied the low mass spirals NGC2976, NGC6689, NGC5949, NGC4605, and NGC5963, where they found a large scatter in the inner slope α: for NGC2976, they obtained α 0.01, compatible with a cored profile, while for NGC5963 they got a cuspy one, α 1.28. The other three galaxies had α 0.80 (NGC6689), α 0.88 (NGC5949), and α 0.88 (NGC4605).
In Figure 4, the top left panel plots the DM halo RC of NGC5963 (black dots with error-bars), together with the RC obtained from a fitted NFW profile (cyan line), from a fitted pseudo-isothermal profile (ISO, short-dashed magenta line), and from the model of [166] (yellow dashed line) taking account of baryonic physics. The top right, and bottom panels display the same data for the cases of NGC5949, and NGC2976, respectively, with the same kind of fitted profiles for NGC5949, and just a flat power law (black line) and the same model [166] (dashed line) for NGC2976. The cuspy density profile of NGC5963 is well approximated by an NFW profile, which reflects in the RCs. NGC5949's RC is fitted equally well by RCs from an ISO or an NFW profile, while NGC2976 displays a very flat inner density (α 0.01. All the three RCs are well approximated by the [168] model.
In other terms, if a large part of dwarfs are well described by cored profiles, others are not. Different results have even been obtained using similar techniques for the same object. For example, the dark matter profile inner slope of NGC2976 is bracketed between −0.17 < α < −0.01, according to [192], while [193] got α = −0.90 ± 0.15, and [194] found, considering tracers being stars, α = −0.53 ± 0.14, or gas, α = −0.30 ± 0.18.
Somehow in agreement with the previous discussion, [195] found that the shapes of observed rotation velocity of galaxies display a much larger variation than in simulations.
The discussion above highlights that the determination of the inner slope of galaxies, even for dwarfs, is no easy task. Moreover, it points out that the CC problem must be defined in terms of the inner mass of galaxies rather than of the inner slope of RCs or density profile. The result from the studies discussed above, and of several others, shows the existence of a range of profiles, and that, even with the improvements of nowadays kinematic maps, there is no agreement on the exact dark matter slopes distribution based on morphologies [36,189,194].  [196]'s simulations, that are well fitted by the Einasto profile [39], in the cases of dwarfs (red line), galaxies (green line), and clusters (blue line). The top right corner is the same as top left, but for the Moore profile. The bottom left and right corners display the residuals (Figure reproduced from [196]).
The situation is even more flagrant going to larger masses (e.g., spiral galaxies) dominated by stars, or especially to smaller masses (e.g., dwarf spheroidals (dSphs)) where biases that enter in the system modelling lead to opposite results.
Several techniques have been used to understand and evaluate this problem on dSphs. The spherical Jeans equation gives results highly dependent on the assumptions, since mass and anisotropy of the stellar orbits are degenerate in such model [197]. Maximum likelihood in the parameter space approach applied to Jeans modelling [198][199][200] is plagued by similar such degeneracies. Schwarzschild modelling has been applied to, e.g., Sculptor and Fornax dSphs, finding cored profiles [201][202][203][204]. Methods based on multiple stellar populations concluded that Fornax (slope measured at 1 kpc) and Sculptor (slope measured at 500 pc) have a cored profile [205][206][207][208]. However, a cusp is found in Draco using a Schwarzschild model [204]. This latter results show that in reality there is no accepted conclusion on a unique inner structures of dSphs. A similar problem is also present in galaxy clusters. Sand et al. [209] combining weak lensing, strong lensing, and velocity dispersion studies of the stars of the BCG (Brightest Central Galaxy) found that out of the clusters MACS1206, MS2137-23, RXJ1133, A383, A1201, A963, only RXJ1133 had a profile compatible with the NFW model, and similar studies of Newman et al. [210] (for A611), Newman et al. [211] (for A383), and Newman et al. [212] (for MS2137, A963, A383, A611, A2537, A2667, A2390) (see also [213]) also found flatter profiles than other studies. For example, Donnaruma et al. [214] found a cuspy profile for A611 combining strong lensing and X-ray observations, among other discrepancies from Newman's [212], which covered seven relaxed, massive clusters with flat and cuspy profiles and an average slope α = 0.50 ± 0.1.
While early observations obtained conflicting results concerning the inner structure of the density profiles, high resolution observations, on average, agree on profiles flatter than the NFW's. At the same time, the new observations show a diversity in the inner structure from galaxy to galaxy, as also shown in [195] simulations.
Even if [195]'s simulation results (shaded green band in Figure 5) are in agreement with the RC of galaxy IC2574 at radius > 6 kpc, their behaviour in the inner parts is completely different. This discrepancy points out that convergence in the inner slope of RC between simulations does not mean that they are correctly describing the whole behaviour of the RC. In fact, the deficit of the mass in the inner part of the profile better characterises the CC problem [195]. The key issue is not the shape of the density profile but the excess amount of DM predicted by CDM in the central kpcs of the galaxy. The tension is already evident at scales at which the circular velocity reaches its asymptotic value [224].   [189]. The dotted cyan represents the RC obtained from a fitted NFW profile, the short-dashed magenta, from the pseudo-isothermal profile fit (ISO), and the yellow dashed line, from the model of [168]. Top right panel: same as the top left panel, for NGC5949. Bottom panel: DM halo RC of NGC2976 (black points with error-bars) obtained by [189]. The solid line is the RC computed from a power-law fit to the corresponding density (slope 0.01), and the dashed line, the RC obtained by the model of [168] (from whom this figure is adapted).
Unfortunately, almost all the observation papers in the literature estimate the inner slope of galaxies through α, except for [206] who measured the inner slope through Γ ≡ ∆ log M/∆ log r < 3 − α for Fornax and Sculptor. The latter slope, an integrated quantity, is more easily evaluated than the local quantity α, but only provides some constraints, not a precise value for the inner slope.
On the other side of the mass spectrum, a clear determination of the inner structure of dSphs, cored or cuspy, is very important, as objects with smaller masses are more likely to display a similar inner profile than that of dissipationless N-body simulations predictions (cuspy). Hydro sims: LG-MR + EAGLE-HR, Figure 5. Comparison of the observed IC2574 RC (filled circles) with the two simulated galaxies DG1, and DG2 of [58]. The green line, and shaded region, represent respectively the median rotation velocity, and scatter, of the galaxies simulated by [195] (from which this figure is reproduced).

Early Solutions
Many solutions have been proposed to solve the CC problem and in general the SSPΛCDMs. A decade ago some authors (e.g., Refs. [175,225]) turned against observations, claiming that the inconsistencies could be due to poor resolution or to an improper account of systematic effects. Non-circular motions, beam smearing, off-centring, slit-misplacement, which tend to systematically lower slopes, were charged with the discrepancies.
In HI observations, the finite beam size produces a smear out of HI emission, giving rise to larger disks. The effect depends on the size of the beam, the HI distribution, the inclination angle of the galaxy, and the intrinsic velocity gradients. The problem can be solved with high spatial resolution observations, <1 kpc (see the following).
In Hα observations a slit misplacement may lead to missing the dynamic centre of the galaxy, with the result of having flatter profiles. The problem can be solved in different ways (e.g., 3d spectroscopy [180,188]). The gas is usually assumed to move on circular orbits, so non-circular motions produce an underestimation of the slope. Those motions are, however, of the order of a few km/s, as shown by [183]. Nowadays high resolution observations can distinguish cored and cuspy haloes by deriving their asymptotic inner slopes from rotation curve data [163].
Several authors [178,196,226] suggested ways to reconcile simulations with observations, claiming their simulations were in good agreement with observations, since they become progressively shallower from the virial radius inwards, and that the discrepancy came from projection effects. Indeed, DM haloes are 3-D, but the practice is to compare spherically averaged circular velocity of DM haloes with the rotation speeds of gaseous disks. In other words, the observational disagreement would be with the fitting formulae, rather than with simulated haloes [178].
Nowadays, this proposal is easily rejected: high resolution DM-only simulations have a minimum inner slopes −0.8 [158], while the inner slope of galaxies observed with high resolution techniques are much smaller.
Another possibility took seriously the failure of the CDM model, claiming the problems were with the simulations [161,227,228]. However, modern simulations do not suffer from their past problems: lack of resolution, relaxation, and over-merging. Already then, convergence tests performed by [229] showed that N-body simulations correctly determine the CDM density profiles. The N-body simulations used in the past were dissipationless, meaning they only took account of DM, while baryons are not negligible in the inner regions of galaxies (inner kpc), and dominate over DM in the central 10 kpc of clusters [209][210][211][212]. Nowadays, high resolution cosmological hydrodynamic simulations are available and we will discuss their important role in the next sections, dealing with the baryon solutions of the SSPΛCDM model.

Baryonic Solutions to the CC Problem
As discussed above, the CC problem could be solved with "cosmological solutions" that do not preserve the ΛCDM model, also known as the concordance model. Such modification could, however, alter the successful predictions of the concordance model, that explains many of the observations and features of our Universe.
Thus, before throwing away such a model, it would be wise to verify if some piece of poorly understood or neglected local physics could be connected to the small scale problems.
ΛCDM solutions of the CC problem are based on "astrophysical solutions", for which some mechanism, "heating" the DM, would produce an inner flatter density profile, such as a. the effect of a rotating bar; b. transferring angular momentum (e.g., from baryons to DM) through dynamical friction [166,247,248]; c. AGN feedback, gas bulk motions generated by supernova (SN) explosions [58,75,76]; d. the presence of a central black hole giving rise to a shallower cusp, as claimed by some authors [249][250][251]; e. the role of angular momentum in structure formation.
El-Zant et al. [247,248] showed how clumps of baryons lose energy, transferred through dynamical friction to the DM component of the system, flattening or erasing the natural DM cusp, both in dwarf galaxies and clusters of galaxies. Other authors studied adiabatic contraction of DM haloes [270,271] 5 , that conversely produces a steepening of the density profiles.
Del Popolo (see also [273]) took simultaneously into account these effects [166]: • ordered angular momentum acquired by the proto-structure through tidal torques; • random angular momentum; • energy and angular momentum exchange between baryons and DM through dynamical friction; • and adiabatic contraction.  [238], and turned into an active research field by Starobinsky [239], f (R) theories are defined by a different function of the Ricci scalar in their Lagrangian [246]. Inspired by the Teleparallel Equivalent of GR, the f (T) theories have been introduced to explain Universe acceleration without using dark energy (see Ref. [240]). Finally, the Modified Newtonian Dynamics, was introduced in 1983 by Milgrom [244,245] as a way to model rotation curves of galaxies. 3 General Relativity secondary infall models have been presented by a group around Mimoso and Le Delliou [265][266][267][268]. 4 The peak height of a proto-structure is defined as ν = δ(0)/σ, where δ(0) is the central peak overdensity, and σ is mass variance (see [269]). ν is larger for more massive objects. 5 Calculated through iterative techniques (e.g., [272]). Angular momentum, and dynamical friction in galaxies and clusters not only tend to flatten their density profiles but also to change their global structure [274][275][276][277].
Ref. [166] showed that comparing dissipationless simulations with real structures containing baryons is not correct: the role of baryons in the inner part of the proto-structure is not negligible, explaining the discrepancy between N-body simulations and observed density profiles. In Figure 6, the evolution of haloes of 10 9 , and 10 14 M are shown. Ref. [166] inscribed itself in the "Dynamical Friction from Baryonic Clumps" scenario (hereafter DFBC scenario) discussed in the following.
In Del Popolo [168], the model was applied to dwarfs galaxies showing the influence of the formation history, the content of baryons, and the environment on their density profiles.
A different process was for the first time proposed by Navarro et al. [278], based on supernovae feedback (see the following), that was able to flatten the DM profile.
Before discussing those mechanisms, it is useful to recall how baryons and DM can "interact".
Smoothly distributed baryons with DM give rise to the adiabatic contraction (AC) of DM and baryons collapse towards the galactic centre, which steepens the DM profile and thus increases the central density of the structure. However, the effects of AC can be counteracted if energy is transferred from baryons to DM. This can happen if a. The orbital energy of incoming clumps is transferred to DM through dynamical friction.
If baryonic matter is expelled from, or even moved in, the galaxy (bulk motions produced by supernovae explosions, [75,76]), this produces a temporary flattening of the gravitational potential, moving DM particles outwards, and flattening the cusp.

Supernovae Feedback Flattening
Since the suggestion from Ref. [37], stressed in many following works, of the importance of baryons in solving the CC problem, the first mechanism envisaged was connected to supernovae feedback [58,75,76,119,165,[278][279][280]289].
Ref. [278] showed that the sudden expulsion of baryons into the halo in a single event could flatten the profile. The process is most efficient for galaxies with shallow potential, such as dwarfs. However, Ref. [290] showed that, while a single explosive event did not move sufficient energy to form a core, repeated moderately violent explosions could reach the goal (however see [291] for a different point of view).
Gelato and Sommer-Larsen [279] studied more in detail the [278] scenario, trying to reproduce DDO154's RC starting from NFW profiles: they tried to reproduce a gas outflow event by abruptly changing the disk potential. They found that it was necessary to expel at least 75% of the disk mass in order to reproduce the RC .
Read and Gilmore [280] showed that repeated outflows followed by gas re-accretion, could give rise to a core even in larger galaxies.
Refs. [75,76] showed that random bulk motions of gas driven by SN explosions in primordial galaxies could form a core. Similar results were obtained in [58] simulations. Refs. [36,164] compared the average slope of THINGS dwarfs with the simulations by [58]. Governato et al. [165] ran simulations to study galaxies larger than in [58], and compared the results with observations. They found, for M * > 10 7 M galaxies, a correlation between the stellar mass M * and the inner slope.
Governato's papers used the code GASOLINE [292], a N-Body+SPH code to simulate galaxies. By means of the "zoom" technique [293], the resolution for gas particles was brought down to M p,gas = 3 × 10 3 M , while the DM particles had M p,DM = 1.6 × 10 4 M , and the softening retained at 86 kpc. The authors performed a two runs with different star forming thresholds: one in which stars formed if the hydrogen density was >100/cm 3 (High Threshold run, HT), and another with hydrogen density threshold >0.1/cm 3 (Low Threshold run, LT). These simulations, similarly to [157], implemented the blast wave SN feedback mechanism [294], and/or early stellar feedback [295]. There, the interstellar medium (ISM) received 10 51 ergs of energy from >8M stars. Energy ejected from SN ejected energy was coupled with the coefficient esf to the ISM. The MaGICC simulations [295] adopted the fiducial esf = 0.1.
Similarly, [296] showed that final result of cusp flattening was generated by combining bursty star formation together with supernovae feedback, resulting in fast oscillations of the inner (1 kpc) galaxy potential, and expanding gas bubbles. This process only starts after the galactic centre accumulated cold gas density reaches >100/cm 3 and stars form 6 . Smaller densities (e.g., 0.1/cm 3 ) do not produce any visible changes in the DM inner density profiles. Governato et al. [165], repeated the calculation for larger mass galaxies (see Figure 7, top left and right panels). Teyssier et al. [289] used the adaptive mesh refinement code RAMSES together with a new SN feedback scheme, finding results in agreement with [296] (see Figure 7, bottom left panel), which showed that M * > 10 7 M galaxies have a flat inner profile. Onorbe et al. [298] found similar results, but for M * > 10 6 M galaxies (Figure 7 bottom right panel) 7 .
This view was recently criticised by [195]. According to their simulations, in systems having V max < 60 km/s, baryons have little effect on the rotational curve even in the inner regions of the galaxy.
According to [195], the cores formed in [296] are fundamentally related to an ad hoc choice of parameters, while in their own simulations had no evidence of core formation. 6 Consistent with [58] assumptions; the threshold >10/cm 3 [297] marks the limit for bulk gas flows. 7 In general, the supernovae feedback mechanism is not able to transform cusps into cores in galaxies with M * < 10 6 M .  [289]. The bottom right panel shows the density profile evolution in the [298] hydrodynamic simulations for three different dwarfs: early (all stars form in early times), medium, late.
Those results are in agreement with that of Ref [299]. The latter group simulated 7 high-resolution dwarfs, living in 1 − 2 × 10 10 M mass halos, with different assembly history. They found no case of flattening of the inner core. Their lowest inner slope (at 0.01-0.02 R vir ) was −0.8 and corresponded to a dwarf formed in a halo with a very extended assembly history, which also implies a more extended star formation rate (SFR) history. In addition, [139] got realistic galaxies in simulations, but formed no cores.
Besides the contradicting results on core formation obtained in different high resolution hydrodynamical simulations, the results of [58,300] (and simulations using their same methods and assumptions) have been criticised on several fronts: a. the energetics of the core formation [301] (galaxies with M * < 10 7 M have too few stars to generate the requested energy to flatten the cusp) and the required baryonic mass, marginally exceeding the baryon content of the dSphs [291]. Figure 8 illustrates that problem in its left panel, from Penarrubia [301], while the right panel is reproduced from Maxwell [302]'s study that arrives at opposite results. Moreover, while the solution to the CC problem with the SNFF model needs a large number of SNs, and thus a large star formation efficiency (SFE), the solution of the TBTF problem, places an opposite demand on the star formation efficiency (SFE); The dotted black line corresponds to a core of size r c = 1 kpc. The right vertical axis displays the luminosity, in stellar mass, derived from a star formation efficiency F * = F * (Mvir) constrained by the number of luminous satellites in our Galaxy. Luminosities are converted into SNeII energy using Equation (6) of Ref. [301] and adopting a strong energy coupling DM = 0.4. The outputs of energy found for SN explosion compatible with the "missing satellite" problem in two different studies [303,304]. The tension between the "core/cusp" and "missing satellite" problems becomes obvious in haloes with M vir < 10M (panel reproduced from [301] b. too high a value of energy coupling, SN 0.4, compared to 0.05, a value deduced by [307]; c. a very high star formation threshold [195,308] required to obtain the results of [58]; d. they present difficulties in solving the TBTF problem [291,309,310]; e. they lack resolution to follow the feedback processes which should transform the cusp into a core [135,[311][312][313].
Refs. [36,164] compared the average slope of THINGS dwarves with the simulations by [58], and [165] made a similar comparison for larger objects, and found a correlation among M * and the inner slope for galaxies having M * > 10 6 M 8 .
Conversely, for M * < 10 6 M , hydrodynamic simulations predict cuspy profiles . This result is in conflict with Ref. [206]'s results for Fornax and Sculptor inner structure. They showed that both galaxies are compatible with a core, using the slope of the mass profile Γ ≡ d log M d log r < 3 − α as it gives a more reliable inner slope of DM haloes. Ref. [314] claimed simulations in agreement with [206], injecting however 50% more energy from SNs.
Gnedin and Zhao [290], considering SN feedback in a Semi-Analytic model based on peculiar assumptions, claimed it cannot produce cores. Without discussing their assumptions, even supposing that they were considering all the details of stellar feedback, their objections are lessened as they did not consider the role of the baryonic clumps on core formation as done in Ref. [166].
The only SNFF simulations forming cores in dwarf galaxies with masses < 10 6 M used the GIZMO code in P-SPH mode [298]. This is much more naturally obtained in the DFBC scenario, that we are going to discuss in the next section.

Gas Clumps Merging
As already reported, the other mechanism able to transform cusps into cores is that proposed by Refs. [247,248]. There are three ways, previously discussed, by which DM and baryons can interact. In the case when baryons form clumps of masses 0.01% of the system mass, they can transferthrough their motion -their orbital energy to DM. As a result, and similarly to the SN feedback scenario, DM particles will move towards outer orbits, flattening the inner DM profile. The process is more efficient when it occurs earliest, on smallest haloes.
Ref. [247] showed how such mechanism works in galaxies, and [248] how it works in clusters. Ref. [282] used a similar model to study the evolution of the cluster C0337-2522, finding that after the formation of the brightest cluster galaxy (BCG), the inner DM profile has a baryonically induced inner slope, 0.49 < α < 0.90, smaller than the NFW profile.
Ref. [283] used N-body simulations, with a hybrid N-body/SPH code, to study the evolution of galaxies in the DFBC scenario. They compared their results between the case of DM-only systems and of mixed systems containing DM and baryons. They found that baryons subhalos heated up the cusp, forming a 1 kpc core.
The scenario may be summarised as follows. Initially, the proto-structure contains DM and diffuse gas in the linear phase. DM goes non-linear first and forms the potential well in which baryons will fall.
The rotating disc turns unstable and forms clumps when its surface density, Σ, becomes too large, namely when Q σΩ/(πGΣ) < 1 [324], with Ω the disc angular velocity, and σ its orthogonal 1-D velocity dispersion.
The largest clumps reach radii of 1 kpc (e.g., [325]) and masses of a few percent of the disc mass. Galaxies containing baryon mass of 10 10 -10 11 M typically form clumps of mass 10 8 -10 9 M (see [320,322,323]). Those long lived clumps ( 2 × 10 8 Myr) remain in Jeans equilibrium [322], and are rotationally supported. Ref. [325] showed that the clumps to survive when, in agreement with the Kennicutt-Schmidt law, the gas is converted into stars at a rate of a few percent, similarly to local star-forming systems. The gas clump remains bound, converting into stars, and thus migrates to the galaxy centre. During the collapse phase, baryons are compressed (adiabatic contraction, Refs. [270,271], e.g., at z 5 in the 10 9 M galaxy shown in Figure 6), making the DM profile more cuspy. As dynamical friction (DF) between baryons and DM transfers energy and angular momentum to the DM component, the clumps migrate to the galactic centre. The cusp then heats up, and forms a core.
At later stages (e.g., around z = 2), supernovae explosions provide other events where gas expulsion from the supernova explosion decreases the surrounding stellar density. However, as soon as the smallest clumps form stars with a small part of their mass, they get destroyed by such feedback 9 .
Despite some considering that SNFF and DFBC scenarios represent different implementation of the same idea, based on some common features (e.g., gravitational interaction yielding clumps to DM transfer of energy), they essentially differ for two main reasons: Firstly, the epoch of effective profile 9 The process of star formation is not efficient. flattening are markedly different. While the flattening provided by the DFBC (see Figure 6), starts at higher redshifts (z < 5), at z 3 the DM density profile, in the pure SNFF scenario, the flattening remains similar to the NFW profile [298,300,326].
Secondly, the energy sources moving the clumps come from opposite natures: in the SNFF scenario, the energy of supernovae is driving the clump [75,76], while the DFBC clumps just "passively" infall to the galactic centre, in the sense of Ref. [75,76]'s definition.
Refs. [327,328] compared the two scenarios with high resolution data, from Refs. [189,194], LITTLE THINGS [187], THINGS dwarfs [36,164], Sculptor, Fornax and the Milky Way, examining their respective predictions for the slope-stellar mass, and the slope-circular velocity relationships. They found the DFBC scenario to perform slightly better than the SNFF, in addition to predicting, in DFBC and differently from SNFF, the emergence of cores at smaller stellar masses than 10 6 M . However, even the DFBC mechanism cannot produce cores in very small dwarfs (M * ≤ 10 4 M ) in agreement with Ref. [329] results.
Finally, recall that cusp can reform in galaxies with a bulge [171], We recall that in galaxies having a bulge, the cusp can reform as shown by [171], even in dwarf galaxies [312].

Cosmological Solutions
As previously discussed, the SSPΛCDM model could perfectly well hint at the CDM's paradigm failure. In such a case, the nature of dark matter should be changed , or that of gravity modified, in models, both of which have already been widely checked with various degrees of success.
The simplest possibility grants DM with a small velocity dispersion (σ 100 m/s nowadays, also denoted by a small DM relic mass) [230,330], and is usually referred to as warm dark matter (WDM). As past values velocity dispersion should be higher, such smear could affect small scale structures. This idea leads to an effect resembling the baryonic solutions discussed above: as DM particles retain higher velocity than in the CDM paradigm, small scales cluster less, leading both to flatter profiles and fewer low mass haloes. Many simulations checked WDM structure formation (e.g., [326,[331][332][333]). Although tuning the WDM particle mass to the scale of the halos considered can solve several of the CDM problems, it is not able to get the correct rotation curves for all galaxies or in the entire mass range for which CDM has problems [334]. Both pure CDM and WDM models were explicitly tested against disk galaxies observed rotation curves by Ref. [335], finding no match, however, taking baryons into account, hydrodynamical simulations find the correct rotation curves (see the discussion in Section 4 above).
Moreover, WDM produces too few subhaloes compared to, say strong-lensing subhalo fraction, as shown by the m = 2 keV thermal relic of Ref. [331], and by several other authors (e.g., [336][337][338]). A 1 kpc core requires a 0.1 keV thermal candidate, while large scale structure imply m 1-2 keV, corresponding to cores of 10-20 pc (see Figure 9, and Ref. [339]). WDM too sharp spectrum fall off lead [340] to conclude it does not improve on ΛCDM. To make it worse, WDM alters the Lyman-α forest structure [341] and thus cannot consistently solve the small scale problems and the observed structure of the Lyman-α forest.
The next simple possibility endows DM with self-interaction (Self-interacting DM, hereafter SIDM, [236]), with cross-section of same order as for nucleon-nucleon ( (m/g) −1 cm 2 ) 10 . Redistribution of angular momentum and energy results from elastic scattering in the galaxies inner region, reduces tri-axiality and forms a Burkert profile core [342]. Some cosmological simulations [212,[343][344][345] running with 0.1-0.5 cm 2 /g cross sections, consistent with clusters of galaxy merging observations, [346][347][348], claim that this scattering mechanism solves the CC problem in dwarfs, MW sized galaxies, and clusters of galaxies. The result of [344]'s simulations are presented in Figure 10 for halo mass ranging from galaxies to clusters and two values of σ/m. As the cored SIDM subhaloes feel more tidal stripping and disruption than CDM subhaloes, it solves the CC problem and improves on WDM [344,345] since it leaves enough subhaloes (see however [334] for a different point of view). Note also that SIDM have some natural appeal since particle models of the "hidden sector" produce them, e.g. [236,349,350].
SIDM has been declined in several variations with altered properties: negative scattering leads to repulsive DM (RDM) [232], a condensate of massive bosons, which superfluid behaviour in the central part of DM haloes smoothes down the cuspy profiles [351]. Such non-relativistic Bose condensate was recently simulated for structure formation [352], where large scale structures could not distinguish CDM universes from their model, while flat density profiles and reduced small scale substructures resulted from the opposition, at small scales, between gravity and the uncertainty principle. A further implementation of this model using scalar field condensation, dubbed Scalar Field Dark Matter (SFDM), produced galaxies flat inner profiles [353].
Wave-particle duality was summoned for Fuzzy DM (FDM) [231], ultra-light (m 10 −22 eV) scalar particles which galactic core sized Compton wavelength cannot be "squeezed" further, resulting in flatter profiles and less substructure.
The last two aspects of SIDM are very common when related to indirect DM detection: Self-Annihilating DM (SADM) [234] proposes that the self-interaction results in annihilation, with cross section-velocity σv 10 −29 (m/GeV) cm 2 . In dense regions, the annihilation leads to possibly detectable levels of radiation. Annihilation reduces the structure's particle number in the centre, thus reducing gravity and consequently expanding particles orbits and flattening the density profile.
The second aspect, Decaying DM (DDM) [235], considers DM to decay into relativistic particles, also leading to radiation detection. The gravitational effect is similar to SADM in structure, reducing significantly the core's density of galaxies while large scale structures remain unaltered.

Modified Theories of Gravity
The preceding alterations of DM spoil the simplicity of the CDM paradigm, a further possible solution to the small scale problems legitimately questions the existence of DM and proposes changing Figure 10. Comparison of the NFW (black line) and Burkert (blue line) density profiles with simulations of SIDM universes, using the two cross section over DM particle masses σ/m = 1 and 0.1, denoted SIDM 1 (blue stars) and SIDM 0.1 (green triangles) respectively. The arrow indicated the location of the Burkert's profile core radius (figure reproduced from [344]). gravity itself: this leads to the branch of modified theories of gravity (MG). Although MG is an old issue, the discovery of the universe's accelerated expansion [354,355] literally exploded interest in it. Most of the more recent alternatives to GR are cosmologically motivated, attempting to replace or supplement the concordance cosmology postulates of "inflation", "dark matter" and "dark energy" [356]. They all build on the premises that, although agreeing with GR locally in time and space, gravity may be quite different in the early universe or at large scales.
In this context, instead of interpreting the "anomalous" rotation curves of spiral galaxies nowadays as the trace of missing mass (DM), they reveal a lack in the gravity theory.
From the discovery of the universe's accelerated expansion by the supernova surveys, Einstein's cosmological constant was rapidly reinstated, and quintessence was proposed to overcome the problems that such constant entails. Alternatives to GR also attempted to explain such results.
MOND demonstrated particular effectiveness in solving the SSPΛCDM. An originally ad hoc modification of Newton's gravitation law was proposed in 1983 by Milgrom [244,245], for Newton's second law 11 : There, with the universal constant a 0 10 −10 m/s 2 , the gravitational force results nonlinearly in the acceleration, a, as the functional form µ(a/a 0 ) tends to 1 for high values of the acceleration, while small accelerations modify Newton's law with µ ∼ a/a 0 . 2T gravity

Modified Gravity
Higher-order CDT Figure 11. MGs tree diagram, reproduced from [358] (with permission from Tessa Baker. See updated version here).
In both forms, a >> a 0 recovers Newton's second law, while the case a << a 0 yields a force F proportional to the velocity squared: such that a = GMa 0 /r.
For a test particle in circular motion around the galactic centre, for small acceleration, far away from the centre, the rotation velocity reaches a plateau V f = 4 √ GMa 0 and one re-obtains the baryonic Tully-Fisher relation V 4 f = Ga 0 M b . The success of MOND extends beyond fitting the RCs of spirals, and reproducing the Tully-Fisher relation: it provides explanations for several galactic phenomena, from dwarfs to ellipticals (see [359,360]), to Freeman's law [361], that sets the upper limit for spirals surface brightness, to Fish's law [362], determining ellipticals characteristic surface brightness, and to the Faber-Jackson relation. MOND is actually the simplest way to reproduce observed scaling relations, such as the relation between the rotation curve's shape and the baryonic surface density (see Figure 15 of Ref. [360], relevant to the diversity of rotation curve shapes at a given V max scale), with the stellar and dynamical surface densities in the central regions of disk galaxies (recently discussed in [363,364]), or the small scatter of the BTFR discussed in Ref. [122]. These tight relations appear as less natural consequences of baryonic solutions to the CC problem than of MOND, from a formal point of view. The RCs of two different dwarf galaxies, UGC11583, and ESO138-G014 are plotted in Figure 12, together with their MOND fit. Although MOND usually fits galactic RCs well, some cases, such as ESO138-G014 here, escape its grasp. This should be taken with caution as, according to Ref. [365], there are misunderstandings on MOND on the problems they quote in their section 3. See also Ref. [366] for a different point of view.
At cluster scales, nonetheless, MOND proves much less successful (however, see [360] for a different point of view). In addition, MOND historically being a mere Newtonian fit, to be considered a full theory requires, also to be applicable on cosmological scales, a relativistic version. One such TeVeS theory was proposed by Sanders and Bekenstein [368][369][370].
Since the successes and problems of MOND and the ΛCDM model appear on complementary scales (galactic vs large scales), Khouri [371] proposed combining both theories, keeping each on their respective successful scales. A wide survey of MOND's successes and problems, and of its attempted relativistic extensions was presented in [360].
In summary, there exist several apparently valid proposals solving the CC problem, as well as for the other small scale problems, but the challenge is to single out the correct one.

The Missing Satellite Problem
Galactic mass halos were noticed by Klypin et al. and Moore et al. [32,41] to present many more subhaloes predicted by N-body simulations than observed satellite galaxies. The scale invariant CDM primordial fluctuations at small scales leads, through collapse, to a large number of subhaloes, hence creating this MSP.
The MW counts 9 bright dSphs, Sagittarius, the LMC and the SMC, thus much less than the 500 satellites, obtained in simulations, with larger circular velocities than Draco and Ursa-Minor (i.e., bound masses > 10 8 M and tidally limited sizes > kpc, see Refs. [34,35] and Figure 13, top left panel).
The problem was confirmed in subsequent cosmological simulations (Aquarius, Via Lactea, and GHALO simulations: [158,377,378]). In the end, every cosmological simulations predicts Milky Way-like galaxies surrounded with at least one order of magnitude more small subhalos (dwarf galaxies) than observed (e.g. Via Lactea simulation [378]).
The key idea of the solution lies in the distinction between visible satellites and the entire population: if only a subset of the population is visible, the observed vs predicted satellites discrepancy can be reduced. Thus, various suppression mechanisms for the visible population have been proposed: a. Tidal stripping from the satellites' parent: presently observed satellites had the largest masses before accretion (LBA), large enough to retain visible stars, resisting tides when accreted by their parent [378]. b. Re-ionisation stripping satellites gas, thus star formation, hence suppressing visible satellites formation [385,386]: presently observed, earliest forming satellites (EF see Refs. [32,385] and Figure 13, top right, and bottom left) are visible because they acquired gas, and thus form stars, before re-ionisation.
Ref. [376] compared the cumulative number of dSphs and UFDs (with M/L 1000, from SDSS data) with the cumulative number of satellites obtained in the Via Lactea simulation, assuming a z = 9 − 14 reionisation epoch, almost solving the problem (see Figure 13, bottom right panel). c. Photo-ionisation from stellar and supernova feedback (e.g., Refs. [303,387,388]), and generally stripping gas by ram pressure (e.g., Ref. [389]). Because of the photo-ionisation threshold, UFDs' baryons to stars conversion efficiency lies in the range 0.1%-1%, thus making it is not clear whether they are "fossils" from reionisation epoch [390]. d. Transfer of angular momentum from baryons to DM through dynamical friction [166,273], which turns also cuspy profiles to cores. The number of visible satellites decreases since tidal stripping acts more on cored density profiles.
In conclusion, adding baryon physics to the usual dissipationless DM model solves the MSP, for the majority of the Cosmology community, as well as the other SSPΛCDM, as shown by recent hydrodynamic simulations (e.g., [73,119,391]) or semi-analytic models [59].

The Too Big to Fail Problem
On closer inspection, eliminating visible satellites from the faint end of the distribution does not exhaust the model vs observation discrepancies in satellites (MSP): the most massive (luminous) satellites also pose problem. Boylan-Kolchin [34,35] discovered, in the Aquarius and the Via Lactea simulations, a population of 10 subhaloes that were too massive and dense, by a factor 5, to host even the brightest satellites of the MW, and dubbed it the "Too Big to Fail" (TBTF) problem 12 . In general, ΛCDM simulations of the MW predict at least 10 subhaloes with V max > 25 km/s, while 12 < V max < 25 km/s for all the dSphs. This is shown in Figure 14. Both the left, and right panels display RCs with V max < 24 km/s while simulations have much larger values. Rotation curves obtained from NFW-shaped subhaloes with V max = (12,18,24,40) km/s, with a 1 σ scatter, taken from Aquarius simulations, present the simulation side of the TBTF problem in the left panel. They are confronted, for this purpose, with the bright dSphs, the dots with error-bars. The results show that V max < 18 km/s for most of the dSphs, all have V max under 24 km/s and only Draco is consistent with a subhalo modelled with V max 40 km/s. The right panel plots V max vs the visual magnitude, M V , and confirms that the simulations produce much more massive subhaloes than the observed dSphs.
This issue can be compared with the CC problem for haloes: the ΛCDM seems to produce too much mass in haloes and subhaloes. Thus, similarly to the CC problem, two main classes of solutions have been proposed: cosmological and astrophysical. As for the CC problem, cosmological solutions modify the perturbation spectrum or the nature of dark matter particles. Astrophysical solutions are driven by baryon physics, and consider, in analogy with solutions to the CC problem: a. The shape of satellites inner densities, shifting from cuspy to cored [73,74] (hereafter Z12 and B13 respectively) , thus making them more susceptible to tidal stripping and even to tidal destruction [77,374]. This picture would see the present-day dwarf galaxies, more massive in the past, transformed and reduced strong tidal stripping (e.g., [392]); b. The suppression of star formation from (a) Supernova feedback (SF) (b) Photo-ionisation [74,393], and (c) Reionisation. This can prevent small mass DM haloes gas accretion, "quenching" star formation after z 10 [385,386,394].
This would suppress dwarfs formation or could make them invisible; c. The dynamical effects of a baryonic disc [73,74]. Satellites crossing such disc experience disk shocking, strong tidal effects, even more so for cored inner profiles.
Alternative solutions to the TBTF problem have been proposed: the TBTF excess of massive subhalos in simulations on the MW could vanish if Einasto's are the correct satellite density profiles, or if the correct measurement of MW's virial mass reduces from 10 12 M to 8 × 10 11 M [395,396].
The above solutions to the TBTF problem are focused on the MW satellites. However, the TBTF problem also concerns the Local Group, and Local Volume galaxies [309,310,397], galaxies that are too massive to be modified by (reionisation) feedback [35,72,291,301,[398][399][400], contrary to the baryonic solution of Refs. [59,73,74,401]. Nevertheless, such limitations of feedback were challenged for a few galaxies by [137,302,314], and for an entire galaxy population by [402].

A Unified Baryonic Solution to the ΛCDM Small Scale Problems
So far, each SSPΛCDM was solved separately through different recipes. Unified solutions were however proposed, after understanding their inter-relations. For example, transforming cuspy density profiles into cored distributions, through the SNFF or the DFBC scenario, also affects the parent halo distribution and number of substructure/satellites as tidal effects of a parent halo on a satellite depend fundamentally on its shape (e.g., [75][76][77]): cuspy satellite structure will not suffer big changes when entering the main halo, while the parent's tidal field can easily strip a cored profile from its gas, in some cases down to its destruction [77].
A recent sketch of a baryonic solution to the SSPΛCDM was proposed [73,74], based on SF explosions [58,75,76,[278][279][280] removing angular momentum and gas from the proto-galaxy, thus flattening the density profile. Indeed, the discrepancy between observed and predicted number, and density, of satellites (see [74]) can significantly reduce, by applying a correction to satellites circular velocity at 1 kpc, v 1kpc [73], to large N-body simulations results (e.g., Via Lactea II, VL2 hereafter).
A similar result is obtained in the DFBC [166,247,248]. Indeed, the Del Popolo model [166] provided a similar correction to v 1kpc than Zolotov et al. [73], which, applied to the Via Lactea satellites, as in Brooks et al. [74], simultaneously solves the MSP and the TBTF problems [59,403]. The model of Ref. [403] is in fact more complete: it couples satellites interaction with the halo to dynamical friction energy and angular momentum exchange from baryons clumps to DM, which flattens the density profiles. The effects of tidal stripping and heating on the satellites were modelled with a combination of the procedures from Refs. [77,404].
In summary, the method proceeded in two main phases. The first phase calculated the satellite density profile flattening from baryonic physics (in particular, the subhaloes' central mass lowering), ∆Vc,1kpc km/s Vinfall km/s considering it isolated, without interactions with the host halo, in the same fashion as in Ref. [166] (see also [405,406], to have a semi-analytical description of halos growth).
The flattening was translated into the difference in circular velocity, at 1 kpc, between the DM-only (hereafter DMO) satellites equivalent to those considered and containing also baryons (hereafter DMB satellites), ∆v c,1kpc = v c,DMO − v c,DMB . Figure 15 displays the resulting ∆v c,1kpc at z = 0 and the fitted dashed line through the output points of the model, given by ∆v c,1kpc = 0.3v infall − 0.3km/s 10 km/s < v infall < 50 km/s The second phase no longer considered the satellite as isolated, and subjected it to the host halo's tidal field and accretion.
The merger history and interacting satellites' growth was followed, and the substructure evolution was tracked in Ref. [403]'s model, taking into account the mass loss induced from tidal stripping, tidal heating, as well as the disc's stripping enhancement, caused by the host halo.
The destruction rates by tidal stripping induced a second correction to satellites outputted from N-body simulations, requiring the link between satellite mass loss, or remaining, and velocity change (e.g., V max ) during infall in the parent. Such link was shown in Figure 16: DMB satellites were shown to lose more mass than DMOs because a. DMB satellites contain gas, DMOs do not; b. DMB satellites profiles are flatter than DMO's, thus tidal stripping affects them more (e.g., [77]).
The maximum tidal loss for DMB satellites corresponds to those in the host galaxy disc's vicinity. The analytic fits from Equation (8) of [77], linking the v max change to tidal stripping mass loss, where x ≡ mass(z = 0)/mass(z = in f all), was also shown in Figure 16, for the exponent values ζ = 0.40 and η = 0.24, corresponding to central density profile logarithmic slopes γ = 1.5, represented by a dashed line, ζ = 0.40 and η = 0.30, corresponding to γ = 1, with a solid line, and ζ = 0.40 and η = 0.37, yielding γ = 0, with a dotted line, respectively. Then following [74], a destruction criterion (in terms of mass loss) was fixed, determining tidally disrupted satellites from Equation (6): mass losses • above 97% mass (x = 0.03), or • above 90% mass, with v infall > 30 km/s and pericentric distance <20 kpc, for a given satellite, was set as destroyed. The photo-heating induced star formation suppression, from the [393] results, produced a third correction. The last step assigned surviving satellites luminosity. The stellar mass change with v infall is synthesised by the data fit from Equation (7), shown as a solid line (figure reproduced from [403]). Stellar masses of satellites were allocated relating their infall velocity v infall to the stellar mass M * , as plotted in Figure 17. The fit to the v infall -M * data yielded the relation 13 Finally, the V-band magnitude, M V was related to M * by applying [74]'s relation from [73]'s simulations, Figure 18 showed the result of the all corrections, with the raw results from the VL2 simulations at z = 0 displayed on the top panel, to compare with the same satellites corrected with heating, destruction, and velocity corrections, as discussed above, presented on the bottom panel. There, • Red filled symbols mark "observable" objects produced by the VL2 simulation, • Filled black circles mark much fainter satellites than "observable", stripped of their stars (see [74,77]), their mass loss ≥90% still does not grant them destruction from the above criteria, • Empty circles indicate totally dark objects -Simple empty circles merely lost all their baryons and thus did not form stars as their mass was smaller than the minimum to retain them, while -Empty circles crossed with an "x" represent destroyed subhaloes from baryonic effects (e.g., baryonic disc, etc).
In agreement with [74], 3 satellites with v 1kpc > 20 km/s were obtained. On top of a reduction of the number of satellites to reach the levels observed in the MW, noticeable from Figure 18, clearing up the MSP, the model, applying baryonic correction to subhaloes, Figure 18. Baryonic effects on VL2 simulation subhaloes, seen through their changes in v 1kpc vs. M V . VL2 satellites are labelled, as in [74], in the top panel, by their velocities vs. M V at z = 0. The same, after baryonic corrections, are presented in the bottom panel. Satellites stripped of their stars from losing enough mass are represented by filled black circles: at infall, their actual luminosities are much fainter and indicated luminosities are only upper limits. Satellites actually observable at z = 0 appear as filled red circles. Empty circles mark dark subhaloes, those with an x being most probably tidally destroyed (figure reproduced from [403]). also reduces their central velocity, solving the TBTF problem. As evidenced in Figure 19, those corrections from baryonic physics are sufficient to solve the velocity and number counts problems of MW-type satellites produced by the VL2 simulation.

Conclusions
The discussions on vertical motions of stars near the Galactic plane by Öpik [407] were the first evaluation of our neighbourhood's matter content. One century has passed since, and present evidence points towards a large part of matter in the universe consisting of DM. Despite this mounting evidence, direct or indirect experiments have still not detected DM particles, to date [4,5].
As DM dominated structures, dwarf galaxies, small scale structures, and galaxy satellites are at the front of our understanding of the nature of DM. Moreover, observation of those objects have presented significant discrepancies with predictions of the ΛCDM model. This review discussed some of those discrepancies such as the CC problem, and issues with satellites.
The understanding of the nature of dark matter, and/or the physics of galaxy formation certainly can improve a great deal from the study of those discrepancies. Noticed more than twenty years ago, we still have no definitive idea of what causes the CC problem. Recently, baryonic effect on small scale structure problems have been the focus of several researches. Interaction of baryon clumps with DM, from M > 10 5 M dwarfs, in the DFBC scenario, can create cores, and MW-type galaxies acquire cuspy profiles, as shown in [327]. The SNFF scenario had found similar results, but with cores formed from heavier satellites, M * > 10 6 M [165]. Both scenarios also solve the subhaloes abundance and TBTF problems by reducing the galaxies central DM. However, the baryonic effects debates are still open, as the SNFF scenario presents difficulties in some small scale issues, in particular in some of the MW

Our model
All VL Subhalos MW dwarfs 20 30 10 Figure 19. Corrections from baryonic physics to cumulative number of MW satellites in terms of circular velocity in VL2 simulation. The classical MW satellites augmented by the ultra-faint-dwarfs from [376] are shown as filled squares with error bars. The Via Lactea subhaloes abundance [378] produces the solid line with diamonds. The model of baryonic corrections applied to VL2 subhaloes yields the abundance shown by the dotted line (figure reproduced from [403]). classical dwarfs [291,301], while isolated galaxies do not clearly present baryonic physics solutions to the TBTF problem. In particular, the SNFF scenario, with future flat inner profiles dwarfs with M * < 10 6 M , would conclude that DM is not cold 14 . Further investigations of the SIDM model or its variant would in such case be required to understand if they can solve the small scale issues. In any case, new hints on the small scale issues would require the future surveys to discover new satellites with their stars velocity measurements.
Be that as it may, the difficulty to distinguish cusps from cores retains the debate on dSphs structure open, as discussed in the introduction (e.g., [408]). GAIA [409] and the Subaru Hyper-Supreme-Camera [410] were proposed to provide possible solutions, but both instruments can only hope to solve the problems for larger dwarves, such as Sagittarius [411]. The stars velocity line of sight components and stars' 2D projection radius that they could provide are plagued with degeneracy between the anisotropy parameter and the density profile, because the determination of the latter is based on Jeans equations. Improvements were proposed measuring just one of three velocity components, and two of three position components [see 412], to obtain proper motions of the dwarves stars [374], and thus the half-light radius density slope, despite the challenge it poses with GAIA [411]. The knowledge on dwarves inner structures remains of fundamental importance, no matter the technical problems. Yet, the SNFF scenario's tension with the ΛCDM model is not fundamental: as shown, the DFBC scenario provides a way to form cores from cusps before gas forms stars, and is more efficient. Thus, although the SNFF would have serious problems if cored dwarves with M * < 10 6 M were found, so long as M * > 10 5 M , the DFBC scenario within the ΛCDM model would not disagree with observations. The Baryonic Tully-Fisher Relation and, mostly, the satellite planes problem remain open problems for the ΛCDM model to tackle, although the first one seems to be on the right path.
Alternate tests of the ΛCDM could proceed from checking • The number of small subhaloes in galaxies' virial radius,predicted to be large: gravitational lenses flux anomalies barely found agreement with those predictions [413]; • The Galaxy's cold tidal streams perturbations by subhaloes (see [414]).
Obviously, collider, direct, and indirect DM particles detection would solve the problem of DM nature and existence. Unfortunately, no evidence of super-symmetry (SUSY) has appeared so far at the LHC, so WIMPS particles DM are yet elusive 15 , while a di-photon excess decay at 750 GeV, rumoured at 3.9 σ significance is still debated. The claimed DAMA/LIBRA/CoGeNT annual modulation continues to be the only, although controversial, news from direct searches and similarly no incontrovertible indirect evidence of DM has been presented, to date.

Conflicts of Interest:
The authors declare no conflict of interest.