Unbalanced St\"{u}ckelberg Holographic Superconductors with Backreaction

We numerically investigate some properties of unbalanced St\"{u}ckelberg holographic superconductors, by considering backreaction effects of fields on the background geometry. More precisely, we study the impacts of the chemical potential mismatch and St\"{u}ckelberg mechanism on the condensation and conductivity types (electrical, spin, mixed, thermo-electric, thermo-spin and thermal conductivity). Our results show that the St\"{u}ckelberg's model parameters $C_{\alpha}$ and $\alpha$ not only have significant impacts on the phase transition, but also affect the conductivity pseudo-gap and the strength of conductivity fluctuations. Moreover, the effects of these parameters on a system will be gradually reduced as the imbalance grows. We also find that the influence of $\alpha$ on the amplitude of conductivity fluctuations depends on the magnitude of the both $C_{\alpha}$ and $\delta\mu/\mu$ in the electric and thermal conductivity cases. This results in that increasing $\alpha$ can damp the conductivity fluctuations of an unbalanced system in contrast to balanced ones.


Introduction
The gauge-gravity duality [1] based on the holographic principle, establishes a relationship between gravitational theory in the bulk with d + 1 dimensions and a quantum field theory on the boundary with d dimensions. This duality can deal with lots of unsolved questions in strongly coupled field theories. One of the main achievements of this duality is the establishment of the holographic superconductors [2][3][4][5].
More percisely, the standard BCS theory [6,7], which can describe the properties of low temperature superconductors, is not capable of explaining the properties of some unconventional superconductors whereas the gauge-gravity duality may helps us to handle strongly coupled systems and understand some features of high temperature superconductors. This duality relies on the mechanism of spontaneously breaking of the U (1) symmetry in the dual field theory. This holographic model undergoes a phase transition from the black hole with no hair to the scalar hair at low temperatures [8,9]. There exists several studies on the holographic superconductors to describe their different aspects [10][11][12][13][14][15][16][17][18][19]. One of the interesting features of the Abelian Higgs model for high themperture superconductor is the second order phase transition as it shown by Landau-Ginzburg theory and the ratio of pseudo-gap frequency (ω g ) to critical (T c ) is similar to the high-T c superconductors (ω g /T c ≈ 8) [20].
It is also interesting to take an effective field theory approach and consider the existence of a U (1) symmetry breaking via the Stückelberg mechanism [21][22][23]. This model depends on a general function F of the scalar field. One of the main features of this phenomenological model is provision of a large group of phase transitions which are first order. There are -1 -some other works having regard to this model [24,25]. Moreover, in conductivity case it might be existed additional resonances at non zero frequencies for some choices of F. One can interpret these poles as sing of the existence of quasiparticles in the superconductor. The similar behavior can be observed once the scalar field mass approaches the BF bound [26].
By ignoring the details of the quantum field theory side and only consider broken symmetries, one can also apply the Stückelberg mechanism to an unbalanced holographic superconductor [27,28]. This model is based on emerging superconducting phase around a quantum critical point [29]. The mechanism of this model is that the superconductive phase happens where the two fermionic species contribute with unbalanced populations or unbalanced chemical potentials. This is a relevant subject both in condensed matter systems and finite density QCD [30]. The unbalanced chemical potential can be produced by magnetic impurities in a system or by an existence of external magnetic field inducing Zeeman splitting of single-electron energy levels. Moreover, in Ref. [31,32], Larkin, Ovchinnikov, Fulde, and Ferrel showed that, except for the normal/superconductor phase transition, system also experience a new state called LOFF phase. This inhomogeneous phases with spatially modulated condensates lead to spontaneously non-trivial spatial modulations in superconducting condensate. Adding a non-trivial charged field on the gravity side leads to the breaking of a U (1) A "charge" symmetry which is the characterizing of superconductivity [2,3,9]. The chemical potential mismatch is also a potential for a U (1) B "spin" symmetry when the scalar field is uncharged [33]. These two gauge fields correspond to two currents in the boundary theory which provides us with the strong-coupling generalization of the two-current model proposed by Mott [34]. Furthermore, The mixing effects of the two currents lead to obtaining the spintronic features. One can, therefore, investigate the mixed spin-electric linear response properties by using the holographic method.
In this paper, we study the unbalanced Stückelberg holographic superconductors where the effect of backreaction of matter Lagrangian on the geometry has been considered. In other words, we want to see that, in the presence of imbalance, how valid the behavior of holographic Stückelberg superconductor obtained in [21,22] is. Or equally, we look for how behaviors of unbalanced systems obtained in [27,28] are changed by considering the Stückelberg mechanism. The effects of this mechanism are characterized by form of function F(ψ) which one can back to the Higgs model by setting F(ψ) = ψ 2 . Therefore, in order to trace these effects, we need to construct the conductivity matrix describing the linear response of the system to a small electric field and a small temperature gradient. In most cases, the results show that the imbalance makes the effect of Stückelberg mechanism weaker. However, in some situations, diagrams illustrate complicated behaviors.
The paper is organized as follows. In Section 2, we introduce the Lagrangian for our model. We also represent results of numerically calculated condensation and phase transition for different cases of δµ/µ (the ratio of chemical potential mismatch to chemical potential where indicates the amount of imbalance) and the function F(ψ). In Section 3, we briefly introduce the process of calculation for all types of conductivities. Then, we study the behavior of conductivities of different unbalanced systems by considering various forms of the function F(ψ). Finally, conclusions are presented in Section 4.

The Model
We consider an extension of the generalized Stückelberg model introduced in [21] which an extra U (1) gauge field B has been added to the action to produce the source of spin current on the boundary. Note that the scalar field ψ is uncharged under the additional gauge field B. Therefore, the bulk action for such an unbalanced Stückelberg model in (3+1)-dimensions is defined as: in which F = dA and Y = dB are the two field strengths associated with the two gauge fields. The Maxwell equation makes the phase of ψ constant, so we take it to be null in order to have real ψ. In addition, this theory is invariant under the local gauge symmetry A → A + ∂Ω(x) and p → p + Ω(x) [21]. Therefore, we can utilize the gauge freedom to fix p = 0. We also set L = 1 and 2κ 2 4 = 1. Moreover, the function F(ψ) can be written in the general form as: It should be noted that the choice of F causes to change the properties of the CFT at the boundary [22]. Moreover, in the effective field theory context, changing F corresponds to a sort of "non normalizable deformation" or equivalently a change in the theory. Because of the positivity of the kinetic term for p, we must take the above function to be positive. It is obvious that for F(ψ) = ψ 2 , our model will reduce to the unbalanced model in Ref. [27,28]. The plane-symmetric black hole with the backreaction effects is described by the metric: We also consider the following ansatz for scalar and vector fields: Furthermore, the temperature of such black holes with the horizon at r = r h is defined as: By varying the action with respect to the metric and scalar and vector fields, we arrive at the following equations of motions.
where prime denotes derivative with respect to r and dot denotes derivative with respect to ψ. For solving the above equations we need to apply numerical method when one imposes two suitable boundary conditions at the horizon r = r H and the conformal boundary, r = ∞. We also take the standard choice of mass as m 2 = −2 [35,36] and restrict the potential to V (ψ) = m 2 ψ 2 containing just the mass term. For our case, in which m 2 = −2 > −4, the Breitenlohner-Freedman (BF) bound [37] is respected. The asymptotic behaviors of the scalar and gauge fields near the boundary are: where ψ 1 (ψ 2 ) can be regarded as the source of the dual condensation operator, O 1 (O 2 ). Since we need the U (1) symmetry to be broken spontaneously, we turn one of the sources off, i.e. ψ 1 = 0 and then set O 2 = √ 2 ψ 2 . According to the gauge/gravity duality, the leading terms of φ(r) (v(r)) are interpreted as chemical potential (chemical potential mismatch) and charge density (charge density mismatch) in the dual theory, respectively. Working in the grand-canonical ensemble, we fix chemical potential (and chemical potential mismatch) and live the charge density (and charge density mismatch) valuable.
At the AdS boundary, we also should set χ → 0 and impose the asymptotic behavior g(r) = r 2 − 2r + ..., (2.15) where is the mass of black hole interpreted as the energy density of the dual field theory [3]. The other boundary conditions are those which should be imposed at the horizon, r = r h . In this region, both g(r) and the temporal components of the gauge fields should vanish; therefore we have By substituting Taylor expansion of fields at horizon in (2.7) and making use of the Einstein equation (2.10), the black hole temperature could be rewritten as where the subindexes h0 and h1 indicate the coefficients of field's expansion about r = r h . Both bulk and the boundary theories have the same time coordinate and, consequently, they have the same complex time continuation and temperature. We numerically solve the equations of motion ((2.8)-(2.12)) by integrating from the horizon out to the infinity with respect to introduced standard boundary conditions. We mostly consider the interval 0 ≤ δµ/µ ≤ 2 where µ is always fixed at µ = 1   Table 1: Value of critical temperature T c as a function of α for given F(ψ) = ψ α and fixed δµ = 1.

Condensation and phase transition
In this part, we are looking for the properties of the phase transition by studying the condensation of the scalar operator for some considered forms of function (2.4). At first, -5 -by checking the second order phase transition diagrams in (T c , δµ) plane for some F = ψ 2 + C α ψ α , with α > 2, and fixed µ = 1. Form Fig. (1), we find out that the critical temperature is not affected by the parameters in (2.4) since we have α > 2. Of course, it could be predictable since at limit ψ → 0 (near the normal phase), the ordinary form ψ 2 for function F(ψ) is dominated. While, if we choose, for instant, F(ψ) = ψ α , like [21] the critical temperature will be affected by α. We check numerically this assertion by plotting T c with respect to α for the function F(ψ) = ψ α and various δµ in Fig. (2). These curves explicitly show the T c dependence on α as well as δµ. However, for our model in which function (2.4) with α > 2, the T c is only affected by δµ/µ. We also represent some data in table (1) which indicate the dependence of T c on α for chosen function F(ψ) = ψ α when we fix the chemical potential mismatch in the way δµ/µ = 1. In next subsections, we choose some specific forms of the function F and investigate the phase transition type for those. We start with the special case of F(ψ) = ψ 2 + C 4 ψ 4 to see how the phase transition behaves in the interval 0 ≤ δµ/µ ≤ 4. Figs. (3) illustrates the change of phase transition order caused by raising C 4 . Moreover, Fig. (4) demonstrates that this influence of reducing C 4 over phase transition is rather more remarkable in the less unbalanced systems. The results are detailed as follows: • Figs (3) and (4) show that the change of phase transition order caused by increasing -7 -    (d) illustrate that our most unbalanced systems, i.e. δµ/µ = 2 and 4, do not experience the first phase transition even at C 4 = 3. As a result, increasing the imbalance in a system makes it harder to switch the order of phase transition from second to first by C 4 .
• We numerically check that the condensations approach zero as with mean field critical exponent β = 1/2 for the second order phase transitions. Thus, β is independent on chosen δµ/µ or C 4 in this case.
As clearly shown in Fig. (5), in this case the first order phase transition occurs for nonvanishing positive C 3 , and δµ/µ has no effect on the order of phase transition. Since all the phase transitions are first order, there is not a relation like (2.18) in this case.
We are interested to investigate the effect of α, for C α < 0, on critical exponent β and search for non mean field behavior. We check that the relation obtained in [21] for interval remains unchanged even in unbalanced systems. As indicated in Fig. (6), the relation has been checked for some different values of δµ/µ and 3 ≤ α < 4.
All the data of various unbalanced systems are almost located on each others in Fig. (6) (c) and (d). This clearly implies that imbalance has nothing to do with the gradient of condensation diagram near the critical temperature. Relation (2.19) means that the β becomes larger than the mean field critical exponent for 3 ≤ α < 4. Such non mean field behavior follows suppression of fluctuations and stability of condensation (as observed in the Gross-Neveu model for massless fermions [38]). Likely, it manifests the existence of long-range interaction and chiral symmetry in the boundary theory [39,40] (a) (b) (c) (d) Figure 6: Condensation near the critical temperature for the function F(ψ) = ψ 2 −ψ α +ψ 4 . Each plot of the first array indicates the condensation for fixed δµ = 0.5, 1.

Conductivity
In this section, we study the properties of system conductivities. In addition to consider mixed spin-electric linear response to external gauge fields fluctuation, let us add the thermal effects, namely the thermo-electric and thermo-spin linear response of the system. Therefore, we can write explicitly the conductivity matrix as follows: which encodes the whole system response. The diagonal components σ A , σ B , and κT stand for "electric", "spin", and "thermal" conductivities respectively. The off-diagonal components obviously indicate mixed effects; i.e. γ, αT , and βT accounts for the "mixed", "thermoelectric", and "thermo-spin" response respectively. The time-reversal invariant equilibrium states of system implies the symmetry of the matrix.
To study the transport behavior of our thermodynamical system we have to consider small variations of the sources and the consequent current flows. Therefore, considering the fluctuations of fields A and B in direction x with time dependence of form e −iωt , we switch on external field perturbations in the bulk in order to obtain the conductivities in the dual field theory as a function of frequency, ω. Afterwards, by substituting derived Einstein equation in the two Maxwell equations on the background, and eliminating metric fluctuations, one arrives at the two following linear differential equations mixing fields A x and B x : Note that the backreaction effect leads to couple the different gauges. Therefore this event causes to appear the mixed spin-electric transport properties of the system [27]. We can also consider near-horizon behavior ansatz which also impose ingoing boundary conditions at horizon. In addition, the asymptotic behavior of fields around boundary r → ∞ are Using introduced method in [27], we can finally get and thermo-electric and thermo-spin conductivities are obtained as follows: Finally, we find that the non-canonical thermal conductivity is given by where we have considered pressure p = /2, as in [27], in order to account for contact terms not directly implemented by the previous computations (see Herzog's review in [4]). To find more details about equations (3.9), (3.10), and (3.11), see [27]. By numerical solution of equations (3.2) and (3.3) and utilizing (3.9), (3.10), and (3.11) we able to study the effects of model parameters and imbalance on all the conductivity types.

Diagrams and behaviors
We restrict ourselves to the superconducting phase of system at specific temperature T = 0.3T c . Like before, chemical potential is always fixed at µ = 1. Note that the imaginary part of conductivities has a pole at ω = 0, which, according to the Kramers-Kronning relation, translates in a delta function at the same point in the real parts, because of the system translation invariance. Since we are working with the fully backreacted solution, translational invariance is preserved due to the lack of dissipation appearing in probe approximation. Furthermore, since the Ferrell-Glover-Tinkham sum rule implies that the area under the curves must be constant at different temperatures, we have a depletion at small frequencies to compensate the development of the delta function at ω = 0 [3]. According to the terminology used in [27], we use the term "pseudo-gap" to describe the depletion at small frequency since the real part of the electric conductivity appears exponentially small with respect to T but is not zero even at T = 0.

Conductivity Behavior with respect to the variation of δµ/µ
We first study conductivity behaviors of a Stückelberg superconductor in the presence of imbalance. We are actually looking for recovering results of [27] for our model. We can therefore fix the function F(ψ) and temperature to study the conductivities related to systems with different imbalances, i.e. δµ/µ = 0, 0.5, 1, 1.5. Fig. (7) illustrates the decline in the pseudo-gap of electric conductivity as system becomes more and more unbalanced, which is also reported in [27]. For example, Fig.  (7) (a) demonstrates that the pseudo-gap of system with δµ/µ = 1.5 almost vanishes. Moreover, one can easily see that the difference between pseudo-gap values for each δµ becomes negligible at high enough value of C 4 .
In order to clarify this point, we also depict ω g /T c as a function of C 4 for different values of δµ/µ in Fig (8). One can demonstrate that the differences between ω g /T c of various unbalanced systems almost vanishe for large C 4 .
-11 -   (7) also indicates that the increasing coefficient parameter C 4 may make the coherent peak turn to a delta function ( See Fig. (7) and (11) (a)). We return to this point at the next subsection. Note that this happen leads to increase the pseudo-gap value by jumping ω g to the next peak frequency. Therefore, the differences between ω g /T c of various unbalanced are still the same for high enough values of C 4 till we ignore this jumping.  In addition, imbalance disturbs the constant values of spin and mixed conductivity. Note that the optical spin and mixed conductivities of a balanced system are the constant values of 1 and 0, respectively. In our model, the optical spin and mixed conductivities of unbalanced systems will be saturated to these values at large ω after some fluctuations.
Opposite to electric conductivity, as it is clear from Fig. (9), the optical spin conductivity becomes more and more depleted at small frequencies as imbalance grows [27]. The appearance of pseudo-gap is also observable in highly unbalanced systems. This opposite behavior of electric and spin conductivities with increasing δµ/µ is usually interpreted as a separation of the dynamics of charge and spin degrees of freedom [27]. Hence, from the pseudo-gap of spin conductivity view point, effect of increasing imbalance is in the same direction as the effect of increasing C 4 .
The real part of mixed conductivities are reported in Fig. (10) (a). It seems that -13 - increase of imbalance causes some fluctuations which are not only intensified but also shifted to larger frequencies.
Diagrams of real part of thermo-electric conductivity for fixed function F(ψ) = ψ 2 +3ψ 4 are represented in Fig. (10) (c). They show some fluctuations at small frequencies before relaxing to −1 at large frequencies. More unbalanced systems tend to generate a positive peak in conductivity before providing some negative fluctuations. This behavior therefore kills the pseudo-gap in more unbalanced systems. Moreover, increase of parameter C 4 not only does not disturb the general behavior with imbalance, but also highlights it (we will return to this point in the next subsection.).
Imbalance also turns thermo-spin conductivity on and makes it stronger such that for more unbalanced system we have larger negative conductivity as shown in Fig. (10) (d).
At the end, the behavior of thermal conductivity may seem like electric conductivity in the sense that imbalance reduces pseudo-gap (Fig. (10) (d)). Whereas, except balanced systems, there is not so different between pseudo-gap of systems with different imbalances. Therefore, similar to the unbalanced holographic superconductors [27], it seems that there is not mono-tonic behavior by changing imbalanced at small frequencies in our model. It -14 - is obvious that the real part of thermal conductivities of systems which are different in imbalance do not rest to the same value at large frequencies, like thermo-spin ones.

Conductivity Behavior with respect to the variation of C 4
In this section, we focus on the special case of F(ψ) = ψ 2 + C 4 ψ 4 in order to identify influences of C 4 . By increasing C 4 the pseudo-gap and the coherent peak become narrower and stronger at low frequencies, respectively. Comparing plots of Fig. (7), one can also observes that parameter C 4 gradually loses its control over coherent peak as system becomes more unbalanced, while this parameter keeps making the pseudo-gap wider in even strictly unbalanced systems. As mentioned, increasing coefficient parameter C 4 (and also α and C α ) may develops an extra delta function in pseudo-gap region. Formation of these kinds of delta function, which is the direct consequence of mapping poles of the imaginary part by a Kramers-Kroning relation just as the case with the delta at ω = 0, actually signs to extra resonances. We should expect such resonances since we have the "vertex" ψ α (∂p − A) 2 with α ≥ 3 providing inelastic scattering [41]. For the spin conductivity, C 4 does not control conductivity fluctuations as much as it does in the case of electric conductivity. Fig. (12) indicates optical spin conductivities for two unbalanced systems with δµ/µ = 1, 1.5. Moreover, increasing both δµ/µ and C 4 results in making the depletion at small frequencies bigger .
For mixed conductivities, increase of C 4 shifts fluctuations to larger ω s. This happens because of suppression of negative fluctuations at small ω and amplification of positive ones at larger ω, see Fig. (13). This makes fluctuations to shift to positive conductivities as well, which is more noticeable in the less unbalanced systems.
In the case of thermo-electric conductivity, Fig. (14), C 4 shows a rather remarkable control over fluctuations. In more unbalanced systems, increasing C 4 intensifies fluctuations not only in the negative direction of conductivity but also in the positive direction (at smaller ω).
Although the creation of the positive fluctuation in more unbalanced system kills the pseudo-gap, it survives in less unbalanced ones (e.g. systems with δµ/µ = 0.5) and also becomes wider by raising C 4 . Fig. (15) shows the appearance of slight fluctuations, caused by increasing C 4 , at middle frequencies in thermo-spin conductivity. Fluctuations are also suppressed as long as the system becomes more unbalanced.
In Fig. (16) for the real part of thermal conductivity, the fluctuations are also dominated by increasing C 4 values. Furthermore, coherent peak gets stronger and is shifted to larger frequencies once C 4 grows. But, as evident in Fig. (16), it seems that these behaviors gently vanish in highly unbalanced systems.

Conductivity behavior with respect to the variation of α
By fixing parameter C α and temperature, we can investigate how conductivities behave with variation in parameter α for unbalanced systems. From (17) part (c) to (h), one can find that pseudo-gap of an unbalanced system becomes smaller by growth of α. Therefore, for unbalanced systems, a larger α makes the coupling weaker.
Moreover, part (a) of Fig. (17) obviously shows that, in balanced systems (δµ → 0), the exponent both α and C 4 control the strength of fluctuations [22]. But, although growth of the α parameter makes coherent peaks of the optical electric conductivity sharper and stronger in balanced systems, it is not what always happen in the case of unbalanced superconductors. For instance, in Figs. (17) (c) and (e)-(h), one can see suppression of fluctuations by increasing in α amounts. It is interesting that for δµ/µ = 0.5 one can find both the prevalent and the opposite behavior for C α = 2 and C α = 6, respectively. Therefore, intensity of peaks depends on the magnitude of C α . Furthermore, for δµ/µ = 1 and 1.5, Figs. (17) (f) and (h) present the opposite behavior even for much larger C α . Indeed, we have failed to reach the large enough C α making the fluctuations intensified by increase of C α (it does not happen even at C α = 14). Note that, in this case, we can also change amounts of the parameters C α and α to make ω g /T c approach 8. For example, in δµ/mu = 1.5 and C α = 8, we can reduce the ω g /T c to about 8 by setting α = 6 (Fig (17) (h)). Fig. (18) also implies that the effect of α on optical spin conductivity leads to exist a slight reduction of the pseudo-gap and reinforcement of the fluctuations. It means that the spin conductivity is not as much sensitive as other conductivity types to parameter α.
However, the mixed conductivity shows intensive fluctuations when the α increases (Fig. (19)). Like the other model parameter C α , the influence of α on unbalanced systems is weaker compared with the balanced system.
As shown in Fig. (20), the general thermo-electric conductivity behavior with α indicates a shift in fluctuations from high frequencies to smaller frequencies. Similar to the electric conductivity, whether increase of α intensifies fluctuations or not depends on the values of both C α and imbalance. From Fig. (20) (b), it is clear that for a less unbalanced systems like δµ = 0.5, parameter C α = 6 is large enough to have fluctuations intensified by α. Nevertheless, for more unbalanced systems δµ = 1 and 1.5, Fig. (20) (a) and (c)-(f), even C α = 8 is not large enough to make fluctuations stronger and, oppositely, they get suppressed by increasing α.
Moreover, about optical thermo-spin conductivity type, one can figure out from Fig.  (21) that the growth of α parameter produces a slight fluctuations in middle frequencies. Similar to the spin conductivity, these fluctuations do not obey an explicit pattern, but they are damped by increasing imbalance. As it is obvious from Fig. (21), we need larger C α to well demonstrate the fluctuations of more unbalanced systems.
The thermal and electric conductivities behave with varying the α parameter in almost the same manner. For the balanced case, the real part of thermal conductivity reduces to the optical electric conductivity. Generally, the variation of the α parameter leads to shift conductivity fluctuations and coherent peak to smaller frequencies, while amplification of -21 -them depends on the value of both imbalance and the coefficient parameter C α . According to Fig. (22) (a) and (b), there exist two opposite behaviors with α for two different coefficient parameter C α associated with the system with relatively low imbalance δµ/µ = 0.5. Nevertheless, like the case of electric conductivity, more unbalanced systems does not reach the intensified fluctuations in our rang of parameter C α (even for C α = 10 in system of δµ/µ = 1). Figure 17: The optical electric conductivities in terms of ω/T c for function F(ψ) = ψ 2 + C α ψ α with different α = 3, 4, 5, 6 (solid curve, dotted curve, dashed curve, and dot-dashed curve). Each row is related to the systems with same imbalance (as before chemical potential fixed at µ = 1), we have δµ = 0, 0.5, 1, 1.5 from up to down. -23 -  Figure 19: The real part of mixed conductivities in terms of ω/T c for function F(ψ) = ψ 2 + C α ψ α with different α = 3, 4, 5, 6 (solid curve, dotted curve, dashed curve, and dotdashed curve) and fixed C α = 6. Each figure is associated with systems with fixed µ = 1 and δµ = 0.5, 1, 1.5 (figures (a), (b), and (c)).  Figure 20: The real part of thermo-electric conductivities in terms of ω/T c for function F(ψ) = ψ 2 + C α ψ α with different α = 3, 4, 5, 6 (solid curve, dotted curve, dashed curve, and dot-dashed curve). Each row is related to the systems with same imbalance (as before chemical potential fixed at µ = 1), we have δµ = 0, 0.5, 1, 1.5 from up to down.

Conclusion
We studied the unbalanced holographic model in combination with Stückelberg mechanism which gives us a highly flexible dual theory. This flexibility provide us to make a suitable theoretical theory matching the experiment. We showed that while model parameter C 4 provides the change of phase transition order from second to first, imbalance makes it harder. In other words, we need larger C 4 for more unbalanced system to change the order of phase transition. Such behavior also can be observed from conductivity diagrams. In most cases, conductivity behavior of more unbalanced systems are not under the influence of the model parameters of Stückelberg mechanism as much as of the less unbalanced one. In other words, Stückelberg mechanism generally lose its effects as system become more unbalanced. Moreover, we have numerically recovered the Eq. (2.19) also for the unbalanced system. Additionally, we have found that imbalance can significantly deviate the system's behavior with model parameters of Stückelberg mechanism. Such deviation may be so strong which invert the behavior in some situations. We can specifically mention the behavior of electric and thermal conductivity with model parameter α. For example, electric conductivity fluctuations of the relatively less unbalanced system with δµ/µ = 0.5 are intensified as C 4 = 6 although they are damped as C 4 = 2. The same has been observed from thermal conductivity too.
At the end, it is worth mentioning that, as future task, we should push more towards the experimental directions and comparisons. By making use of the method introduced in [42], it would be interesting to take advantage of the freedom of F to simultaneously match two phenomenological behavior, i.e. the phase transition and conductivity.