Minima of Classically Scale-Invariant Potentials

We propose a new formalism to analyse the extremum structure of scale-invariant effective potentials. The problem is stated in a compact matrix form, used to derive general expressions for the stationary point equation and the mass matrix of a multi-field RG-improved effective potential. Our method improves on (but is not limited to) the Gildener-Weinberg approximation and identifies a set of conditions that signal the presence of a radiative minimum. When the conditions are satisfied at different scales, or in different subspaces of the field space, the effective potential has more than one radiative minimum. We illustrate the method through simple examples and study in detail a Standard-Model-like scenario where the potential admits two radiative minima. Whereas we mostly concentrate on biquadratic potentials, our results carry over to the general case by using tensor algebra.

The essence of scale invariance is that the interactions between fields are also responsible for the emergence of all of the observed scales [55,56]. Remarkably, this is a natural achievement of quantum theories, where the fluctuations need not respect scale symmetry even when the action contains no explicit scale.
The appearance of divergences in loop diagrams, for instance, forces the introduction of an arbitrary mass scale in the problem, required for the regularisation of the problematic contributions. Because observables cannot depend on an arbitrary parameter, the explicit dependence of a process on this quantity must be compensated by an implicit dependence of the involved physical parameters. This is the essence of the Callan-Symanzik equation, which determines the evolution of couplings, masses and fields with the renormalisation scale, resulting in the loss of scale symmetry at the quantum level. That quantum fluctuations generally violate scale invariance is also demonstrated by the mechanism of dimensional transmutation [55]. In this case, the quantum corrections induce the spontaneous symmetry breaking of the theory by producing a radiatively generated minimum in an otherwise trivial scalar potential. Because the minimisation equation relates the resulting vacuum expectation value (VEV) of the scalar field to the quartic coupling of its potential, the mechanism effectively allows to trade one parameter for the other.
In this paper, we concern ourselves with the impact of quantum corrections on the extremum structure of scale-invariant potentials, proposing a new mathematical framework to gauge their contribution. The central object of our analysis is the effective potential of a generic scale-invariant theory, whose minima determine the vacuum state and the spontaneous breaking of symmetries. The computation of the minima of the effective potential is generally a cumbersome task: On top of the dimensionality of the problem that scales with the number of scalar fields, the analysis is complicated by the nature of the radiative corrections. In the absence of scale invariance, such contributions seldom modify the properties of the theory determined by the dominant tree-level term. Within scale-invariant models, instead, these corrections often constitute the leading order contribution. In multi-scalar models, the analysis of the effective potential is also hindered by potentially large quantum corrections which depend on different mass scales. These enter the problem through logarithms of the ratios of mass terms to a common renormalisation scale. Consequently, minimising all quantum corrections in the region of interest by choosing a suitable value for this parameter is a straightforward matter only when the effective potential exhibits a single scale. More involved scenarios must then rely on different approaches to ensure that perturbativity is retained [57][58][59][60][61][62][63][64][65]. Scale-invariant theories are not immune to this problem, as several mass scales could result from non-degenerate field-dependent scalar, vector or fermion masses. To cope with the difficulties posed by several different scales, we use the renormalisation group (RG) improvement method proposed in [66,67], but our formalism is not limited to this.
In what follows, we thus present a new approach to the problem of radiative minima of a classically scale-invariant potential based on the method of [68] (successfully used in [69]). We put particular emphasis on biquadratic potentials because they comprise the lion's share of models in the literature. That the Standard Model Higgs boson transforms as an SU (2) L doublet quite naturally leads to biquadratic potentials when considering portal couplings to new scalar fields that transform non-trivially under gauge transformations. Other possibilities encompass Z 2 mirror symmetries, O(n)-symmetric and biconical theories. Our strategy is to write a biquadratic potential in terms of the matrix of quartic couplings and the scalar field vector, using a simple formalism that adds to understanding by building on the geometric intuition and is easily implemented with any modern computer algebra system. As a first result, we obtain general expressions for the stationary point equation, its solution and the mass matrix in a compact way that does not depend on the number of scalar fields. The minimum equation is solved via a systematic iteration procedure, whose lowest order yields the Gildener-Weinberg approximation. Not limited to biquadratic potentials, most of our conclusions also apply to generic potentials via tensor algebra.
Our second main result concerns scenarios where the classically scale-invariant effective potential admits multiple radiatively generated minima. Although the possibility was originally mentioned in Ref. [56], so far it has been practically neglected. In particular, we determine conditions that the coupling matrix has to satisfy in a minimum and identify semi-analytical criteria that assess the presence of additional minima in the effective potential. We explicitly demonstrate the possibility by using a toy model that incorporates two complex scalar fields, two Weyl fermions and an SU (2) gauge symmetry. Because our toy model is essentially a bare-bones version of the Standard Model extended by a new scalar and the scalar potentials in most such models are biquadratic, we expect our results to be relevant for future studies of cosmological phase transitions and gravitational wave signals in scale-invariant models.
The paper is organised as follows. We start by surveying general and biquadratic tree-level scalar potentials and their properties in section 2. The effective potential and its renormalisation group improvement is reviewed in section 3. In section 4 we discuss the conditions for the effective potential to develop a minimum. Analytical and numerical examples of minimum solutions are presented in section 5, while in section 6 we present detailed scenarios and discuss the criteria for an effective potential to have multiple nontrivial minima. We conclude in section 7.
2 Tree-level potential The generic renormalisable quartic potential of n real scalar fields φ i , collected in the vector Φ, is given by where λ ijkl is the symmetrised tensor of quartic couplings. 1 For biquadratic potentials, the expression simplifies to where Λ is the symmetric quartic coupling matrix. The • symbol denotes the Hadamard product, defined as the element-wise product of matrices of same dimensions: For example, the elements of the Hadamard square of the field vector, Φ •2 , are given by φ 2 i . Biquadratic potentials (2.2) admit a Z n 2 group, so each scalar φ i is odd under a private (Z 2 ) i symmetry. The radial coordinate ϕ in the field space is given by the norm of Φ, which can be expressed as where e = (1, . . . , 1) T , the vector of ones, is an identity element of the Hadamard product. Note that the relation (2.3) is quadratic in the field Φ but only linear in its square Φ •2 . We can use the vector e to construct identity elements of the Hadamard products of higher dimension, for instance the matrix ee T with elements all equal to unity. The tree-level field-dependent scalar mass matrix for the biquadratic potential (2.2) is given by where diag(v) designates a diagonal matrix with entries given by the vector v. For the first term, we have used ∇ Φ Φ T = I.

Effective potential
Considering one-loop quantum corrections, the effective potential can be written as where V (0) is the tree-level potential and V (1) is the radiative correction. Using dimensional regularisation in the MS scheme, the one-loop contribution is given by where µ is the renormalisation scale. The supertrace ranges over tree-level field-dependent scalar, vector and fermion mass matrices m 2 S,V,F . The matrix C, peculiar to the MS scheme, is a block-diagonal constant matrix with diagonal elements equal to c S,F = 3 2 for scalars and fermions, and c V = 5 6 for gauge bosons. The one-loop correction (3.2) can also be rewritten in the form [66] [66]. By construction, the scalar and vector boson terms in eq. (3.5) are always positive, while the fermion term is negative. Notice that when computed at the minimum of the potential, the tree-level mass of the scale-invariance Goldstone is negative and results in a spurious imaginary part of the logarithm. The exact solution to this problem requires a shift of the Goldstone mass m 2 G → m 2 G + ∆ [73,74] (see also e.g. [75]), but the numerical impact of the spurious part -in particular at the one-loop level -is negligible. In practice, it is then sufficient to consider only the real part of the logarithm, which corresponds to taking absolute values of the eigenvalues of the scalar mass matrix. Regardless, we will also show how to implement the exact solution in our formalism.

Renormalisation group
The renormalisation scale µ is arbitrary and the effective potential must not depend on it. Therefore, changes in µ are to be counterbalanced by changes in the masses, couplings and fields, as formalised in the Callan-Symanzik equation where g i include all mass parameters and coupling constants and γ is the matrix of anomalous dimensions of the scalar fields. The functions β i = µ dg i /dµ describe the running of couplings and masses with renormalisation scale.
In the case of a biquadratic tree-level potential, the β-function of the coupling matrix Λ is given by a matrix β. The Callan-Symanzik equation (3.6) implies that B, β and γ are related to each other via At one-loop level, γ = 0 in the absence of gauge and Yukawa interactions. As shown in appendix A, the matrix β-function is given by + gauge and Yukawa contributions, (3.8) where Diag(Λ) ≡ Λ • I is the diagonal matrix with the same diagonal as Λ. Conditions for Yukawa interactions to preserve the symmetries of biquadratic potentials are discussed in appendix B and a brief summary of gauge contributions is given in appendix C. Our result (3.8) for biquadratic potentials agrees with general results in the literature [76][77][78][79].

Improved potential
The validity of the perturbative expansion can be improved by taking a different renormalisation scale at each point of the field space. A common choice to RG-improve the effective potential sets µ = ϕ, or t = ln(ϕ/ϕ 0 ), which cancels only the second term in the one-loop radiative correction (3.3), yielding 2 Following the analysis in [66], we instead choose a field-dependent renormalisation scale µ * such that the one-loop part V (1) of the effective potential (3.1) vanishes (see appendix F for improvement procedures that use a different renormalisation scale). The effective potential then acquires the form of the tree-level potential, with running couplings in Λ(t * ) given by the solutions to the β-functions. The variable t * quantifies the displacement along the RG flow from the initial conditions set at a scale µ and is formally given by Notice that we could absorb the last term in (3.9) into the quartic coupling matrix by writing Λ eff ≡ That is, all quartic couplings are shifted by A/ϕ 4 . The last term can also be absorbed into the quartic coupling tensor in a general potential (2.1), but contrary to the biquadratic case, no unique prescription can be given.
or, through eq. (3.3), by An expression for t * can be found by solving V (1) = 0 numerically at a given point of field space, yielding µ * = µe t * . (3.13) To a first approximation, t * can be written as where the subscript in A 0 and B 0 indicates that the quantity is evaluated for couplings held at the fixed scale µ, neglecting the running. Because both A 0 and B 0 are homogeneous functions of order four, the first term in (3.14) depends only on the field vector length ϕ, while the second term is purely a function of the direction of the field vector. This distinction holds only at the level of t (0) * , since the second term in t * also depends on ϕ through the running couplings.
With the t (0) * approximation, the method ceases to be valid on the B 0 = 0 hypersurface which, being tangent to the RG flow, is a characteristic hypersurface of the RG equations [66]. Conversely, with t * , the B = 0 hypersurface is not problematic because the running of couplings in B regularises the divergence that occurs in the approximate expression. In order to maintain analytic control, in the following we employ the approximationdenoting henceforth t (0) * with t for simplicity -taking care of defining t (0) * in a subspace where B 0 = 0. In particular, after solving the RG equation for arbitrary initial conditions, we use the RG invariance of the effective potential to initialise the parameter in the vicinity of the radiative minimum, where B 0 > 0 holds as we show in section 4.4. Since B 0 depends on the field values but does not depend on the running couplings, this quantity necessarily remains positive everywhere but at the origin, the denominator in t is finite and the RG flow well defined.
Notice that t (0) * coincides with t * at the scale µ where these parameters are defined. Consequently, with t (0) * specified as above, the approximation holds in the vicinity of the radiative minimum as the quantum corrections in V (1) necessarily remain small in the field space region of interest. Although this choice guarantees the precision required for the computation of the mass spectrum of the theory, we observe that the t (0) * approximation progressively worsens when considering field values away from the radiative minimum.
The formalism of Ref. [66] can also be applied beyond the one-loop order, though the analytical expression of t * becomes rather involved. At two-loop order [80], it is possible to obtain an approximate relation for t * by retaining the leading terms in the form

Anomalous dimensions
If the matrix γ of anomalous dimensions does not vanish, the scalar fields are rescaled as so dΓ(t)/dt = −γ(t). In general, γ can be non-diagonal and may not commute with e Γ , but as a result of the Z n 2 symmetry of a biquadratic potential both γ and e Γ are diagonal and commute.
Our effective potential then takes the form where we used the fact that e Γ is diagonal to pull it through Φ.
If we ignore higher order terms involving two derivatives with respect to t, the formulae for the stationary point equation and the mass matrix, computed in the next section, will hold through the substitutions Λ → e 2Γ Λe 2Γ and β → e 2Γ βe 2Γ .

Minima of the effective potential
In this section, we first analyse the vacuum stability of the potential as it is a necessary condition to have a global minimum at finite field values. Then, we derive the stationary point equation, calculate the mass matrix and show how to solve the equation iteratively. Most of the discussion pertains to biquadratic potentials, but several results can be generalised to arbitrary cases by using tensor algebra.

Vacuum stability
A necessary condition for a finite global minimum to exist is that the scalar potential must be bounded from below, otherwise the theory is unstable and no state of lowest energy exists. We thus begin our investigation by requiring the stability of the scalar potential.
Because the second term in eq. (3.12) depends only on the direction of the field vector, t is a monotonic function of scalar fields in the limit of large field values [66]: The one-loop RG-improved effective potential is thus bounded from below if it obeys the tree-level stability conditions for values of the running quartic couplings given at large scales.
A biquadratic tree-level potential (2.2) is not directly a function of the field Φ, but rather of its Hadamard square Φ •2 , which is a vector of non-negative elements. Moreover, the potential is a homogeneous function of the fields. For these reasons, a necessary and sufficient condition for the biquadratic tree-level potential (2.2) to be bounded from below is that the matrix Λ be copositive [81]. A matrix A is copositive if it is non-negative on vectors with non-negative elements [82], that is x T Ax ≥ 0, ∀x ≥ 0.
Copositive matrices include the usual positive-semidefinite matrices, thus a positivesemidefinite Λ is sufficient, but not necessary, for the stability of the potential. Moreover, although Λ could be positive-semidefinite at a certain scale, this property is generally not preserved by the running of couplings. Indeed, even if the only contributions to the running of Λ arise from scalar interactions, we can regroup the terms in the β-function (3.8) as which is a sum of a matrix with non-negative elements and a positive-semidefinite matrix. Therefore, the running Λ remains copositive but may not be positive-semidefinite if its non-diagonal elements are large compared to the diagonal entries. The copositivity of a matrix can be ascertained, for example, via the Cottle-Habetler-Lemke (CHL) theorem [83]: Suppose that the order n − 1 principal submatrices of a real symmetric matrix A of order n are copositive. Then A is copositive if and only if det(A) ≥ 0 ∨ some element(s) of adj(A) < 0. It is now clear that the copositivity of the running coupling matrix Λ must hold in the limit of large field values for the effective potential to be bounded from below. On the contrary, we will show that a radiatively generated minimum requires, for its existence, a violation of the copositivity conditions in a finite range of field values. This criterion allows one to survey for the presence of an effective potential minimum solely on the basis of the running of quartic couplings. In particular, as we will see in section 4.4, the potential V < 0 in a minimum, implying det(Λ) < 0. Therefore, by the theorem (4.4), all elements of adj(Λ) must be non-negative in order for the minimum to exist.

Stationary point equation
We obtain the stationary point equation for the effective potential (3.10) by setting to zero its gradient with respect to the field vector Φ. Considering a biquadratic tree-level potential (2.2), the equation for the stationary point of the effective potential is given by where we have used the definition of the matrix β-function β = dΛ/dt, the chain rule and the relation (3.7). As anticipated, we neglect anomalous dimensions which can easily be taken into account with the substitutions Λ → e 2Γ Λe 2Γ and β → e 2Γ βe 2Γ in the results below. In the rest, Φ will mostly denote the field vector in the minimum, i.e. the vector of VEVs. The VEV of the radial direction ϕ = √ Φ T Φ is denoted by v ϕ . The field-space gradient of the t variable is given by where we neglect terms proportional to dA/dt and dB/dt consistently with the adopted t (0) * approximation. The two terms in the gradient (4.6) are orthogonal to each other. This is because A 0 and B 0 are both homogeneous functions of order four, so Multiplying the stationary point equation (4.5) from the left by Φ T and using (3.7), we obtain the radial stationary point equation which, in the case of a single scalar field with self-coupling λ and β-function β, reduces to the well-known relation 4λ + β = 0 [84][85][86]. Notice that the radial equation alone does not fully specify the minimum: it defines a (n − 1)-dimensional hypersurface in the field space, but the minimum direction on this surface is determined only by eq. (4.5).
In the next subsection we show that B > 0 in a minimum, hence for eq. (4.7) we must have V < 0 at the same point in field space. 3 Because the potential has to be bounded from below for large field values, it also follows that the potential vanishes at a radial distance beyond the minimum radius v ϕ . Such a point, together with the origin of the field space, identifies the (tree-level) flat direction associated with the minimum, be it of a Coleman-Weinberg or Gildener-Weinberg type [66].
The quantities A and B embodying the quantum corrections enjoy the same discrete symmetry as the biquadratic tree-level potential and, therefore, depend on scalar fields only via Φ •2 . We then apply the chain rule, and write the stationary point equation (4.5) as The first term in the resulting Hadamard product, Φ, gives the trivial vacuum solution: for a non-zero minimum, it can be used to force up to n − 1 VEVs to vanish.

Non-trivial minima
Let us assume, to begin with, that all of the VEVs are finite. Then in (4.9) the factor in parentheses must vanish. Using the radial equation (4.7) to trade B for V , we have Upon multiplying eq. (4.10) from the left by adj(Λ), we obtain (4.11) and observe that the potential value in the minimum is proportional to det(Λ). Using we finally find an implicit solution for the minimum, where the inequality must be satisfied in order for the solution to be physical (with the equality as a limiting case). In section 4.6 we solve eq. (4.13) iteratively, starting from the lowest order solution presented in section 4.5.
The general expressions for ∇ Φ •2 A 0 and ∇ Φ •2 B 0 entering the second term of (4.13) are (4.14) where we used eqs. (3.5) and (3.7). The gradient of the mass matrix is more precisely written as ∇ Φ •2 ⊗ m 2 , but we leave the tensor product implicit everywhere in order to avoid clutter. We denote with B S,F,V the contribution of scalars, fermions and vectors to B, respectively. The first term in second line of eq. (4.15) requires further action and, for purely scalar contributions, reduces to While ∇ Φ •2 B 0 can be written as the product of a matrix times the vector Φ •2 , the gradient ∇ Φ •2 A 0 cannot be written in this form. Therefore, it is impossible to solve the minimisation equation (4.13) analytically. Notice that in the presence of a dominant mass scale m 2 , eq. (3.11) gives and also

Mass matrix and minimum conditions
A stationary point of the potential is a minimum if the eigenvalues of the associated mass matrix (the Hessian matrix) are positive. The full quantum-corrected mass matrix is given by where the tree-level scalar mass matrix m 2 S is given in eq. (2.4), the gradient ∇ Φ t in eq. (4.6) and the gradient of B is where we have neglected second-order terms which contain two derivatives of t (such as dβ/dt). Finally, the double gradient of t is given by The explicit form of the last term of eq. (4.21) is involved, but its most important properties can be readily seen. We know that It can be shown, in addition, that the first two terms in eq. (4.21) are canceled by a part of the last term.
Sandwiching the mass matrix between the VEV vectors, we have Therefore B > 0 is a necessary condition for the extremum to be a minimum (at least as long as couplings are perturbative and dB/dt is negligible). Consequently, the radial equation (4.9) implies that the potential in the minimum respects V < 0. From Φ T m 2 S Φ = 12V we see that V < 0 is equivalent to having a negative eigenvalue in the tree-level mass matrix m 2 S , and therefore m 2 S ≺ 0. The potential can only be negative if the coupling matrix Λ is not copositive in the minimum. By the CHL theorem (4.4), we have that det(Λ) < 0 and all elements of adj(Λ) are non-negative. (All the submatrices of Λ are copositive, because we assumed that there are no minima in lower-dimensional subspaces of the field space.) Thus, the running coupling matrix Λ has a negative eigenvalue.
Hence, in a minimum where the VEVs of all fields are non-zero, we have that all hold and are equivalent. Besides minima, quantum corrections may generate saddle points indicated by V ≥ 0 (copositive Λ) and B ≤ 0. If there is at least one finite minimum, then the origin is a saddle point with V = B = 0. Therefore, if the potential has two local minima with a saddle point between them, then necessarily B = 0 at two points located between the saddle point and each minimum. A possible complication with the conditions will be discussed in the end of subsection 4.8.

Minimum solution to lowest order
If we approximate ∇ Φ •2 t by the first term in eq. (4.12), then we have The solution is given by where the unit vector in the direction of the minimum is When V = 0, the vector n •2 gives the tree-level flat direction. In the minimum, however, V < 0, so the lowest order solution Φ •2 is not a null eigenvector of Λ evaluated at the minimum. However, for practical purposes, the difference between the minimum direction and the flat direction is negligible in most cases. 4 Substituting the solution (4.25) into the radial minimisation equation (4.7) yields giving a condition on the running couplings that must hold in the minimum at the lowest order of approximation. For the value of the effective potential at the minimum, we find whereas the quantum-corrected mass matrix, in this approximation, is In fact, the last term in eq. (4.29) can be dropped, because it will be canceled by a part of the ∇ Φ ∇ T Φ (A 0 /B 0 ) term (see discussion below eq. (4.21)). We can now obtain the usual Gildener-Weinberg potential [56] in the radial direction by making further approximations [87]. First, we insert the lowest order minimum solution (4.25) into the potential (3.9). Then, by approximating the running to the linear order 4 Because the matrix adj(Λ) has only positive elements, it is guaranteed to have a Perron-Frobenius eigenvector v with all elements positive, corresponding to the largest eigenvalue r of the matrix. It is interesting that the vector (4.25), normalised, is usually very close to the Perron-Frobenius eigenvector. The eigenvalues of adj(Λ) are products of the eigenvalues of Λ with one eigenvalue absent in each product. At the flat direction, Λ has a zero eigenvalue, so adj(Λ) is a projection matrix with only one non-zero positive eigenvalue and the vector (4.25) is the Perron-Frobenius eigenvector. In the minimum, the other eigenvalues of adj(Λ) are non-zero and negative, but small compared to the Perron-Frobenius eigenvalue, so adj(Λ) ≈ rvv T and thus n •2 ∝ v.
in t, we have where we have used the radial equation (4.7) in the last step.

Iterative solution
In order to go beyond the Gildener-Weinberg approximation, we solve the minimum equation (4.13) iteratively. After selecting an RG flow by giving initial conditions for the couplings at an arbitrary scale µ in , we start the iteration procedure from eq. (4.25), which gives the minimum direction at zeroth order. The zeroth order coupling matrix Λ 0 = Λ(t 0 ) is obtained by inserting n •2 0 in the radial equation (4.7) -resulting in eq. (4.27) -and finding the scale µ 0 = µ in e t 0 at which the radial equation is satisfied. In order to improve the accuracy of the t (0) * approximation, we now use the RGE solutions to reparametrise the running of couplings, taking as equivalent initial conditions their values at µ 0 . As a consequence, after the reparametrisation, eq. (3.14) reads t 0 → t(ϕ = v ϕ0 , µ = µ 0 ) = 0 and we can solve it to obtain the norm v ϕ0 of the zeroth order solution The first order solution is then obtained by inserting Φ •2 0 into the RHS of eq. (4.13): where again Λ 1 = Λ(t 1 ) satisfies the radial equation (4.7) at the scale µ 1 = µ 0 e t 1 . Proceeding as before, we obtain the first iteration v ϕ1 for the norm after reparametrising the RGE solutions so that t 1 → t(ϕ = v ϕ1 , µ = µ 1 ) = 0. The procedure can be repeated until sufficient precision is attained.
Notice that, at first order, the second term in Φ •2 1 is orthogonal to e: The method can also implement the solution to the spurious negative mass of the Goldstone boson proposed in [73,74], by shifting -at each iteration step -the fielddependent scalar mass matrix m 2 S by ∆ Φ G Φ T G at that step, where ∆ is the correction to the Goldstone mass and Φ G is the Goldstone direction at that order.
A criterion for the convergence of the iteration is that the spectral radius of the Jacobian -the largest of absolute values of its eigenvalues -of the RHS of eq. (4.13) at the fixed point be less than unity. To assess convergence it is necessary to evaluate the Jacobian at least to the first order: the potentially large term ∇ Φ •2 (A 0 /B 0 ), which could -in principle -spoil the convergence, does not enter at zeroth order. The Jacobian is calculated in detail in the appendix E. In the example analysed below, we find that the iteration converges to a fixed point for a large range of couplings.

Inverse problem
In the previous sections, after specifying the quartic, gauge and Yukawa couplings at a given scale, we have solved the minimisation equations and found the minimum of the potential together with the mass matrix. We consider now the opposite problem: given the full mass matrix M 2 and the desired minimum Φ •2 , is it possible to find a corresponding quartic coupling matrix Λ and the matrix β-function at that point?
In the mass matrix (4. 19), Λ appears explicitly in the tree-level part (2.4): In the first term, we can use the minimum condition (4.5) to substitute 4ΛΦ We can now approximate t to the lowest order, t ≈ (1/2) ln(ϕ 2 /µ 2 ), so that ∇ Φ t and ∇ Φ ∇ T Φ t depend on Φ only. At this stage, with 8B = Φ T M 2 S Φ, the only unknown is the gradient 2∇ Φ B = 4Φ • βΦ •2 . Although in principle the β-functions are all given in terms of Λ, gauge and Yukawa couplings, still the dependence is non-linear.
We solve the issue by considering the lowest order solution (4.29) for M 2 S , where we neglect the dependence on the unknown ∇ Φ B and drop the last term, thereby obtaining that The above formula can be used as a starting point for a "reverse" iteration that leads to the desired Λ. At higher orders, we can use the expression (3.14) for t and also obtain the

Minimum in a field subspace
If up to n − 1 components of the field vector Φ vanish in the minimum, we can bring the coupling matrix, the β-function matrix β and the field vector and ∇ Φ t into a block form upon a suitable permutation of the fields: Consistency with the null VEVs forces the elements of ∇ Φ t corresponding to zero elements of Φ to vanish too, otherwise these VEVs could be pushed to non-zero values. That the elements do vanish follows from the fact that the effective potential enjoys the same symmetries as the tree-level potential. Therefore, we must solve the stationary point equations (4.10) and (4.7) restricted to the Φ 1 subspace: where the running of Λ 11 still depends directly on Λ 12 and, indirectly, on Λ 22 . The full mass matrix also admits a block form where Φ 2 is the subspace of fields with null VEVs. In fact, only the fields in the Φ 1 subspace, whose VEVs break the corresponding Z 2 symmetries, are allowed to mix with each other and M 2 22 then necessarily remains a diagonal matrix. A positive-definite M 2 22 constrains the off-diagonal Λ 12 , which also enters eqs. (4.38) and (4.39) via β 11 . The only direct constraint on Λ 22 is that it be copositive.
As discussed in detail in section 6, the effective potential can have more than one minimum. If two minima are located in the same field subspace (e.g. on the same axis), or in two orthogonal subspaces (e.g. two different axes), the corresponding existence conditions are independent of each other. It is also possible that two stationary points fulfil minimum conditions in part (e.g. have V < 0), while only one of them can be an actual minimum. Such a possibility arises when a stationary point lies in a subspace on the border of the subspace of the other point, for instance with one stationary point on the φ 1 -axis and another on the φ 1 φ 2 -plane. In such cases, it is possible to assess which stationary point is an actual minimum by checking if the associated mass matrix is positive-definite.

Generalisation to generic potentials
In general, a scalar potential can have the form (2.1) and, unlike a biquadratic potential (2.2), contain terms such as φ 3 1 φ 2 . While the minimum solution for the biquadratic case has a specific form, other results carry over with minimal changes. In tensor notation, we can write the RG-improved potential as where Λ is the tensor of quartic couplings. A polynomial f (x) of order m and its coefficient tensor A are related as Ax m = mf (x) and The minimisation equations are given by where β is the β-function tensor. The radial equation is Consequently, the relation Φ T M 2 S Φ = 8B still holds and B > 0 is a necessary condition for the minimum to exist. Therefore V < 0 in the minimum, the Λ tensor must not be positive-definite and one of the tensor eigenvalues λ of Λ must be negative.
The E-eigenvalue equations for the Λ tensor [88,89] are given by but tensors also admit N-eigenvectors and N-eigenvalues, given by the solutions of Ax m−1 = λx •(m−1) . We refer the reader to ref. [90] for a thorough discussion. The number of E-eigenvectors of a symmetric tensor of order m in R n is They can be complex, but only real tensor eigenvectors and eigenvalues are physical solutions of eqs. (4.51) and (4.52). Eliminating φ i from eqs. (4.51) and (4.52) yields the characteristic polynomial φ Λ (λ) of the tensor, its degree is given by eq. (4.53). The multivariate resultant of a system of polynomial equations is a polynomial in their coefficients, which vanishes if and only if the equations have a common root. The free term of φ Λ (λ) -the product of all E-eigenvalues -is given by the resultant res Φ (ΛΦ 3 ). For that reason, the resultant is also called the hyperdeterminant. Having a negative tensor eigenvalue then requires res Φ (ΛΦ 3 ) < 0.
In the limit of large field values, the tensor Λ must be positive-definite in order for the potential to be bounded from below: all of its eigenvalues and those of its principal subtensors -obtained by setting one or more fields to zero -must be positive [91,92].
Note that for a quartic potential of two fields, the resultant is proportional to its discriminant computed with one field set to unity. Unfortunately, for a larger number of variables, the calculation of the resultant is computationally expensive, but the tensor eigenvalues can still be readily calculated numerically. Unlike for biquadratic potentials, in general all the potential coefficients cannot be determined by specifying the desired mass matrix.

Examples
In general, it is not possible to find analytically the direction in field space along which the minimum of the effective potential lies. Two notable exceptions are given by the case where a minimum is on an axis of field space and by the "democratic" setup, where the minimum lies equally distant from all axes. These configurations are, in a sense, complementary: the minimum on an axis is the extreme case of the minimum in a subspace (section 4.8), where all but one fields are set to zero. In the democratic case, instead, all fields are equally nonzero. Despite this evident difference, the scenarios are united by the fact that the second term in ∇ Φ t identically vanishes and, therefore, the lowest order minimum solution (4.25) is exact. To go beyond these simple examples, we will resort to the iterative solution of section 4.6 for more realistic cases.

Minimum on an axis
We start with the simplest case, corresponding to a minimum of the effective potential that lies along a field space coordinate axis. This is also the simplest example of the possibilities analysed in section 4.8, where not all the fields develop a finite VEV. Let the solution with a non-zero v ϕ be in the direction Φ = v ϕ e i , where e i is the unit vector along the i-th axis. Note that e •2 i = e i . From the radial minimisation equation (4.7), we obtain 4e T i Λe i + e T i βe i = 0 (5.9) or 4λ ii + β ii = 0, (5.10) which -if we assume only scalar contributions -takes the form The minimisation equation (5.11) is a quadratic equation in λ ii . Of the two solutions, one is dominated by λ ij and yields a perturbative value for λ ii . We will see that it is a minimum of the effective potential. 5 The tree-level mass matrix (2.4) is given by where Λ i is the ith column of the coupling matrix Λ. Since it is a diagonal matrix, the tree-level masses are simply m 2 i = 12λ ii v 2 ϕ < 0 and m 2 j =i = 2λ ij v 2 ϕ . In order to calculate the quantum-corrected mass matrix, we use 14) and through eq. (4.19) find that the mass matrix is a diagonal matrix with entries More concretely, suppose that the minimum of the potential (5.1) is on the φ 1 -axis. Then λ 1 < 0 at the minimum and β 1 > 0. The block form of the coupling matrix in eq. (4.37) is given by Λ 11 = (λ 1 ), Λ 12 = (λ 12 ), Λ 22 = (λ 2 ) and Φ 1 = (v ϕ ). The copositivity of Λ 22 thus requires λ 2 > 0. We have, neglecting λ 1 and λ 2 , that For perturbative couplings, the requirement M 2 2 > 0 implies λ 12 > 0. An example of potential (5.1) with a minimum on the φ 1 -axis is shown in figure 1, where we used λ 1 = −7.916 × 10 −3 , λ 2 = 0.1, λ 12 = 0.5 and the couplings have been evaluated at the radiative minimum for µ = exp(−3/4)v ϕ . The field values are given in units of v ϕ . The minimum (indicated by a red dot) lies in a basin of negative potential values. We have M 1 = 0.113 v ϕ and M 2 = 1.006 v ϕ . The second term in (5.18) amounts to a 0.6% correction to M 2 .

Democratic minimum
The other case that can be analytically solved results from a "democratic" choice of the couplings: we set all self-couplings to a common value λ ii = λ s and all portal couplings equal to each other λ ij = λ p , i = j. Hence, where ee T − I is a matrix with zeros on the main diagonal and ones everywhere else. Provided that possible gauge and Yukawa couplings are also the same for all scalars, the symmetry of the problem determines the following solution: Obviously, given the assumptions, the β-functions β s for the self-couplings are equal, and so are those for the portals, β p . Eq. (5.21) thus takes the form The tree-level mass matrix (2.4) is given by where we used diag(e) = I. The eigenvalues of the matrix (5.23) are given by the eigenvalues of the last term, shifted by the coefficient of the unit matrix. Since (ee T )e = ne, the eigenvalues of ee T are n, with multiplicity 1, and 0, with multiplicity n − 1. That is, the tree-level masses along the minimum direction and in the orthogonal directions are given, respectively, by The minimum field vector (5.20) is an eigenvector of the tree-level mass matrix with m 2 ϕ v 2 ϕ = Φ T m 2 S Φ and an eigenvector of the quantum-corrected mass matrix. Inserting the eigenvalue m 2 ⊥ in the eigenvector equation of m 2 S , we find that e T n ⊥ = 0, as required. The vector n ⊥ has unit norm and is orthogonal to the VEV Φ in field space.
The quantum-corrected masses at the minimum are where we used M 2 ϕ = 8B/v 2 ϕ = −16V /v 2 ϕ and the fact that both ΛΦ •2 and βΦ •2 are proportional to e and thus orthogonal to n ⊥ . For the two-field potential (5.1), we have in the democratic case λ 1 = λ 2 and the minimum solution is An example of democratic minimum for the potential (5.1) is shown in figure 2, where we took λ 1 = λ 2 = 0.244 and λ 12 = −0.5 in the minimum with µ = exp(−3/4) √ 2v ϕ . The masses are M ϕ = 0.23 v ϕ and M ⊥ = 1.43 v ϕ . The second term in M ⊥ amounts to a 1.5% correction.

Iterative solution
To exemplify the iterative solution to the stationary point equation, we choose λ 1 = 0.06, λ 2 = 1.125 and λ 12 = −0.5 at µ in = 1 in the potential (5.1). The chosen values of the self-couplings differ by an order of magnitude to strongly depart from the democratic case discussed above.
First, we solve the RGEs to obtain the running couplings. All the elements of adj(Λ) are positive in the neighbourhood of µ in = 1, indicating that the potential may have a minimum lying on the field plane. 6 We establish that the effective potential does have a minimum by studying the running of det(Λ): there is a flat direction at µ = 0.789, below which the determinant is negative.
To begin the iteration procedure, we insert the lowest order solution (4.25) in the radial equation (4.7), which we solve for t 0 and Λ(t 0 ) simultaneously. Then we reparametrise our RGE solutions and solve t(ϕ = v ϕ0 , µ = µ 0 ) = 0 to determine v ϕ0 . The first order iteration is obtained by inserting Φ  Table 1. Iterative solution to the minimum. The first iteration is already very close to the true minimum; to the given accuracy, the third iteration gives the same result as the second.
Λ(t 1 )), and finding v ϕ1 from t 1 after a further reparametrisation. The iteration procedure is repeated until sufficient precision is obtained.
The results after each iteration step are given in table 1, where we present the renormalisation scale µ at the minimum and the corresponding values of the couplings. In place of the field components φ 1 and φ 2 , we give the radial coordinate v ϕ and the minimum vector angle cos θ. In addition, we display the masses and the mixing angle cos α, which is close, but not exactly equal, to cos θ. The first order correction is already within sight of the actual minimum and the second and third order give identical results to the given accuracy. The mass M 1 in the lowest order approximation is overestimated by 4.5%, whereas the mass M 2 is 11% smaller than the true value at the minimum. Table 1 also shows the values of the spectral radius ρ(J ) of the Jacobian given in the appendix E. When it is less than unity at the fixed point, the iteration converges. Figure 3 shows the iteration steps to the minimum overlaid on the contour plot of the potential V . The lowest order solution is indicated by a yellow dot, the first iteration by an orange dot and the second iteration by a red dot. The third iteration already coincides with the second one for all practical purposes. The angle between the fixed point minimum vector and the lowest order solution (4.25) differs by less than a degree, but the length v ϕ0 of the lowest order solution differs from the actual length by 2.8%.
To analytically gauge the size of the corrections, we compare the two terms in the solution (4.13) for Φ •2 at first order, neglecting the Goldstone mass for simplicity. The first term is (5.28) and the second term is Note that ∆Φ •2 is only approximately equal to Φ •2 1 −Φ •2 0 , as the values of quartic couplings used in Φ •2 0 and Φ •2 1 are different. Because ∆Φ •2 is orthogonal to e, it is convenient to compare the projections Φ •2 0 and ∆Φ •2 on e ⊥ = (1, −1) T . Furthermore, we expect the iteration to converge when the ratio of norms is less than unity. For the example at hand, e T ⊥ ∆Φ •2 /e T ⊥ Φ •2 0 is 3.4% and ∆Φ •2 / Φ •2 0 is 1.7% at the minimum. Note that as ∆Φ •2 ∝ |λ 1 − λ 2 |, the correction is zero when the self-couplings are equal, which is exactly the case of the democratic minimum (section 5.2). Now that we have found the minimum, let us test our method for the inverse problem of section 4.7 and find the coupling matrix Λ, starting with the minimum vector Φ and the mass matrix M 2 S . The lowest order approximation (4.36) to the coupling matrix gives λ 1 = 0.0579, λ 12 = −0.4904 and λ 2 = 1.152, so the relative error with the true coupling matrix is Λ − Λ 0 / Λ = 0.25. The renormalisation scale in the minimum is obtained from solving the equation t(ϕ = v ϕ = √ Φ T Φ, µ) = 0 for µ = 0.728v ϕ . Beginning with these values, we calculate B and t and use eq. (4.35) to obtain Λ via iteration. It takes 10 iterations to bring the relative error Λ i+1 −Λ i / Λ i+1 below 1%, after which the relative error with the true coupling matrix is Λ − Λ 10 / Λ is about 3% (further iteration does not improve the precision).

Multiple minima
The lion's share of literature considers only the presence of one minimum in the effective potential. Although it has been pointed out that more minima could co-exist [67], to our best knowledge, a detailed discussion of the implications of such a possibility is still absent.
The simplest way to construct an effective potential with several separate minima is to arrange them in different same-dimensional subspaces of the field space, e.g. on field space axes (along the lines of the discussion in section 4.8). For example, we discussed a minimum on a field space axis in section 5.1 within the two-field model with the potential (5.1). In this scenario, in order to have a minimum on each axis, it is clear that both self-couplings λ 1 and λ 2 must be negative below certain scales, resulting in separate minima that need not be degenerate (v ϕ1 = v ϕ2 ). Such a situation is depicted in figure 4 for the values where v ϕ1 is the VEV of the φ 1 field. In this case, it is the positive portal coupling λ 12 that separately drives both the self-couplings through zero. The value of the potential in the φ 1 minimum is V = −7.92 × 10 −4 v 4 ϕ1 and the particle masses are M 1 = 0.113v ϕ1 and M 2 = 1.006v ϕ1 . In the φ 2 minimum with v ϕ2 = 1.165 v ϕ1 (at µ = 0.550v ϕ1 ), the value of the potential is V = −8.04 × 10 −4 v 4 ϕ2 = −1.46 × 10 −3 v 4 ϕ1 and the particle masses are M 1 = 0.132v ϕ2 = 0.153v ϕ1 and M 2 = 0.870v ϕ2 = 1.01v ϕ1 . In a similar manner, it is possible to arrange two or more minima situated on different field planes and so on. This way for constructing effective potentials with multiple minima is robust and can effortlessly be implemented in models with more than two fields.
In more involved cases, several minima lie in the same subspace of field space (an axis or a plane, for instance). In order to have two minima on the same field space axis, the self-coupling must increase, decrease and increase again. Its β-function must therefore first become negative and then positive again, crossing zero twice. Similarly, if two or more minima are on a field plane or in a larger subspace, B has to run through zero between every pair of stationary points, as shown in section 4.4.

Two minima on an axis
We start by discussing how the toy model effective potential can admit two minima on the same axis of the scalar field space.
In order to generate two minima on the φ 1 -axis, the running self-coupling λ 1 must be negative in two different intervals. Let us set λ 1 negative at low scales (t = 0). Its β-function β λ 1 must then be positive at low scales, negative at intermediate scales, and positive at high scales to ensure vacuum stability. In short, it must run through zero twice. For negligible λ 1 contributions, we have that The positive contributions from scalars and gauge bosons must dominate at low scales and high scales, while at intermediate scales, instead, the negative Yukawa contributions must prevail. In order to keep the minimum on the axis we need a positive portal coupling λ 12 .
As the portal coupling grows with the running, it will dominate the large scale region. The non-Abelian gauge coupling g, on the contrary, keeps β λ 1 positive at low scales but  Figure 5. Running of couplings that generates two minima on a field space coordinate axis. Left panel: Running of the scalar self-couplings λ 1 and λ 2 , the portal coupling λ 12 , the gauge coupling g and the Yukawa coupling y. Middle panel: The self-coupling λ 1 is negative in two regions and the effective potential has a minimum in each. Right panel: Such running of λ 1 is due to a β λ1 that is positive at low and high scales and negative at intermediate scales.
decreases with the scale. An almost constant Yukawa coupling then ensures that β λ 1 is negative at intermediate scales.
Since we already know the minimum direction, we find the scales of both minima from the radial eq. (4.7). We choose the values of the couplings in the first minimum with unit VEV v ϕI (at the renormalisation scale µ = 0.264v ϕI ) as λ 1 = −1.962 × 10 −4 , λ 2 = 0.01832, λ 12 = 0.4003, g = 0.7039 and y = 0.6287. The running of couplings is shown as a function of t = ln(µ/0.01v ϕI ) in the left panel of figure 5. The middle panel shows the running of the self-coupling λ 1 and the right panel the running of its β-function β λ 1 .
The value of the potential in the first minimum is V = −4.906 × 10 −5 v 4 ϕI . The mass matrix is diagonal. The mass of the pseudo-Goldstone of the scale invariance is M 1 = 0.0280v ϕI . Since the global U (1) associated by the singlet Φ 2 is not broken, both its components have the same mass M 2 = 0.447v ϕI . The field is not rescaled in the first minimum, because there we choose our initial conditions, but in the second minimum we have to take into account anomalous dimensions. We find that the renormalisation scale µ = 7.098 × 10 18 v ϕI in the second minimum, so φ 1 is scaled by exp(Γ 1 (t)) = 1.0541, while φ 2 does not scale. The VEV in the second minimum is therefore v ϕII = exp(−Γ 1 (t))1.282 × 10 19 v ϕI = 1.216 × 10 19 v ϕI , the value of the potential is V = −7.170 × 10 −5 v 4 ϕII = −1.582 × 10 72 v 4 ϕI and the particle masses are M 1 = 0.0322v ϕII = 3.925 × 10 17 v ϕI and M 2 = 0.6611v ϕII = 8.057 × 10 18 v ϕI . As expected, a saddle point separates the two minima.

Two minima in a plane
Another scenario we can explore with our toy model has two minima in the field plane, i.e. two minima, where both scalars have non-zero VEVs. We can use a similar configuration of couplings as in the previous example. The exception is that now the portal coupling λ 12 must be negative. However, a negative λ 12 tends to run slower, so we need a larger λ 2 to make its absolute value large at high scales. For that reason, the second minimum is generated right before the couplings hit a Landau pole. Whereas this issue should be addressed in realistic models, it is irrelevant for the purpose of our example.
In order to forbid minima on the axes, we require λ 1 > 0 and λ 2 > 0. To check whether two minima appear, we track the running of det(Λ). This quantity must run thrice through zero, similarly to the λ 1 in the case of two minima on the φ 1 -axis treated in the previous subsection. We also need adj(Λ) > 0, which holds on the whole range considered, since λ 1 > 0, λ 2 > 0 and −λ 12 > 0 for all t. Figure 6 shows the running of couplings, of det(Λ) and of β det(Λ) as a function of t = ln(µ/0.01v ϕI ).
We find the first minimum via the iteration procedure in section 4.6. We give everything in units of its radial VEV v ϕI . The values of couplings in the first minimum are λ 1 = 0.1195, λ 2 = 0.1718, λ 12 = −0.2872, g = 0.7208 and y = 0.7169 at µ = 0.231v ϕI . The first minimum is at φ 1 = 0.738v ϕI , φ 2 = 0.674v ϕI . The value of the potential is V = −3.873 × 10 −5 v 4 ϕI . The masses of the singlet-doublet mixed mass eigenstates are M 1 = 0.0248v ϕI and M 2 = 0.530v ϕI . The imaginary part of Φ 2 is the Goldstone of the breaking of a global U (1).
The iteration does not converge for the second minimum due to large quantum corrections in the vicinity of the Landau pole. Therefore, we use the full expression for t * to improve the effective potential and find the minimum numerically. The second minimum is at µ = 3.27 × 10 16 v ϕI , which scales the φ 1 component by exp(Γ 1 (t)) = 1.0267, while φ 2 does not scale. The second minimum is then at φ 1 = 0.919v ϕII , φ 2 = 0.333v ϕII with v ϕII = 3.17 × 10 16 v ϕI . The value of the potential at the second minimum is V = −5.943 × 10 −5 v 4 ϕII = −6.001 × 10 61 v 4 ϕI and the particle masses are M 1 = 0.809v ϕII = 2.56 × 10 16 v ϕI and M 2 = 0.0620v ϕII = 1.96 × 10 15 v ϕI . Again, a saddle point separates the two minima.

One minimum on an axis, one on the plane
Is it also possible to arrange one minimum on an axis and another on the field plane? With only two scalars it is impossible, though such a configuration may be achievable considering three or more fields. With two scalars, in fact, the running portal coupling λ 12 keeps its sign. To have a minimum on the plane, we need a negative portal coupling, but this is in contradiction with the λ 12 > 0 requirement imposed to keep the orthogonal field mass (5.18) positive for the minimum on an axis. Even in the case of a nearly vanishing portal, other couplings cannot produce a positive mass via quantum corrections (as in the second term of eq. (5.17)) without non-perturbative values.

Conclusions
We have considered the minimisation of effective potentials with multiple scalar fields by way of a novel matrix formalism in which calculations are presented in a compact and intuitive manner. Our approach is easy to implement in any modern computer algebra system. The matrix formalism could also be useful for studying the symmetry breaking due to negative mass-squared terms.
To demonstrate the formalism in a clear way, we consider an RG-improving method such that the one-loop corrections vanish due to a specific field-dependent choice of the renormalisation scale [66]. The resulting effective potential has the tree-level form with running couplings evaluated at that scale. In the main, we are concerned with popular biquadratic potentials, but many of our results apply to more general potentials via tensor algebra. The approach is not limited to the particular improvement procedure adopted and can be adapted to other choices of the renormalisation scale, such as µ = ϕ, through the straightforward changes detailed in appendix F. Dimensionful terms arising from thermal corrections or from a non-minimal coupling to gravity in a curved background can also be treated in our matrix formalism because they respect the biquadratic symmetry of the scalar potential. In the latter case, for example, the renormalisation scale can depend on the curved background [94]. 8 Our first main result concerns the stationary point equation (4.13) for the scalar potential, derived within a matrix formulation. At the lowest order (4.25), the minimum solution has the same functional form as the corresponding flat direction [68], but differs from it due to quantum corrections. Neglecting those, the lowest order solution reproduces the Gildener-Weinberg approximation. In general, the minimum equation can only be solved iteratively. We show that the iteration converges and discuss in detail the first order correction with a semi-analytical treatment that goes beyond purely numerical considerations. The lowest order solution is exact in two special cases, which we solve analytically: a minimum on a field space axis and a "democratic" minimum, where all scalar fields have equal self-couplings and identical portal couplings.
We stress that the formalism is not limited to Gildener-Weinberg minima. In fact, given a generic effective potential with a radiatively generated minimum, since at large field values we must have V > 0 and in the minimum necessarily it is V < 0, there is always a point in field space with V = 0 at some intermediate value of the radial coordinate. For Gildener-Weinberg solutions, this flat direction and the minimum are aligned. Differently, when the tree-level and one-loop contributions are comparable (as typical of Coleman-Weinberg solutions), the "flat direction" and the minimum may not be aligned at all. As our criterion for the existence of a minimum given in (4.23) depends only on the quartic couplings, our formalism can address both mentioned cases.
Once the direction in field space of the effective potential minimum is determined, one can calculate the quantum-corrected mass matrix (4.19). Its positive-definiteness requires the potential and the quartic coupling matrix to satisfy a set of conditions in the minimum (4.23). We also solve the inverse problem, obtaining expressions for the quartic couplings and the β-functions at the minimum given the desired minimum field vector of VEVs, physical masses and mixing angles for the fields.
The second main result concerns the exploration of effective potentials with several physical minima. In the literature, only a single Gildener-Weinberg type of minimum is generally considered, although more minima can be present. For example, a two-scalar biquadratic effective potential can admit two minima located on different axes of the field space (figure 4), two minima on the same axis (figure 5), or two minima in a common field plane at different distances from the origin (figure 6). We provide semi-analytical criteria to determine whether an effective potential has more than one minimum. In such cases, the running quartic couplings must satisfy the same set of conditions (4.23) at two different scales or in orthogonal field subspaces. Examples are given in a minimal toy model with two complex scalar fields, two Weyl fermions and an SU (2) gauge symmetry meant to resemble the SM.
Because most of the models in the literature are concerned with biquadratic potentials, we expect that our methods will facilitate the exploration of such effective potentials. Considering classically scale-invariant effective potentials with more than one minimum may open new directions in phenomenological studies of phase transitions and cosmic gravitational wave signals.

A β-function matrix
We derive the β-function for the quartic coupling matrix Λ of the biquadratic potential (2.2). In the absence of gauge or Yukawa interactions, the anomalous dimensions at one-loop level are zero, and eq. (3.6) implies that

B Biquadratic symmetry-preserving Yukawa couplings
The Yukawa Lagrangian for left-handed complex Weyl fermions ψ i is given by where = iσ 2 is the fermion metric. In general, Yukawa couplings do not respect the Z n 2 symmetry of a biquadratic potential.
In order for a Yukawa coupling involving the field φ a to respect the symmetry, the fermion bilinear ψ j ψ k must be odd under (Z 2 ) a . Hence, one of the fermions must be odd under (Z 2 ) a and both fermions may also be odd under a subset of the other Z 2 symmetries.
An even parity corresponds to the Z 2 charge 0, since +1 = e i0π ; similarly an odd parity corresponds to the Z 2 charge 1, since −1 = e i1π . It is most convenient to start counting the Yukawa matrix indices from zero and number a scalar or fermion by treating its representation as a binary number 9 -with the order of bits reversed, so that the numbering of (Z 2 ) a begins from the left: 0 = (0, . . . , 0), 1 = (1, 0 . . . , 0), 2 = (0, 1, 0, . . . , 0), 3 = (1, 1, 0, . . . , 0), . . . , 2 n − 1 = (1, . . . , 1). There can be several fermions in the same representation -for instance if the fermion carries colour -whose Yukawa couplings we sum together for simplicity. The n scalars are naturally numbered from 2 0 to 2 n−1 , where only one bit equals 1 and others are 0, so the scalar φ a has number 2 a . A Yukawa coupling between φ a , ψ i and ψ j can thus be non-zero only if 2 a ⊕ i ⊕ j = 0, where ⊕ is the bitwise 9 We note that a mixed radix number system could be used to number representations of more complicated discrete symmetries whose factors are not only Z2, but also higher ZN symmetries.
The field-dependent fermion mass matrix is given by Then we have, for example, The bitwise xor is commutative and associative and satisfies Eqs. (B.6) therefore imply and, in turn, that 2 a ⊕ 2 b = 2 c ⊕ 2 d , 2 a ⊕ 2 d = 2 b ⊕ 2 c and 2 a ⊕ 2 c = 2 b ⊕ 2 d , that is a ⊕ b = c ⊕ d, a ⊕ d = b ⊕ c and a ⊕ c = b ⊕ d. We see that assuming e.g. a = b implies c = d, and so on. On the other hand, if we set e.g. a = b = c = d, the equations (B.6) have no solution. Let us show that there are no non-biquadratic terms in tr m 4 F . Having a non-biquadratic term requires that at least two φ i differ, e.g. a = b without loss of generality. Then if either a or b equals c or d, we saw that the term must be biquadratic. Therefore all a, b, c and d must be different from each other to have a chance to produce a non-biquadratic term. For that, we need n ≥ 4.

C Gauge contributions to the β-function
The scalar fields are coupled to gauge bosons via their covariant derivatives, given by where the representations θ A of the group generators are Hermitian matrices. Because we have decomposed complex scalars in terms of their real components, the θ A matrices are purely imaginary and antisymmetric. The field-dependent vector boson mass matrix is where κ is 1/2 for Weyl and 1 for Dirac fermions [79] and where C 2 (S) a is the quadratic Casimir operator C 2 (S). In the presence of several gauge groups g 2 C 2 (S) → k g 2 k C k 2 (S). We see from eq. (D.3) that the gauge contribution is always diagonal, while in general the Yukawa contribution is not. We will now specialise to Yukawa couplings that obey the Z n 2 symmetry of the potential, as detailed in appendix B. We have, for Yukawa couplings of Weyl fermions given by eq. (B.1), that where we used δ (2 a ⊕i⊕j)0 δ (2 b ⊕i⊕j)0 = δ ab δ (2 a ⊕i⊕j)0 . Hence, the matrix of anomalous dimensions for a model with Z n 2 symmetry is diagonal and given by 16π 2 γ ab = i,j y a ij 2 δ (2 a ⊕i⊕j)0 − 3g 2 C 2 (S) δ ab . (D.5)

E Convergence of the iterative solution
A criterion for the convergence of the iteration is that the spectral radius of the Jacobian of the RHS of eq. (4.13) at the fixed point be less than unity. The spectral radius of a matrix is given by the largest absolute value of its eigenvalues. The Jacobian is where we have used the minimisation equation 0 = ∇ Φ •2 V and eq. (A.6). In addition, the radial minimum condition (4.7) yields and we have that Thus at zeroth order, for example, we have where V (0) (t) is the tree-level potential with field-dependent effective couplings, and the running parameter is given by The stationary point equation is Using the chain rule (4.8), we can write Eq. (F.6) as We solve eq. (F.8) for Φ = 0. An implicit equation for the solution is given by where we used the radial eq. (F.7).
To lowest order, we neglect the last term in eq. (F.9) and take V ≈ V (0) , yielding where we used eq. (4.28) with e → e M . For the first iteration, we insert the lowest-order solution (F.10) in (F.9).
The mass matrix around the minimum is then given by where m 2 S is the tree-level scalar mass matrix. Notice that the term proportional to B is canceled by a similar term resulting from ∇ Φ ∇ T Φ A. For a minimum in a field subspace, like in section 4.8, the mass matrix still obtains a block-diagonal form due to the discrete symmetry of the potential.