Bifurcations in the RG-flow of QCD

Bifurcation analysis is used to study an effective model of QCD$_4$ with four-fermi interactions. Our analysis supports the scenario of a fixed point merger at the lower edge of the conformal window. This indicates square root scaling of the anomalous scaling dimensions of the fermion fields just above the lower edge and exponential scaling just below. We also predict existence of new fixed points in this model whose (dis)appearance may indicate transitions of the flow within the conformal window. Furthermore, we make new predictions for the critical value $(N_{f}/N_{c})_{\textrm{crit}}$ at the lower edge. We also obtain exotic spiraling flows that are generated by complex scaling dimensions of the effective four-fermi interactions. Finally, we extend the model by adding a scalar field that couples with a Yukawa interaction term and study the modifications it causes to RG-flows.


Introduction
Renormalization group describes the flow of a theory in the space spanned by the coupling constants called the "theory space". In quantum field theory renormalization group flows are governed by a set of beta functions that result in ordinary differential equations. The flow induced by these beta functions provide information about the change of coupling constants and anomalous field dimensions with the energy scale. An interesting aspect of the RG flows are the fixed points where the beta functions vanish. In particular one is interested in determining the positions of these fixed points and their stability properties. In addition one would like to obtain possible renormalized trajectories generated by the beta functions. Since the RG flows are generated by ordinary differential equations, it is natural to study them with the tools developed for the theory of dynamical systems. In this article, we use analytical and numerical 1 methods for bifurcation analysis to study bifurcations in renormalization group flows, as recently proposed by Gukov [3].
Bifurcation theory, which we summarize in Appendix A, categorizes topological changes in a dynamical system and is therefore an ideal tool to study the appearance and disappearance of fixed points. Using this categorization, we can not only determine easily when fixed points appear and disappear on the RG flow, but also obtain information about the nature of these transitions where the fixed points (dis)appear. As discussed in [4], there are three typical ways in which fixed points of an RG-flow (dis)appear: • Merger of a fixed point with non-trivial coupling constant and a fixed point with a trivial coupling constant. At this merger the stability of the non-trivial fixed point changes.
• Merger of two fixed point with non-trivial coupling constants. At this merger both fixed points disappear.
• A divergent fixed-point as one or multiple coupling constants run off to infinity.
The first and the second scenarios are examples of standard bifurcations. The second is naturally described by the well known saddle-node bifurcation, while the first is called a transcritical bifurcation. It should be noted that in this scenario the fixed points usually do not merge as was stated in [4] but collide and exchange stability. However, at this transition the coupling constant at the non-trivial fixed point generally becomes negative, hence moves into a physical uninteresting region. Furthermore, it is an easy exercise to show that the first scenario always induces a linear scaling of the anomalous dimension at the fixed point [3], while the second scenario implies square root [3] and exponential (also known as BKT or Miransky scaling) scaling of the anomalous dimensions [4] at the two sides of the merger. Finally, the third scenario above with diverging fixed points is much harder to describe using the techniques from the bifurcation theory. However, this scenario is mostly encountered in a perturbative analysis of the RG flows which becomes invalid anyway at the divergence. It turns out that in some cases this last scenario is replaced by the second one in an effective field theory.
In this paper, we apply bifurcation analysis to an effective model of QCD as suggested in [3]. We consider QCD 4 with gauge group SU (N c ) with an arbitrary number of colors N c and flavors, N f , as the bifurcation parameters. The Lagrangian is given by where G A are the gluon field strengths and ψ i are the fermions of flavor i. Furthermore a, b, c are used to label the color in the fundamental representation, and A, B, C is used for the color in the adjoint representation. We have We will consider the massless case, m i = 0 ∀i, where the Lagrangian becomes invariant under a global SU (N f ) L × SU (N f ) R × U (1) symmetry which can be broken to SU (N f ) V by formation of a chiral condensate.
The two-loop beta function for the gauge coupling is given by [5] β [2] g ≡ Λ with α g := g 2 (4π) 2 and are the scheme independent beta-function coefficients. This beta function has a (degenerate) trivial root at α g = 0, and a non-trivial root at α g = − b0 b1 . The non-trivial fixed point becomes negative and therefore unphysical for and where the latter bound runs between 3.4 and 2.6 for 1 < N c < ∞, thus it is always below the first bound for any N c . The region between these bounds, where a non-trivial fixed point exists, is known as the conformal window (cf. figure 1), and its existence was first discussed in [5,6].
Within the conformal window the non-trivial fixed point is attractive in the IR limit, while the trivial fixed point is attractive in the UV. Above the conformal window, the trivial fixed point becomes an IR attractor. Since α nt g becomes small close to x = 11/2 higher order contributions to the beta function (2) can be neglected hence the value of upper bound is robust. However, close to the lower edge α g becomes very large, invalidates the above estimate and (6) becomes unreliable. In particular, other fixed points resulting from higher order terms might influence this critical value.
The type of phase transition that occurs at the lower edge of the conformal window, that we call x crit , and its value has been a central problem in the literature. Use of Schwinger-Dyson equations [17,19,21,22] generically suggest that the theory is in the confining phase with chiral symmetry breaking for x below x crit , see figure 1. As for the value of x crit , estimates for N c = 3 are typically in the range 6 x crit 12 [7][8][9][10][11][12][13][14][15][16][17][18][19][20][21][22][23], and follow from the various techniques (cf. Table 6 of [3] for a recent overview of various studies and their prediction). A more recent study [24,25] uses Potts models to discuss the behaviour at the lower edge of the conformal window. The transition at the lower edge could be first order [26,27] or infinite order [4,22], and can be analyzed through several methods. This is an important open problem as characterization of this phase transition will lead to a better understanding of the strong interactions and might be relevant for model building beyond the standard model [27]. Furthermore, functional RG and holographic approaches indicate that the conformal window is continuously connected to the quark gluon plasma phase at finite temperature [28][29][30]. Therefore studies of the conformal window may provide crucial information on the quark gluon plasma and vice versa.
Several works in the literature [9,24,25,27,28] point toward the possibility that the IR stable fixed point α g = −b 0 /b 1 annihilates with another fixed point of the RG flow at x crit in the form of a saddle-node bifurcation (second option in the list above). This would in turn indicate that a UV irrelevant operator, such as a 4-fermi interaction, crosses marginality precisely at x crit [3] as already anticipated in the early literature [20,21,23]. In this paper, we investigate this possibility by considering an effective model for QCD 4 that includes four-fermi interactions and determine all fixed points within this model in and around the conformal window. Furthermore, we discuss how they appear and disappear as the bifurcation parameters N c and N f are varied. Finally we study the effects of adding a scalar operator with a Yukawa type interaction to our results in this effective theory. The paper is organized as follows. In the next section we introduce the effective model with the four-fermi interactions, derive the beta function equations both for finite values of N c and N f and in the Veneziano limit, and study the fixed point structure using bifurcation theory in detail. In section 4, we extend the model including the scalar operator. We conclude and discuss our results and list possible future work in the last section. Appendix A provides a short introduction to bifurcation theory in dynamical systems and Appendix B details the calculation of the betafunctions in perturbation theory.
2 An effective model for QCD 2

.1 Beta functions in an effective model for QCD
We study the theory in the Wilsonian picture. For this, we write down the effective action, which is given by the Legendre transform of the Lagrangian (1) with m i = 0. We impose the symmetries obeyed by the Lagrangian on the effective action and expand the effective action in the quark fields ψ ai up to the 4-quark interactions, as was done in [15]. This yields where L 4f denote the four-fermi interactions and we redefined the gauge fields One can write down four independent four-fermi operators [31][32][33][34] respecting the symmetries. We choose a basis as in [15,34]: where η = γ ψ is the anomalous dimension of the fermion field and One can calculate the beta functions of the four-fermi interactions using the Wetterich equation given in Appendix B.3. The beta function for the gauge coupling α g cannot be derived from the effective action (7) easily as one needs an expansion in the gluon fields in the effective action [15] which complicates the calculation substantially. Instead, we can start from the perturbative equation (2) and add the effective four-fermi interactions to it. We provide this calculation in appendix B.2 and obtain the same result as in [15]. The complete set of beta-functions is then given by 2 We find slightly different beta functions for the four-fermi interactions from the ones presented in [15]. They agree in the Veneziano limit however.
where we defined the rescaled variables Let us further define and perform the rescaling which makes the set (10) amenable to the Veneziano limit where N → ∞ and x remains finite. After the rescaling we obtain

Veneziano Limit
In this section, we study the Veneziano limit of the model (13) as was done in [3] and [15]. The Veneziano limit, where we have defined t = − ln(Λ/Λ 0 ). It is important to note that all higher order contributions in the perturbative theory are either suppressed or give constant contributions in the Veneziano limit with the rescaled variables. We now observe that the first three equations decouple from the last two. Therefore, we can reduce analysis to the system 3     α The domains of the variables and the parameters are given by x ∈ R + , α ∈ R + 0 and g S , g V ∈ R. Due to the fact that the model (10) is perturbative plus a few higher order corrections, we expect better results for small parameter values. As a rough bound we will use α g < 1, |g i | < 1, and notice when this bound is exceeded.
The beta function for α g has a double root for α t g = 0 and another root for α nt . Therefore, the fixed points of the RG-flow lie on these manifolds. The manifold defined by α g = 0 is an invariant set of the system, sinceα g vanishes for α g = 0.
The behavior of the fixed points projected on the (x, α g ), (x, g S ), (x, g V ) planes and (x, α g , g S ), (x, α g , g V ) and (x, g S , g V ) spaces are shown in Figure 2. We use the following color coding in this figure: • Solid red line: stable node, 3 negative eigenvalues (IR attractor).
The equilibria on the α g = 0 manifold all have a trivial eigenvalue, and therefore do not fit the color coding as described above. In order to use the color coding, we define the sign of a trivial eigenvalue on this manifold by approaching the equilibrium along a line of constant g S and g V from positive but small α g . The trivial eigenvalue will then approach 0 from either positive or negative values. If the trivial eigenvalue is slightly positive for 0 < α g 1, we define the trivial eigenvalue to be positive and vice versa. 4 In Figure 2, 4 up to 8 equilibria are shown, depending on the value of x. Furthermore, a few bifurcations are found. The parameter and coordinate values of these bifurcations are reported in Table 1. The Branching points all correspond to transcritical bifurcations. The transcritical bifurcation at x = 1 is located in the subspace α g = 0, while the other transcritical bifurcations are located on the intersection of the subspaces α g = 0 and α g = α nt g . In the latter case, both a trivial and a non-trivial equilibrium 5 are involved. At all the non-trivial equilibria, the equilibrium has a positive eigenvalue, when α g < 0, that becomes negative when α g > 0, while the trivial equilibrium has a negative eigenvalue that becomes positive. At the branching point at x = 5.5 a stable node with α g > 0 is created. This indicates the upper edge of the conformal window as discussed in the introduction of this section. Finally, at x = 4.05 two non-trivial equilibria collide and disappear through a saddle-node bifurcation (limit point). This saddle-node bifurcation could indicate the lower edge of the conformal window through a fixed point merger, as was discussed in [15] and [3]. The two other non-trivial fixed points diverge for x 2.5, close to the divergence of the large N c limit in equation (6). In Figure 3a, a phase portrait on the invariant plane α g = 0 is shown at x = 5. We find 4 fixed points, and 4 red critical heteroclinic orbits, connecting those fixed points, and marking the boundaries of an invariant set. All orbits inside this set have a finite UV and finite IR fixed point, while Figure 2: The fixed points of the RG-flow of (15) projected on the (x, α g )-plane, (x, g S )-plane, (x, g V )-plane, (x, α g , g S )-space, (x, α g , g V )-space and the (x, g S , g V )-space. LP: Limit Point. BP: Branching Point. orbits outside the invariant set have a diverging IR limit or diverging UV limit. The structure of this flow does not change, when x is varied in the region x > 1. Only the location of the fixed points and connecting orbits will change continuously. At x = 1, the two fixed points with smallest value of g V collide and exchange stability, via a transcritical bifurcation, yielding a similar phase portrait.
In Figure 3b, a phase portrait of the flow in the 3-dimensional space is shown at x = 5. Here, all 8 fixed points are shown, together with a few heteroclinic orbits connecting them. The red heteroclinic orbits form a skeleton of a 3-dimensional invariant set. This invariant set is bound by six invariant planes. The solid red lines lie on the intersections of these planes. The dashed red lines indicate the presence of other intersections, of which the exact locations haven't been found numerically. All orbits inside this set have a finite UV and IR fixed point, while orbits outside the invariant set have a diverging IR limit, a diverging UV limit or both. More precisely, the bulk of the invariant set consists of structurally stable orbits going from the equilibrium at (α g , g S , g V ) = (0, 0, −0.333) to the equilibrium at (α g , g S , g V ) = (0.032, 0.003, 0.000), and the boundary planes consist of nonstable heteroclinic orbits between the other fixed points. Within the invariant set the theory is in a chirally symmetric phase, while outside the set the chiral symmetry is broken. The structure of this flow does not change when x is varied in the region 4.05 < x < 5.5. Only the location of the fixed points and connecting orbits will change continuously. At x = 4.05 the two fixed points shown at (α g , g S , g V ) = (0.032, 0.003, 0.000) and (α g , g S , g V ) = (0.073, 0.322, −0.084) disappear through a saddle-node bifurcation.
In Figure 3c, a phase portrait in the whole space is shown at x = 4. Here, only 6 fixed points are left. The six fixed points together span two unstable bounded 2-dimensional invariant sets. The solid red lines border these sets, and the dashed red line indicates the location, where the last border is expected, which hasn't been found numerically. The 3-dimensional invariant set, which was present at x = 5, has been broken down at the saddle-node bifurcation to two 2-dimensional invariant planes. The structure of the flow does not change when x is varied in the region 0 < x < 4.05, but the two non-trivial fixed points diverge when x is lowered towards x = 2.5. The orbits outside of the invariant planes all diverge in the IR limit. The eigenvector with 0 eigenvalue at the saddle-node bifurcation 6 is mainly directed in the g S direction, as can be found in Table 2. Therefore, the divergence is the strongest in the g S direction, indicating that the chiral symmetry is broken due to the scalar part of the four-fermi interaction, as discussed in [35]. Furthermore, we notice that there still exist chirally symmetric trajectories with a zero UV and finite IR limit, but these are unstable, and therefore need a very precise fine-tuning within this model.

Fixed points and Codim-1 bifurcations
In the model out of the Veneziano limit there are 2 bifurcation parameters N c and x = N f /N c , which allows for more interesting types of bifurcations. Again, we are interested in the renormalized trajectories of the theory and the invariant sets enclosed by the renormalized trajectories, indicating the existence of chirally symmetric theories. We have numerically found that the model contains up to 20 different real fixed points with α g ≥ 0 for integer values of (N c , N f ) ∈ [1, 10] × [1, 100]. The beta function for α g yields two roots where the trivial root has multiplicity two. The (α g = 0)-hyperplane is an invariant set of the model, which allows to study the model without gauge interactions separately. Setting α g = 0, we find up to 12 fixed points, which can be continued using MatCont. Two of those disappear through a saddle-node bifurcation at x > 7.5, which is above the conformal window. The 10 remaining fixed points play an important role for the behavior for α g = 0. For every of these ten fixed points, we find a related fixed point on the manifold α g = α nt g . For large values of x all the non-trivial fixed points have α g < 0, but for smaller values of x they have α g > 0, and become physically relevant. On the intersection of the two hyperplanes, where α g = α nt g = 0, i.e. when x = 11 2+3g V , all non-trivial fixed points collide with a trivial fixed point, through a transcritical bifurcation. The lines at which these transcritical bifurcations occur are shown in Figure 4a. Using the transcritical bifurcations, we relate all the non-trivial fixed points to a trivial fixed point.
The 10 non-trivial fixed points with positive α g exist for large values of x. However, for small values, they either disappear through a saddle-node bifurcation or diverge in one or multiple couplings. In Figure 4b the saddle-node bifurcation lines are shown. Here, we see that two pairs (black lines) of fixed points will disappear through a saddle-node bifurcation which is connected to the saddle-node bifurcation that was found in the Veneziano limit 7 . The two bifurcation curves are close together for large N c , but split up at small N c . One of the curves starts to increase steeply around N c = 2, while the other starts decreasing around N c = 6, and diverges in the coupling constants at N c = 2. Therefore the saddle-node bifurcation of this pair does not exist for N c = 1, and the fixed points diverge instead. Then there is a third pair of fixed points (upper blue line) that disappears through a saddle-node bifurcation for all N c around x = 4.8. However this value starts to increase around N c = 4. A fourth pair disappears through a saddle-node bifurcation when 1 < N c < 16, diverges for N c ≥ 16, while at N c = 1 the fixed points will not exist for positive α g . The last two fixed points do not disappear through a saddle-node bifurcation, but will diverge in multiple couplings at small x.
(a) Transcritical bifurcations on the intersection of the (αg = 0) and (αg = α nt g ) subspaces. Black line: transcritical bifurcation marking the upper edge of the conformal window bifurcation curves.
(b) Saddle-node bifurcations in the (αg = α nt g ≥ 0)-subspace. Black lines: saddle-node bifurcation, that was also found in the Veneziano limit. Using these saddle-node bifurcation curves, we can partition the 10 non-trivial fixed points into 5 pairs. Since all 10 non-trivial fixed points are related to a trivial fixed point through a transcritical bifurcation curve, we can partition the 20 fixed points into groups of 4 consisting of 2 trivial and 2 non-trivial fixed points. It turns out that the heteroclinic orbits between the fixed points within each of these sets have a similar topology, which is as shown in Figure 5: at first, for large x, the two non-trivial fixed points lie below the α g = 0 hyperplane, and the α g = 0 hyperplane is locally attractive in the IR limit, then, at a specific value of x, one of the non-trivial and one of the trivial fixed points collide through a transcritical bifurcation and the non-trivial fixed point gets a positive α g creating a locally attractive IR limit set bounded by orbits between the 3 fixed points with α g ≥ 0. Next, at smaller specific x, the other non-trivial fixed point will collide with the other trivial fixed point and also get positive α g value. After this, at small x, the non-trivial fixed points will either disappear through a saddle-node bifurcation or diverge, making the α g = 0 hyperplane locally repulsive in the IR limit. We can also visualize the structure of the flow on the invariant (α g = 0)-space. A graphical representation of the fixed points and a few critical orbits is shown in Figure 6 for (N c , N f ) = (3,15) and (N c , N f ) = (3, 3). We notice that two points (labeled by 5 and 7) do exist for (N c , N f ) = (3, 15), but do not exist for (N c , N f ) = (3,3). This is due to a saddle-node bifurcation in the (α g = 0)-space at (N c , N f ) ≈ (3, 5.9). Furthermore two fixed points (labeled by 6 and 9) have exchanged stability through a transcritical bifurcation at (N c , N f ) ≈ (3, 7.7). The curves corresponding to these two bifurcations are shown in Figure 7e and 7a respectively. In the remainder of this section, we'll use the numbers as in Figure 6 to label the fixed points. Since we have related every trivial fixed point uniquely to a non-trivial fixed point, we give the non-trivial fixed point the same number as its trivial counterpart. However, we add the label 'a' to the trivial fixed points and the label 'b' to the non-trivial fixed points. We have numerically found all local bifurcations 8 involving these 20 fixed points for α g ≥ 0, and these are shown in Figure 7. This figure also includes the saddle-node and transcritical bifurcation curves shown in Figure 4.   Figure 6.
We use the bifurcation diagram 7 to distinguish and visualize several phases within the RG-flow generated by our model (10). For this, we set N c = 3, and take N f ∈ {11, 12, 13, 15}. The topology on the invariant (α g = 0)-space will be the same for these four values of x, since there are no bifurcations in this region. Only the location of the fixed points will change continuously. Projections of the fixed points and a few heteroclinic orbits in this 4-dimensional invariant subspace, at α g = 0, are shown in Figure 8. The fixed points and a few heteroclinic orbits in the α g = α nt g subspace are shown in Figures  At N f = 15, we see 10 fixed points and a few heteroclinic orbits connecting them, both in the α g = 0 and the α g = α nt g hyperplane. The red lines shown in the figures form part of a skeleton of a closed set in the 4-dimensional subspaces α g = 0 and α g = α nt g . These two 4-dimensional sets on α g = 0 and α g = α nt g form a 4-dimensional boundary of a larger 5-dimensional invariant subspace of the full 5-dimensional space. The other 4-dimensional subspaces that bound the full 5-dimensional invariant set are spanned by orbits going from α g = 0 towards α g = α nt g . Dashed red lines indicate parts of the skeleton that haven't been found numerically, but are expected to exist. The flow within the 5-dimensional invariant set is everywhere directed from the α g = 0 towards the α g = α nt g > 0 subspace. This set consists of structurally stable renormalized trajectories that have finite UV (with trivial α g ) and IR limits (with strictly positive α g ) indicating the theory to be in a symmetric phase. Outside the set, orbits have diverging IR and/or UV limit, and theories will have a broken chiral symmetry. The boundaries consist of unstable orbits between other fixed points, and therefore indicate unstable symmetric theories.
When the flavor number is reduced just below N f = 15, at N f = 14.7 and N f = 14.6, two pairs of fixed points disappear through a saddle-node bifurcation. Because of these bifurcations, the invariant set of orbits indicating symmetric theories is reduced. However, the most attractive IR fixed point with dim(W s ) = 5 still exists. This indicates that the invariant set is locally attractive in the IR limit. Reducing the flavor number even further to N f = 12 and N f = 11 more fixed points disappear through a saddle-node bifurcations. The invariant set of orbits at α g = α nt g is reduced at every step, and therefore the set of asymptotically stable orbits with finite UV and IR limit becomes smaller, making the set of possible symmetric theories smaller. For N f ≤ 12.2, there is no fixed point left with dim(W s ) = 5. Therefore, the remaining invariant set is not an IR attractor anymore. Furthermore, we expect that the dimension of the invariant set with finite trajectories reduces, when more fixed points have disappeared through a saddle-node bifurcation. When the invariant set splits up into multiple invariant sets with dimension smaller than 5, symmetric theories only exist within this model under extreme fine-tuning. This is a strong indication that the lower edge of the conformal window has been crossed.
In Figure 8, there is another peculiarity visible. Here we see a heteroclinic orbit with spiraling behavior indicating the existence of renormalized trajectories with complex scaling behavior in and close to the α g = 0 hyperplane. Equilibrium 5 and 7 are associated to this behavior.  (3,15) and α g = 0. Solid line: heteroclinic connections. Dashed lines: heteroclinic connections that are expected to exist, but haven't been found numerically. Left: projection on the (g S , g V , g V2 )-space. Right: projection on the (g V , g V1 , g V2 )-space. and α g = α nt g . Solid line: heteroclinic connections. Dashed lines: heteroclinic connections that are expected to exist, but haven't been found numerically. Left: projection on the (g S , g V , g V2 )-space. Right: projection on the (g V , g V1 , g V2 )-space.    (13) at (N c , N f ) = (3, 13) and α g = α nt g . Solid line: heteroclinic connections. Dashed lines: heteroclinic connections that are expected to exist, but haven't been found numerically. Left: projection on the (g S , g V , g V2 )space. Right: projection on the (g V , g V1 , g V2 )-space.  (13) at (N c , N f ) = (3, 12) and α g = α nt g . Left: projection on the (g S , g V , g V2 )-space. Right: projection on the (g V , g V1 , g V2 )-space. Figure 13: Fixed points and a critical heteroclinic orbit of the QCD 4 model (13) at (N c , N f ) = (3, 11) and α g = α nt g . Left: projection on the (g S , g V , g V2 )-space. Right: projection on the (g V , g V1 , g V2 )-space.

Operators crossing marginality and degenerate bifurcations
We have found 4 saddle-node bifurcation curves (cf. figure 4b, 7) with non-trivial value of α g . All four indicate the disappearance of two fixed points when the value of x is decreased. Along these curves we can identify the eigenvector corresponding to a zero eigenvalue as the relevant operator that crosses marginality. These critical eigenvectors at N c = 3 are reported in Table 3. We have found that these operators are very dependent on x, but not on N c . From the table we conclude that 3 of the saddle-node bifurcations have 2 complex eigenvalues at N c = 3, while one has only real eigenvalues. We notice that the normalized eigenvectors corresponding to the complex eigenvalues are all directed along the four-fermi couplings, and have very small contributions (|η αg | < 0.007) along the e αg direction, meaning that this complex scaling behavior involves operators, which are linear combinations of the effective four-fermi interactions. Furthermore, we see that the operator crossing marginality is a linear combination of all operators. However the scalar four-fermi interaction is most relevant in three of the saddle-node bifurcations, while in the saddle-node bifurcation of equilibria 5b and 7b the V 1 interaction is most relevant. In the complete model we have encountered nine bifurcation points with higher codimension. An overview can be found in Table 4. Six of these points lie in the physical uninteresting region N c ∈ (0, 1), while three of them (CP, SNT1 and IP1) may have physical relevance, since N c ≥ 1. Those points will be discussed below. The point labeled with 'CP' corresponds to a cusp bifurcation and the points labeled with 'SNT' correspond to saddle-node-transcritical bifurcations, as discussed in Appendix A. The point 'TCT' labels an intersection of two independent transcritical bifurcation curves and corresponds to the normal form The points 'IP' label more degenerate bifurcation points. Bifurcation At the cusp bifurcation 2 saddle-node bifurcation curves with non-trivial value of α g meet. The first curve is associated to the pair of equilibria (6b,8b) and the second to (6b,9b).
The point SNT1, describes a saddle-node-transcritical bifurcation, of a saddle-node curve that lies in the subspace α g = 0 and involves the equilibria (2a,3a), and two transcritical curves that also lie in the α g = 0 subspace and involve the equilibria (2a,4a) and (3a,4a). IP1 is a highly degenerate bifurcation point. The model reduced on the invariant subspace α g = 0 has a pitchfork bifurcation at this point. However, in the 5-dimensional model the branches of the equilibrium curves meeting at this pitchfork bifurcation are all transcritical bifurcation curves of the equilibria (6a,6b), (8a,8b), (9a,9b) and (6a,9b). In addition, to this pitchfork bifurcation we find a saddle-node-transcritical bifurcation at this point. The curves associated to this saddlenode-transcritical bifurcation all lie in the α g = 0 subspace. The saddle-node curve involves the equilibria (8a,9a) and the branches of the transcritical bifurcation involve (6a,8a) and (6a,9a). 9 On top of this, there is a saddle-node bifurcation curve with non-trivial α g intersecting at this point. This curve involves the equilibria (6b,8b). Together with the branches of transcritical bifurcation curves of (6a,6b) and (8a,8b), that were contained in the pitchfork bifurcation, we could also see this point as saddle-node-transcritical bifurcation of these points. In short this point is highly degenerate and can be unfolded into two saddle-node-transcritical bifurcation and a pitchfork bifurcation of transcritical bifurcation curves.

Discussion of the Results
In the complete model we can distinguish 10 pairs of fixed points, that may be relevant for the physics in and below the conformal window. These points can be divided into three different sets: {1, 2, 3, 4}, {5, 7} and {6, 8, 9, 10}. Here, we used the numbering as in Figure 6. The first set contains the most repulsive fixed points, while the last set contains the most attractive fixed points.
The set of pairs of fixed points {6, 8, 9, 10} consists of points that were also found in the Veneziano limit. The non-trivial points of the pairs 8 and 10 disappear through a saddle-node bifurcation at x ≈ 4 for large values of N c like in the Veneziano limit. Furthermore, the non-trivial points of the pairs 6 and 9 diverge for N c ≥ 16, and the trivial points of the pairs 6 and 9 have a transcritical bifurcation at x = 1 for large values of N c . This behavior was also found in the Veneziano limit. Equilibrium 10 is the most attractive and the point with non-trivial α g in this pair describes the IR limit of the structurally stable renormalized trajectories that exist within the invariant set. We see that like in the Veneziano limit the non-trivial IR point disappears through a saddle-node bifurcation with the non-trivial equilibrium of pair 8 at (N c , x) = (3, 4.06) or (N c , N f ) = (3, 12.2). This saddle-node bifurcation is a good candidate for the lower edge of the conformal window in this model. This prediction is higher than but comparable to the predictions summarized in [3]. Furthermore, it is close to a recent prediction on a similar model including effective interactions [34]. In addition, we see that the transcritical bifurcation of the equilibria of pair 10, where the most attractive equilibrium gets a positive α g is constant at the ratio x = 5.5. This line indicates the upper edge of the conformal window, and the same value was found in the Veneziano limit. This is consistent with earlier predictions, as was expected, since α g is small, so the two-loop perturbative beta function of α g without four-fermi interactions should give a good approximation. However, we have also found a few other non-trivial fixed points which already exist above the ratio x = 5.5 in our model. This might indicate that the non-perturbative beta function for α g contains fixed points for ratios above x = 5.5, which are not present in the 2-loop beta function.
If we take the saddle-node bifurcation of the non-trivial equilibria of the pairs 8 and 10 as the lower edge of the conformal window, we see in Figure 7a that the the ratio x, indicating the lower edge of the conformal window, increases for N c ≤ 2. Furthermore, we see that the two fixed points that were diverging in the Veneziano limit for small values of x will disappear through a saddlenode bifurcation for values 1 < N c < 16. This saddle-node bifurcation takes place for ratios below the conformal window for N c > 5, and for 1 < N c ≤ 5 it takes place in or above the conformal window, while at N c = 1 these fixed points do not exist for strictly positive α g . This saddle-node bifurcation curve affects the closed invariant set of structurally stable orbits, and could therefore be relevant for the physical behavior of QCD in the conformal window: for N c ≥ 16, we expect similar behavior as in the conformal window, where these two fixed points diverge around x ≈ 2.5. Here, the fixed points form extremal points of a 5-dimensional closed invariant set, which is discontinuously broken down to a smaller invariant set, when the saddle-node bifurcation of the pairs 8 and 10 takes place. On the other hand, for 5 < N c < 16, the points no longer diverge, but also disappear through a saddle-node bifurcation. However, this happens for a value of x which is below the lower edge of the conformal window. Therefore, the behavior within the conformal window won't be changed much by this saddle-node bifurcation. In the region 3 ≤ N c ≤ 5, this saddle-node bifurcation takes place within the conformal window, resulting in a discontinuous change in size of the closed invariant set of structurally stable orbits (or renormalized trajectories with finite IR and UV limit) at the bifurcation. I.e. it does not affect the range of the conformal window, but it induces an discontinuous change in the set of symmetric theories, which could indicate a transition within the conformal window in the sense that a non-empty set of finite renormalized trajectories, which are present near the upper edge of the conformal window, are not present close to the lower edge and have diverging IR limits instead. At N c = 2, the saddle-node bifurcation takes place above the conformal window, and therefore this discontinuous transition won't take place within the conformal window. However, the set of symmetric theories in the conformal window will be smaller than at larger values of N c . At N c = 1, this saddle-node bifurcation only takes place for negative α g , but not for positive α g . Furthermore, the lower edge of the conformal window has increased at N c = 1 making the conformal window smaller than at higher values.
The set of pairs of fixed points {1, 2, 3, 4} contains the four pairs of fixed points that are the most repulsive of the 10, and are therefore most relevant for the UV physics. Pair 1 is the most repulsive pair and therefore describes the UV limit of the renormalized trajectories within the invariant set of structurally stable orbits. At a curve x > 8 the equilibria within this pair undergo a transcritical bifurcation, as can be seen in Figures 7c. Below this curve the most repulsive of the pair is the trivial equilibrium making the theory asymptotically free. Above this curve the α g = 0 hyperplane is repulsive, and therefore the theory always diverges in the UV limit. The non-trivial fixed points of the pairs 1 and 2 both diverge when x is decreased. Furthermore, we see that two of the non-trivial equilibria, the ones from the pairs 3 and 4, disappear through a saddlenode bifurcation for N c ≥ 2 this bifurcation happens for smaller values of x than the saddle-node bifurcation which we took as the lower edge of the conformal window. We note that this bifurcation curve also converges to x ≈ 4 in the Veneziano limit, where the saddle-node bifurcation was found in the Veneziano limit. We find that both this curve and the saddle-node bifurcation curve of the non-trivial fixed points 8 and 10 correspond to the saddle-node bifurcation that was found in the Veneziano limit, since the beta functions for g V1 and g V2 have two fixed points at this point, which shows that there are actually 2 saddle-node bifurcations in the Veneziano limit, that differ in the values g V1 and g V2 . This saddle-node bifurcation curve of the non-trivial fixed points of the pairs 3 and 4 is another candidate for the lower edge of the conformal window. If we take the saddle-node bifurcation curve of points 3 and 4 as the lower edge of the conformal window, the saddle-node bifurcation curve would indicate a discontinuous change in the invariant set of structurally stable orbits indicating the symmetric theories. Furthermore, this would imply that the lower edge of the conformal window is just below N f = 12 at N c = 3, as was found in lattice studies such as [11], and from Schwinger-Dyson equations with ladder resummations in [19][20][21][22][23], where the lower edge was predicted at Finally, we notice that there is a last set of two pairs of fixed points, {5, 7}, which wasn't found in the Veneziano limit. These fixed points have complex eigenvalues with relatively large complex part compared to the real part. This induces a spiraling behavior in the system. These points disappear through a saddle-node bifurcation both on the α g = α nt g and the α g = 0 subspace. These fixed points might have a physical interpretation or might be a consequence of the approximate nature of our model. Furthermore, the other fixed points (except for the pairs 8 and 10) can have complex eigenvalues, but for those the complex part is smaller than the real part, making the spiraling behavior less visible. Also, the eigenvectors corresponding to the complex eigenvalues are largely directed into the direction of the effective four-fermi couplings. The complex eigenvalues indicate that when the orbits move away (towards) the fixed point there is both a scaling, indicated by the real part of the eigenvalue and a rotation, indicated by the complex part of the eigenvalue. The scaling determines the velocity of flowing away from (towards) the fixed point, while the rotation indicates that the relevant operator, which is a linear combination of the basis operators rotates in theory space. In our case this rotation only happens in the subspace spanned by effective operators, since the operators with complex eigenvectors are mostly directed along the effective operators.
In conclusion, we find two saddle-node bifurcations that could indicate the lower edge of the conformal window, in the sense that both involve the merger of a non-trivial fixed that is created through a transcritical bifurcation at x = 5.5, and thus corresponds to the pertubative Banks-Zaks fixed point. These two fixed points merge with other non-trivial fixed points, which are created at values x > 5.5 through transcritical bifurcations. These other fixed points could become apparent in higher order perturbation theory, when the beta function becomes polynomials of higher order in α g . The two possible lower edges of the conformal window are then given by (N c , N f ) = (3.00, 12.2) and (N c , N f ) = (3.00, 11.3). Apart from these 2 pairs non-trivial fixed points we find 3 more pairs of non-trivial fixed points that are all created through transcritical bifurcations at x > 5.5 and either disappear through a saddle node bifurcation or diverge. Furthermore, we find that most bifurcations occur at almost constant x(N ), for values of N ≥ 3.

QCD 4 with a Scalar Field
In this section we extend our previous model with a meson-like scalar field Φ that couples to the fermion fields with the coupling constant y. The field Φ is invariant under SU (N c ) and transforms in the adjoint representation of This extension is based on an extension mentioned in [34]. The Lagrangian is extended with a massless scalar field and a Yukawa coupling with the fermion fields yielding extra terms in the Lagrangian: The set of beta functions is then extended to a set of 6 differential equations in 6 variables and 2 parameters, which is calculated in appendix C: where we have defined α y = y 2 (4π) 2 . We rescale the couplings as in (12), together with the rescaling N c α y → α y , and we define the ratio x := N f Nc , and N := N c . We find

The Veneziano Limit
We take the Veneziano limit (N → ∞), which yields The first 4 equations decouple from the last two and we can reduce analysis to The domain of the variables and parameters is given by x ∈ R + , α g , α y ∈ R + 0 and g S , g V ∈ R. Due to the fact that the model (19) is perturbative plus a few higher order corrections, we expect better results for small parameter values. As a rough bound we will use α i < 1, |g i | < 1, and notice when this bound is exceeded.
The beta function for α g has a double root for α t g = 0 and another root for α nt g = 11−2x+3x 2 αy−3xg V 13x−34 . Therefore, the fixed points of the RG-flow lie on these manifolds. The manifold defined by α g = 0 is an invariant set of the system, sinceα g vanishes for α g = 0.
In addition, the beta function for α y has a root for α t y = 0 and another root for α nt y = 3αg+2g S 2+x . Therefore, the fixed points of the RG-flow lie on these manifolds. The manifold defined by α y = 0 is an invariant set of the system, sinceα y vanishes for α y = 0.
The location of the fixed points projected on the (x, α g )-plane, (x, α y )-plane, (x, g S )-plane and the (x, g V )-plane is shown in Figure 15, and the projections on the (x, α g , α y )-space, (x, α g , g S )space, (x, α g , g V )-space, (x, α y , g S )-space, (x, α g , g V )-space and the (x, g S , g V )-space in Figure 16. Here, we use the following color coding: • Solid red line: stable node, 4 negative eigenvalues (IR attractor). The equilibria on the α g = 0 manifold all have a trivial eigenvalue, and therefore do not fit the color coding as described above. In order to use the color coding, we define the sign of a trivial eigenvalue on this manifold by approaching the equilibrium along a line of constant α y , g S and g V from positive but small α g . The trivial eigenvalue will then approach 0 from either positive or negative values. If the trivial eigenvalue is slightly positive for 0 < α g 1, we define the trivial eigenvalue to be positive and vice versa. 10 Since β αy = 0 if α y = 0, we find the same fixed points and bifurcations as in the previous model without scalar field (cf. Figure 2). In addition, we now find five fixed points with non-trivial value of α y . We find two additional branching points and one additional limit point, as can be seen in Figure 15: The equilibria of the RG-flow of (21) projected on the (x, α g )-plane, (x, α y )-plane, (x, g S )-plane and the (x, g V )-plane. LP: Limit Point. BP: Branching Point. Figure 16: The equilibria of the RG-flow of (21) projected on the (x, α g , α y )-space, (x, α g , g S )space, (x, α g , g V )-space, (x, α y , g S )-space, (x, α y , g V )-space and the (x, g S , g V )-space. LP: Limit Point. BP: Branching Point. Table 5. Furthermore, at 3 transcritical bifurcations, that were found in the α y = 0 model, we now find three equilibria that intersect and exchange stability, corresponding to the normal form ẋ = αx + x 2 , y = αy + y 2 .
One of the additional fixed points branches from the transcritical bifurcation at x = 1.000 that was already found in the model without a scalar, and diverges rapidly in α y . In addition, one new fixed point branches from one of the equilibria with non-trivial value of α g at x = 3.926, which also diverges rapidly in α y . Furthermore, at the branching point at x = 8.173, an extra fixed point branches of with a non-trivial value of α g and α y , which diverges when x is lowered towards x ≈ 4. Another new fixed point with non-trivial α y exists in the scalar model. For this point α y → 0 if x → ∞ and α y → ∞ when x is lowered towards x ≈ 2. From this point another fixed point branches at x = 11.647. This point disappears trough a saddle-node bifurcation at x = 5.497. The other point that disappears at this saddle-node bifurcation branches at x = 5.5 from the trivial fixed point (α g , α y , g S , g V ) = (0, 0, 0, 0). This point is the stable node in the system. This stable node disappears through a saddle-node bifurcation at 5.497. We therefore find a conformal window with non-trivial for α y at x ∈ (5.497, 5.501), which is smaller than a prediction on a model without effective four-fermi interactions, where x ∈ (5.24, 5.5) [34]. Table 5: Bifurcations found in the Veneziano limit of QCD 4 with a scalar field (21).

Small N c regime
In this section, we analyze the complete model (19) outside the Veneziano limit. The beta function for α g has two roots: where the trivial root has multiplicity 2. The trivial root makes the α g = 0 hyperplane an invariant set of the model (19). The beta function for α y also has two roots: The trivial root makes the α y = 0 hyperplane an invariant set of the model (19). Since α y = 0 is an invariant hyperplane, we find the same bifurcations as in previous section. However, on top of those, there exist fixed points and bifurcations in the α y = 0 space. We find up to six fixed points with α g , α y = 0 and up to six fixed points with α g = 0, α y = 0, depending on the values of (N c , x). As in the previous section, we find that fixed points with non-trivial α g always collide on the intersection α nt g = 0 with a fixed point with trivial value α g = 0 through a transcritical bifurcation. This happens at where the second condition is a reduction to the model without scalar field. Similarly, fixed points with non-trivial α y always collide on the intersection α nt y = 0 with a fixed point with trivial value α y = 0 through a transcritical bifurcation. This happens at where only the positive root is a physically relevant solution. In this way, one can relate most fixed points with (α g > 0, α y > 0) to a fixed point with (α g = 0, α y > 0), a fixed point with (α g > 0, α y = 0) and a fixed point with (α g = 0, α y = 0). We find that the orbits between these points are always as indicated in Figure 17. Furthermore, a non-trivial fixed point 11 with positive α g and α y can appear through a transcritical bifurcation at (α g = α nt g = 0, α y > 0), at (α g > 0, α y = α nt y = 0) or at (α g = α nt g = 0, α y = α nt y = 0). These transcritical bifurcation curves are shown in Figure 18 by dash-dotted red lines, dash-dotted black lines and solid red lines respectively. When x decreases non-trivial fixed points will either diverge in one or multiple couplings or disappear through a saddle-node bifurcation with another fixed point. We find three saddle-node bifurcation curves with (α g = 0, α y = 0) as shown in Figure 18a, 18b and 18c. They're shown together with the related 12 saddle-node bifurcation curves that were found for α y = 0 in the previous section. In addition, a saddle-node bifurcation curve with (α g = 0, α y = 0) related to the saddle-node bifurcation in Figure 18c is shown. Figure 17: Topology of the RG-flow between 4 fixed points. Figure 18a shows the bifurcations of the sets of fixed points which were indicated as pairs 8 and 10 in the previous section. Here, we see the transcritical and saddle-node bifurcations that were found in the previous section. In addition, there is an extra transcritical bifurcation at the solid red line, where a non-trivial fixed point branches of. The solid red line converges to the branching point at x = 5.5 that was found in the Veneziano limit. Furthermore, we see a transcritical bifurcation at non-trivial value of α y (dash-dotted red line), where a non-trivial fixed point splits of. This transcritical bifurcation converges to the branching point found at x = 11.6 in the Veneziano limit. The two non-trivial fixed points disappear through a saddle-node bifurcation at the solid blue line, which is related to the blue dashed line representing the saddle-node bifurcation that was found in the α y = 0 space. The solid blue line converges to limit point that was found at x = 5.5 in the Veneziano limit. Figure 18b shows the bifurcations of the sets of fixed points which were indicated by the pairs 3 and 4 in the previous section. Here, we see the transcritical and saddle-node bifurcations that were found in the previous section. In addition, there is an extra transcritical bifurcation at the solid red line, where a non-trivial fixed point branches of. The solid red line converges to the branching point at x = 5.5 that was found in the Veneziano limit. Furthermore, we see a trancritical bifurcation at non-trivial value of α y (dash-dotted red line), where a non-trivial fixed point splits of. This transcritical bifurcation converges to the branching point found at x = 11.6 in the Veneziano limit. The two non-trivial fixed points disappear through a saddle-node bifurcation at the solid blue line, which is related to the blue dashed line representing the saddle-node bifurcation that was found in the α y = 0 space. The solid blue line converges to limit point that was found at x = 5.5 in the Veneziano limit. Figure 18c shows the bifurcations of the sets of fixed points which were indicated by the pairs 6 and 9 in the previous section. Here, we see the transcritical and saddle-node bifurcations that were found in the previous section. In addition there is an extra transcritical bifurcation at the solid red line, where a non-trivial fixed point branches of. The solid red line converges to the branching point at x = 7.4 that was found in the Veneziano limit. Furthermore, we see a transcritical bifurcation at non-trivial value of α g (dash-dotted black line), where a non-trivial fixed point splits of. This transcritical bifurcation diverges in the Veneziano limit. The two non-trivial fixed points disappear through a saddle-node bifurcation at the solid blue line, which is related to the blue dashed line representing the saddle-node bifurcation that was found in the α y = 0 space. The solid blue line intersects the black dash-dotted line at SNT6. Therefore only in a small interval N c ∈ (1.2, 1.5) the non-trivial fixed points disappear through a saddle-node bifurcation, while at larger N c the non-trivial fixed points that branch of from the solid red line and the black dash-dotted line diverge for small values of x. In addition, for these fixed points we find a saddle-node bifurcation curve not only in the α y = 0 space, as was the case in the previous figures, but also in the α g = 0 space (dash-dotted blue line). Figure 18d shows a few more transcritical bifurcation curves lying in the α nt g = 0 and/or α nt y = 0 subspace. Non-trivial fixed points that are related to these transcritical bifurcation curves, all diverge in one or multiple couplings for small x instead of disappearing through a saddle-node bifurcation. We notice the the black, solid red and dashed red line involve the sets of equilibria 1 and 2, as labeled in the previous section. The eigenvectors that cross zero at the saddle-node bifurcations shown in Figure 18a and 18b are shown in Table 6 for N c = 3 along with the eigenvalues of the other eigenvectors. We see that in this model the bifurcations are strongly triggered by the α g interaction. Furthermore, as in the previous, model we find two complex eigenvalues for the saddle-node bifurcation of point the nontrivial fixed points 3 and 4. Eigenvectors corresponding to the complex eigenvalues are directed along the four-fermion couplings, and have very small contributions (|η αg | < 0.004, |η αy | < 0.016) along the e αg and e αy direction, meaning that this complex scaling behavior involves operators, which are linear combinations of the effective four-fermi interactions. Furthermore, we find multiple bifurcation points in Figure 18, which are reported in Table 7. The first cusp bifurcation and the points SNT5, IP1 and IP2 were already found in the previous model. In addition, we find an extra cusp point in the saddle-node bifurcation curve in the α g = 0 space. This curve has another cusp bifurcation at IP1, making IP1 even more degenerate. The points TCT1 up to TCT6 all describe the intersection of two transcritical bifurcation curves, and can be described by the normal form Finally, IP3 describes a point where 2 saddle-node bifurcation curves and 2 branches of a transcritical bifurcation intersect, and can therefore be described as a saddle-node-transcritical bifurcation, where an additional saddle-node bifurcation intersects. Four fixed points are involved in this bifurcation.  (19) shown in Figure 18.

Discussion of the Results
In this section, we have seen how the addition of a scalar field alters the bifurcations found in the previous section. In particular, we would like to know whether we find similar behavior in the (α g > 0, α y > 0) space as we found in the (α g > 0, α y = 0) space. We find the transcritical bifurcation at x = 5.5, which marks the upper edge of the conformal window. This transcritical bifurcation now indicates the intersection of three fixed points: one on the (α g = 0, α y = 0) subspace, one on the (α g > 0, α y = 0) subspace, and one on the (α g > 0, α y > 0). Furthermore, at this line there are two different transcritical bifurcation curves for different values of (g S , g V , g V1 , g V2 ) as was also found in the previous model. One of these bifurcations is shown in Figure 18a and the other in Figure 18b. We see that both saddle-node bifurcations that could indicate the lower edge of the conformal window in the model with α y = 0 have a related saddle-node bifurcation in the (α g > 0, α y > 0) space, indicating that the phase transition at the lower edge of the conformal window is similar to the previous section and involves the merger of fixed points, followed by a walking behavior of the coupling constant, that changes continuously to a running behavior. However, the non-trivial saddle-node bifurcations are found at a larger value of x than the ones in the α y = 0 subspace, and are triggered by the α g interaction instead of the g S interaction. This curve is very close but slightly below x = 5.5 at large values of N c and slightly decreases towards smaller N c . Therefore, the conformal window is much smaller in the theory that contains a scalar field in the adjoint representation, although fixed points and renormalized trajectories between those points remain to exist for smaller values of x, but these are very unstable and require a large degree of fine-tuning within this model.
Furthermore, we find that the fixed points, which were labeled by 5 and 7 have no related fixed points with non-trivial value of α y . Therefore, we do not have the spiraling behavior related to those fixed points as we had in the previous section. For the saddle-node bifurcation of the pair labeled by 6 and 9 there is only a related bifurcation with non-trivial and positive α g and α y in a small interval N c ∈ (1.2, 1.5), which makes it irrelevant for the physical model, since we require integer values of N c . Therefore, we find that the non-trivial fixed points of the sets 1, 2, 6 and 9 all diverge for small x.

Conclusions and Outlook
We have evaluated an effective model for QCD 4 , and studied the fixed points and their bifurcations in and close to the conformal window. We have found a fixed point merger (saddle-node bifurcation) in the Veneziano limit at N f /N c = 4.049 indicating the lower edge of the conformal window. This prediction is similar to the ones by [3,15], where the same model was analyzed. Furthermore, we have continued this bifurcation to the small N c regime, and found that this curve is more or less constant up to (N c = 3.00, N f = 12.17). For N c ∈ {1, 2} the window is choked off. On the other hand, there exists another fixed point merger (corresponding to the same merger in the Veneziano limit), which is more or less constant up to (N c = 3.00, N f = 11.26). For N c = 1, these fixed points diverge instead of disappearing through a saddle-node bifurcation, while N c = 2 is an intermediate value where they merger at large coupling constants. Both mergers could indicate the lower edge of the conformal window. In both cases the transition at the lower edge of the conformal window is an infinite order phase transition, generated by the effective scalar four-fermi interaction for N c > 2. Furthermore, we have found a few extra sets of fixed points, that either disappear through saddle-node bifurcations or diverge.
The existence of multiple fixed points and multiple saddle-node bifurcations in this model, opens up many new questions about the lower edge of the conformal window. It might very well be the case that the infinite order perturbative beta function contains more fixed points apart from the well known Banks Zaks fixed point. This could result in a one or more transitions within the conformal window, where pairs of fixed points disappear through mergers. The lower edge is then reached when the last fixed point merger takes place. It would be interesting to try to relate such transitions to phase transitions that have been predicted in the quark gluon plasma. Furthermore, it would be interesting to extend the model with higher order effective interactions, and see which of the fixed points and saddle-node bifurcations remain to exist.
In this paper, we already extended the effective QCD model with an scalar field coupling through a Yukawa interaction to the fermions. This interaction turns out to destabilize the model. In this model we found that the 2 candidates for the conformal window both exist in this model as well. Although the lower edge of the conformal window is located at larger values of N f /N c for such a model, as was already predicted in [34]. The other fixed points on the other hand are pushed into a regime α g < 0 and/or α y < 0.
In addition, we have found that the model allows for complex scaling dimensions of the effective interactions. These complex scaling dimension contain both a scaling and a rotation of the associated operator. The rotation indicates that the (ir)relevant operators associated with the complex eigenvalue rotates in theory space, while the orbits move away from/towards the fixed point. Such behavour was early found for example in [36,37]. For one pair of fixed points in the effective QCD model these complex scaling dimensions induced a clearly visible spiraling behavior in the RG-flow. This pair disappeared the the model extended with a scalar field and Yukawa coupling.
It would be interesting to study these kind of exotic orbits in more detail, and maybe relate them to other exotic RG-flows that have been found, such as the bouncing solutions recently found in holographic studies [38] Finally, we have shown that the use of bifurcation analysis in the study of RG-flows can be very useful tool in mapping out the fixed points and their transition in complicated model with multiple variables and parameters as was suggested in [3]. It would be interesting to use numerical tools such as the MATLAB package MatCont [1,2] in other studies of renormalization group flows. In particular for non-perturbative beta functions this can lead to new insights. An alternative to our numerical approach would be to use Conley index theory and describe the exit sets for RG-flows as was suggested in [3].

A Dynamical systems and bifurcation theory A.1 Dynamical Systems
This appendix presents a short overview of relevant terminology from the field of dynamical systems, partially following [39]. A dynamical system is a tuple (T, X, Φ t ), where X is a state space, for time t ∈ T , Φ t : X → X is a function defined on the state space, and typically T = N 0 or T = N for discrete dynamical systems, and T = R + 0 or T = R for continuous dynamical systems. The collection of all maps {Φ t } t∈T is called the flow of the dynamical system and should satisfy Therefore, the flow forms a group, if Φ is invertible, and a semi-group otherwise.
For a point x ∈ X, one can define the orbit O Φ (x) := {Φ t (x)|t ∈ T )}. Furthermore, one can define the notion of an invariant set of a dynamical system.
Definition. Given a dynamical system (T, X, Φ t ), • a set that is both forward and backward invariant is called an invariant set of the dynamical system.
If the map Φ t is non-invertible and T is only defined on positive numbers, every forward invariant set is called an invariant set of the dynamical system. It is assumed above that Φ t (x) is defined for all (t, x) ∈ T × X. If it is not the case, i.e., given x ∈ X, Φ t (x) is only defined for a subset of T that includes t = 0, the dynamical system is called local, and all definitions should be modified accordingly.
A natural way to generate a dynamical system is by a set of ordinary differential equations (continuous-time dynamical system) or a set of difference equations (discrete-time dynamical system). In this article, we consider dynamical systems generated by a set of differential equations, since those can be related to the beta functions of the renormalization group. Given the space X = R n , one can consider the system of autonomous ordinary differential equationṡ where x ∈ R n and f : R n → R n , which is assumed to be (sufficiently) smooth. If one defines the map Φ t : R n → R n by Φ t (x 0 ) = x(t, x 0 ), then the tuple (R, R n , Φ) is a (local) continuous-time dynamical system. The solutions for specified initial conditions, x(0, x 0 ) = x 0 , define the orbits and the visualization of the flow is called the phase portrait. In the study of dynamical systems equilibria 13 and periodic, homoclinic, and heteroclinic orbits are of particular interest.
Definition. Given a continuous-time dynamical system (R, R n , Φ t ), then a point x ∈ R n is called Definition. Given a continuous-time dynamical system (R, R n , Φ t ), a point x ∈ R n , and an orbit O Φ (x), then the orbit is called • a periodic orbit if x is not an equilibrium and Φ T (y) = y ∀y ∈ O Φ (x) for some T > 0; minimal such T is called the period; The topology of the phase portrait is characterized by its invariant sets and their stability. An invariant set is asymptotically stable, if all orbits starting in a neighborhood of the invariant set stay in this neighborhood and converge towards the invariant set if t → ∞ and unstable otherwise. For generic equilibria the stability is determined by the eigenvalues of Jacobian matrix at the equilibrium. At the equilibrium, one can define the notions of a stable, unstable and critical eigenspace.
An equilibrium is called hyperbolic if dim(T c ) = 0, stable if dim(T s ) = n and unstable if dim(T u ) ≥ 1. Furthermore, an equilibrium with dim(T s ) = n or dim(T u ) = n is a node and an equilibrium with dim(T s ) ≥ 1 and dim(T u ) ≥ 1 is a saddle. The classification of non-hyperbolic equilibria not satisfying these conditions is more involved, and won't be discussed further in this article. For a generic equilibrium in a continuous-time dynamical system, one can also define the notions of a stable and unstable manifold of the equilibrium.
Definition. Given a continuous-time dynamical system (R, X, Φ t ), and an equilibrium point x ∈ X of Φ, then • the stable manifold of x is defined as W s (Φ, x) = {y ∈ X : lim t→∞ Φ t (y) = x}, and • the unstable manifold as W u (Φ, x) = {y ∈ X : lim t→−∞ Φ t (y) = x}. For a dynamical system defined by (24), the stable manifold is at the equilibrium tangent to the stable eigenspace, while the unstable manifold is at the equilibrium tangent to the unstable eigenspace. For nonhyperbolic equilibria, one can introduce the concept of a center manifold W c loc (Φ, x), which is n c -dimensional and tangent to the critical eigenspace at the equilibrium. The existence of this locally defined manifold is implied by the Center Manifold Theorem.

A.2 Bifurcation Analysis
A continuous-time dynamical system defined by a set of ordinary differential equations can be analyzed by solving the ODE's numerically. However, if the differential equations also depend on a set of parameters, f α (x), α ∈ R m , one is often not interested in the exact solutions for specific values of α, but in the topological properties of the phase portrait such as the number of equilibria, limit cycles and their stability or the existence of homo-and heteroclinic orbits and chaotic attractors. In general, changes of the phase portrait, when α is varied, occur smoothly except for some specific values of α, where the topology of the phase portrait changes. Changes in the topology of the phase portrait are called bifurcations. Topology changes often occur near or within invariant sets. Bifurcations can be divided in 2 subclasses: • local bifurcations: bifurcations that can be detected by considering small but non-shrinking neighborhoods of equilibria or limit cycles.
• global bifurcations: bifurcations that cannot be detected in this way. Topology changes typically involve multiple invariant sets.
In addition, bifurcations may be called subcritical or supercritical, depending on the stability of the equilibria or limit cycles that appear or disappear at the bifurcation. Furthermore, bifurcations can be characterized by their codimension. The codimension of a bifurcation in a generic system is given by the number of independent conditions that have to be satisfied for a bifurcation to happen, and for local bifurcations is often related to the number of critical eigenvalues of the Jacobian matrix at the equilibrium. Consequently, the codimension of a bifurcation is always smaller or equal to the dimension of the parameter space. In general, many types of bifurcations can be found. For an overview of possible bifurcations and their conditions, we refer to [40]. In this appendix, we discuss bifurcations and related terminology that is relevant for this article. In doing so, we closely follow [40]. Bifurcations in low dimensional systems and with low codimension can often be found using analytical techniques, but for higher dimensional systems and in particular for higher codimension of the bifurcation one often has to rely on numerical techniques. Several computer programs have been developed for bifurcation analysis. In this article we make use of the Matlab package Matcont [1,2]. Bifurcation software relies on the conditions that have to be satisfied in order for a bifurcation to happen.
In different dynamical systems one can expect different bifurcations, which are not immediately defined by the same conditions. However, one would like to have a way of categorizing bifurcation in various systems. For this, we need the notion of topological equivalence. Consider two continuoustime dynamical systems:ẋ = f (x, α), x ∈ R n , α ∈ R m ; (25) y = g(y, β), y ∈ R n , β ∈ R m .
Definition. Dynamical systems (25) and (26) are topologically equivalent if • there exists a homeomorphism of the parameter space p : R m → R m s.t. β = p(α); • there exists a parameter-dependent homeomorphism of the phase space h α : R n → R n s.t. y = h α (x), mapping orbits of the system (25) at parameter values α onto orbits of the system (26) at parameter values β = p(α), preserving the direction of time.
Since many bifurcations are locally defined, one is often only interested in local topological equivalence near specific points.
Definition. Dynamical systems (25) and (26) are locally topologically equivalent near the origin, if there exists a map (x, α) → (h α (x), p(α)), defined in a small neighborhood of (x, α) = (0, 0) ∈ R n × R m such that • p : R m → R m is a homeomorphism defined in a small neighbourhood of α = 0, p(0) = 0; • h α : R n → R n is a parameter-dependent homeomorphism defined in a small neighborhood U of x = 0, h 0 (0) = 0, and mapping orbits of the system (25) in U onto orbits of the system (26) in h α (U ), preserving the direction of time.
By coordinate translations one can easily generalize this definition to local topological equivalence near arbitrary points. Having these definitions, one can categorize bifurcations by topological equivalence to normal forms. Consider a polynomial continuous-time dynamical systeṁ which has an equilibrium z = 0 at β = 0, satisfying k bifurcation conditions. Here, σ is a vector of the coefficients of the polynomial g(z, β; σ). (27) is a topological normal form for the bifurcation if any generic system (25) with equilibrium x 0 satisfying the bifurcation conditions α 0 is locally topologically equivalent near (x 0 , α 0 ) to (27) for some values of the coefficients σ i .

Definition. System
We notice that the system (25) can have higher dimension than is needed for the bifurcation to occur: x ∈ R n , while the normal form has z c ∈ R nc with n c < n. In such a case, there is a continuation of the critical center manifold. For nearby parameter values (25) is locally topologically equivalent toż where g(z c , β; σ) is a normal form on the center manifold and z s ∈ R ns , z u ∈ R nu correspond to the stable and unstable eigenspaces.
We now give a short discussion of local bifurcations which are encountered throughout the article. We begin with bifurcations in scalar one-parameter ODEṡ In the n-dimensional case, each normal form below is the normal form (28) on the one-dimensional center manifold.

Saddle-Node Bifurcation
This bifurcation is a codimension one bifurcation at which two equilibria collide and disappear, and is also known as fold or limit point bifurcation. The bifurcation occurs at (x,ᾱ), if the following conditions hold for the system (31): where a is the quadratic coefficient of the system (31). The system is then locally topologically equivalent near (x,ᾱ) to the normal forṁ near (y, β) = (0, 0) (see Figure 19). β y Figure 19: Saddle-Node bifurcation in the systemẏ = β − y 2 .

Transcritical Bifurcation
This bifurcation is a codimension one bifurcation in the class of systems that always have a trivial equilibrium, and is also a branching point: At the bifurcation, two equilibria collide and exchange their stability. The bifurcation occurs at (x,ᾱ), if the following conditions hold for the system (31): The system is then locally topologically equivalent near (x,ᾱ) to the normal forṁ y = βy ± y 2 , y, β ∈ R near (y, β) = (0, 0) (see Figure 20). Notice that we can transform the normal form of the tanscritical β y Figure 20: Transcritical bifurcation in the systemẏ = βy − y 2 .
bifurcation to the normal form of the saddle-node bifurcation by the non-invertible coordinate transformation This is called a nonversal unfolding of the saddle-node bifurcation.

Pitchfork Bifurcation
This bifurcation is a codimension one bifurcation for systems with reflectional symmetry, and is also a branching point: At the bifurcation one equilibrium splits into three equilibria. The bifurcation occurs at (x,ᾱ), if the following conditions hold for the system (31): • f xα (x,ᾱ) = 0.
The system is then locally topologically equivalent near (x,ᾱ) to the normal forṁ y = βy ± y 3 , y, β ∈ R near (y, β) = (0, 0) (see Figure 21). The − sign gives the supercritical case where a stable equilibrium splits into two stable equilibria and an unstable equilibrium, and the + sign is the subcritical case, where an unstable equilibrium splits into two unstable equilibria and a stable equilibrium. y β Figure 21: Supercritical pitchfork bifurcation in the systemẏ = βy − y 3 .
We continue with two bifurcations in scalar two-parameter ODEṡ In the n-dimensional case, each normal form below is the normal form (28) on the corresponding one-dimensional center manifold.

Cusp Bifurcation
This bifurcation is a codimension two bifurcation in generic systems, where two branches of saddle-node bifurcation meet tangentially. Nearby, the system can have three equilibria which collide and disappear at the saddle-node bifurcations. The bifurcation occurs at (x,ᾱ), if the following conditions hold for the system (35): where c is the cubic coefficient of the system (35). The system is then locally topologically equivalent near (x,ᾱ) to the normal forṁ near (y, β) = (0, 0) (see Figure 22).

Saddle-Node-Transcritical Bifurcation
This bifurcation is a codimension two bifurcation in the class of systems that always have a trivial equilibrium, where two branches of a transcritical bifurcation curve meet on a saddle-node bifurcation curve. There exist two types of this bifurcation. One can be seen as a nonversal unfolding of a Cusp bifurcation, while the other can be obtained as an unfolding of a degenerate Bogdanov-Takens bifurcation. Here, we discuss only the first case.   Figure 23: Saddle-node-transcritical bifurcation in the systemẏ = β 1 y + β 2 y 2 + y 3 . Nearby, the system can have three equilibria: One pair has a transcritical bifurcation at one of the branches and another pair has a transcritical bifurcation at another branch. The last possible pair of equilibria collides and disappears at the saddle-node bifurcation curve. If the bifurcation occurs at (x,ᾱ), then the system (35) is locally topologically equivalent near (x,ᾱ) to the normal form given byẏ = β 1 y + β 2 y 2 ± y 3 , y, β 1 , β 2 ∈ R near (y, β) = (0, 0) (see Figure 23). We then find a transcritical bifurcation along the line β 1 = 0, and a saddle-node bifurcation along β 1 = ± β 2 2 4 . The normal form can be transformed to the normal form of the cusp bifurcation through the transformation Finally, we discuss for completeness the last codimension one local bifurcation in generic ODEs (25), i.e. Andronov-Hopf bifurcation. Note that it does not occur in our RG-flows but might appear in other models. For this bifurcation, the center manifold is two-dimensional, so it is sufficient to consider (25) with n = 2 and m = 1.

Andronov-Hopf Bifurcation
This bifurcation is a codimension one bifurcation, and implies the birth of a periodic orbit (limit cycle). The bifurcation occurs at an equilibrium (x,ᾱ), if a pair of two conjugate complex eigenvalues of the Jacobian matrix crosses the imaginary axis. Generically, at the bifurcation the system is locally topologically equivalent near (x,ᾱ) to the normal forṁ y 1 = βy 1 − y 2 ± y 1 y 2 1 + y 2 2 , y 2 = y 1 + βy 2 ± y 2 y 2 1 + y 2 2 , near (y 1 , y 2 , β) = (0, 0, 0). The plus sign corresponds to the subcritical case, where an unstable limit cycles is created and the minus sign to the supercritical case, where a stable limit cycle is created (see Figure 24).
B Beta functions in QCD 4

B.1 Feynman Rules for QCD 4 with a Four-Fermi Interaction
The theory as discussed in section 2 has three propagators: 14 The fermion-gluon interaction (gauge vertex) is given by: Furthermore, there are four four-fermi interactions: Notice that the vertices have been split up in two parts. This is done to make clear along which lines the color charge is conserved, and will make the calculations in appendix B.3 more tractable, but has no further physical meaning.

B.2 Beta Function for the Gauge Vertex
In this appendix, we reproduce the method used in [15] to calculate the beta function for the gauge vertex. We're only interested in gauge invariant contributions to the beta function, but the Wetterich equation does not necessarily respect gauge invariance. Therefore, we do not use the exact renormalization approach. Instead we start from the well known 2-loop perturbative beta function for QCD 4 with N c colors and N f massless flavors [5]: Next, we would like to include effective interactions induced by the four-fermi couplings. These can be added perturbatively to the beta functions. The beta function is then given by where the last term is represents the perturbative corrections induced by the effective couplings. These perturbations should be chosen such that they're gauge invariant and such that they represent perturbations, that are included in lim n→∞ β g , but not in β [2] g . There is one one-loop correction induced by the effective four-fermi interaction, which is represented by = However, this correction turns out to be gauge dependent [31], and therefore we discard it. A possible two-loop correction is given by = For the effective four-fermi interaction one could take O i , i ∈ {S, V, V 1 , V 2 }. However, since the interactions should represent perturbations, which are included in lim n→∞ β [n] g the chirality should be the same in the whole diagram. This excludes the interactions O S and O V1 . Furthermore, due to the tracelessness of the generators of SU (N ) the contribution due to O V2 vanishes. Therefore, only O V contributes, and its contribution is in the sharp cutoff limit given by yielding a contribution to the beta function of the form
Using the Feynman rules from appendix B.1, we find the following contributing 1PI diagrams: = + + There are 29 different diagrams of the first type, 12 of the second type and 4 of the third. Five diagrams of the first type contribute to the beta function of G S . One of those is proportional to g 2 S : One diagram is proportional to g S g V : which is evaluated as 15 Two diagrams are proportional to g S g V1 : which is evaluated as G S G V2 2π 2 Λ 2 O S . Five diagrams of the first type contribute to the beta function of G V . One of those is proportional to g 2 S : Two diagrams are proportional to g 2 V : Two diagrams are proportional to g V g V2 : Nine diagrams of the first type contribute to the beta function of G V1 . One of those is proportional to g 2 S : Two diagrams are proportional to g 2 V1 : One diagram is proportional to g S g V2 : . Two diagrams are proportional to g V g V1 : . Ten diagrams of the first type contribute to the beta function of G V2 . Two of those are proportional to g 2 V : which are evaluated as One diagram is proportional to g 2 V1 : which is evaluated as Four diagrams are proportional to g 2 V2 : which are evaluated as One diagram is proportional g S g V1 : . Three diagrams of the second type are proportional to g S α g and contribute to the beta functions of G S and G V1 : The second and the third both evaluate to 0 and the first is evaluated as Three diagrams of the second type are proportional to g V1 α g and contribute to the beta functions of G S and G V1 : The first and second both evaluate to 0 and the third evaluates to Three diagrams of the second type are proportional to g V α g and contribute to the beta functions of G V and G V2 : The first and third evaluate to 0 and the second is evaluated as Three diagrams of the second type are proportional to g V2 α g and contribute to the beta functions of G V and G V2 :L ai L aiL bj L dj +L ai L cī L bj L dj L bj L dj L ai L ci + (L → R).
The first and third evaluate to 0, while the second is evaluated as Furthermore, the third type of correction is proportional to α 2 g and consists of 2 diagrams contributing to the beta functions of G S and G V1 : which are evaluated as and two contributing to the beta functions of G V and G V2 : which are evaluated as 9 32

C Beta Functions for QCD 4 with a scalar field
When the scalar terms (18) are added to the model (10), the set of beta functions is changed. As a consequence, new Feynman rules are added: First, we consider the function β αg . Due to the addition of the scalar field terms, the two loop beta function for β αg is changed to [44] β [2] g ≡ Λ dα g dΛ = −2b 0 α 2 g − 2b 1 α 3 g − 2N 2 f α y α 2 g .
The additional terms due to the four-fermion interactions can then be added as was done in appendix B.2. The beta function β αy is given by [45] Λ dα y dΛ = 2α y γ φ + γψ ψ .
The anomalous dimension of the scalar field is given by 2N c α y [44]. Using the four-fermi interactions, one can find a non-perturbative approximation to the anomalous dimension of the mass of the fermion field [34] such that we get three contributing diagrams + + .
The first two diagrams contribute N f α y and −6 N 2 c −1 2Nc α g , and two diagrams can be found that give a contribution to the second diagram: which contribute −2N c g S and 8g V1 .
Next, we consider the beta functions four the four-fermi interactions, which are a slight variation of the ones found for the model without scalars in appendix B.3. The anomalous dimension of the fermion fields is changed due to the presence of the scalar fields to η = α y [34]. Furthermore, we find additional contributions to represented by three types of diagrams + + .
Diagrams represented by the third type all evaluate to 0. For the first type we find four diagrams that each contribute to one of the beta functions of the four-fermi interactions: contributing respectively 4g V α y to β g S , g S α y to β g V , −2g V2 α y to β g V 1 and −2g V1 α y to β g V 2 .
Furthermore, for the second type there are 2 contributing diagrams: contributing 2α 2 y to β g V 1 and 2α 2 y to β g V 2 respectively.