Symmetry breaking patterns for two coupled complex scalar fields at finite temperature and in an external magnetic field

A model of two coupled complex scalar fields is studied at finite temperature and under an external magnetic field. The results are obtained in the context of the nonperturbative method of the optimized perturbation theory and contrasted with those obtained in perturbation theory and in the one-loop approximation. The emergence of phenomena related to inverse symmetry breaking and symmetry nonrestoration are analyzed.


I. INTRODUCTION
Phase transitions are ubiquitous in nature. In physics, phase transition phenomena range from those possibly happening at high energies in the early Universe, e.g., at the electroweak [1] and quantum chromodynamics (QCD) scales [2,3], to the ones in lower-energy systems, like in condensed matter [4] and in atomic gas systems [5]. The first studies along this topic give a good qualitative and quantitative idea of how these phenomena might occur. The general expectation is that a system initially in a symmetry broken state will have its symmetry restored at sufficiently high temperatures, beyond some critical temperature. An unambiguous example is the rotational symmetry of a ferromagnetic system, which is broken at temperatures below a critical value, the Curie critical temperature, and then gets restored as the temperature is raised beyond that critical value.
The above picture gets less clear when more than one symmetry is at play and multiple critical points can be present. This is a common situation in condensedmatter systems, where reentrant (intermediate) phases are present [6]. In this case, the intermediate phase might, e.g., be a state which is less symmetrical than the previous one and later phases. This is a case we can characterize as being of an inverse symmetry breaking (ISB). The opposite case is also possible, where the intermediate phase would display a more symmetrical state and then, later, at higher temperatures, be displaying less symmetry. This would be an example of what we can call symmetry nonrestoration (SNR). Both ISB and SNR are also possible states that can emerge in high-energy physics. Systems exhibiting such forms of phenomena can be constructed in terms of multiple scalar fields with both self-couplings and intercouplings and were first investigated by Weinberg in Ref. [7] in the context of a two coupled scalar field system at finite temperature in quantum field theory.
In this work we are motivated in investigating the interplay of both temperature and external fields, more specifically, an external magnetic field, in the phase transition patterns for a system of coupled complex scalar fields, symmetric under a global U (1) × U (1) symmetry. With the ISB and SNR phenomena gaining recently more interest in different contexts, the investigation of not only temperature as a cause of phase transition in these models, but also when the system becomes subject to other external conditions, becomes pertinent. Phase transitions in the presence of external magnetic fields are relevant for understanding many important physical systems, ranging from superconductor materials in the context of condensed matter models [27], in astrophysical systems [28], in QCD and heavy-ion collision problems [29], and in the early Universe [30]. In particular, it is well known that magnetic fields can lead to phenomena enhancing the symmetry breaking, e.g., the magnetic catalysis [31], or leading to an opposite effect, the inverse magnetic catalysis [32]. Thus, magnetic fields can compete with the temperature in nontrivial ways. The motivation for the present work is to understand this interplay of temperature and magnetic field effects as far as the ISB and SNR phenomena are concerned. To our knowledge, no previous work have studied such combined effects in the context of ISB and SNR. Hence, in our opinion, this makes the present study of relevance, in particular in the context of the aforementioned motivations for studying the combined effects of temperature and magnetic fields in phase transitions in complex systems in general.
To carry out the present study, we will make use of the nonperturbative method of the optimized perturbation theory (OPT) [16,[33][34][35][36][37][38][39] (for a recent review on the OPT method, see also Ref. [40] and references therein). This will help us to gauge the consistency of the ISB and SNR phenomena against the perturbative and one-loop expansions when applied to the present problem.
The remainder of this work is organized as follows. In Sec. II, we introduce the model and its implementation in the context of the nonperturbative OPT method. In Sec. III, the effective potential is explicit derived in the OPT method at first order. The relevant expressions at finite temperature and in an external magnetic field are given. Some general results for the critical points and phase structure of the model are given in Sec. IV. In Sec. V we present several numerical results representative of the effects produced by both temperature and external magnetic field. Finally, our conclusions are given in Sec. VI. One appendix is included to show some of the technical details for the renormalization of the model.

II. THE MODEL AND ITS OPT IMPLEMENTATION
We consider a model with two coupled complex scalar fields, symmetric under a global U (1) × U (1) symmetry. The most general Lagrangian density with renormalizable interactions and preserving the global symmetry is given by The stability of the potential requires the self-couplings to satisfy λ φ > 0 and λ ψ > 0, but the intercoupling λ can be negative, provided that λ φ λ ψ > 9λ 2 .
We work with the model described by Eq. (2.1) in the context of the nonperturbative OPT method. The implementation of the OPT is given by an interpolation of the Lagrangian density in the form where L 0 is the Lagrangian density of the free theory, which is modified by an arbitrary mass parameter η, while δ works as a bookkeeping (dimensionless) parameter used only to keep track of the order at which the OPT is implemented and it is set to 1 at the end. Applying the standard interpolation given by Eq. (2.2) to Eq. (2.1) gives where Ω 2 φ = m 2 φ + η 2 φ , Ω 2 ψ = m 2 ψ + η 2 ψ , and η φ,ψ are the OPT (mass) parameters.
As usual in the OPT method, any quantity evaluated at any finite order in δ will depend explicitly on the mass parameters η φ,ψ . These parameters are then fixed by an appropriate variational principle applied to the physical quantity that is being computed. As in most OPT studies, here we make use of the principle of minimal sensitivity (PMS) [33] 1 . As we are interested in the phase structure of the model, we also follow the prescription of previous works [39,[42][43][44][45] of applying the PMS directly to the effective potential computed in the OPT method. Applying the PMS directly to the effective potential is also quite convenient, since already at first order in δ it is able to produce nontrivial solutions leading to nonperturbative results. As demonstrated in previous works [41], this is also very convenient given the fast convergence of the OPT/PMS method 2 . Let then V (k) eff be the effective potential evaluated at order δ k in the OPT. The PMS principle for the present application then reads The optimum valuesη φ ,η ψ are functions of the original parameters of the theory and it is through them that the OPT produces nonperturbative results.

III. THE EFFECTIVE POTENTIAL FOR THE
From the interpolated Lagrangian density equation (2.3), we can construct all contributions to the effective potential in the OPT method. By expressing the complex scalar fields φ and ψ in terms of their real and imaginary components, φ = (φ 1 + iφ 2 )/ √ 2 and ψ = (ψ 1 + iψ 2 )/ √ 2, respectively, we follow the usual prescription for obtaining the effective potential. We first shift the fields around the respective background expectations values, which can be taken along the real components of the fields without loss of generality, φ 1 → φ 1 = φ 1 + φ 0 and ψ 1 → ψ 1 = ψ 1 + ψ 0 , with φ 1 = φ 2 = ψ 1 = ψ 1 = 0 and φ 1 = φ 0 and ψ 1 = ψ 0 . Thus, at first order in the OPT, i.e., at first order in δ, the effective potential 1 For other ways of fixing the mass parameters of the OPT and how they compared to the PMS, see for instance Ref. [41]. 2 We point out that the fast convergence of the OPT/PMS method has been demonstrated in particular to critical three-dimensional models [46]. For a recent discussion about the convergence properties of the OPT method under different optimization procedures, see, e.g., Ref. [40].
is explicitly given by where the last two terms shown in Eq. (3.1), proportional to 1/ , are the contributions at order δ generated by the mass counterterms needed for renormalization (note that throughout this paper, we work with dimensional regularization in the MS scheme). The tree-level contribution in Eq. (3.1), V 0 (φ 0 , ψ 0 ), is given by with the required counterterms at O(δ) in the OPT scheme given by the terms shown in the last line in Eq. (3.2). Note that at the present order of the OPT, one only requires the mass counterterms ∆m φ and ∆m ψ , for the φ and ψ fields, respectively, and a vacuum counterterm ∆V (vertice counterterms are required when carrying out the OPT at second order). Details of the renormalization scheme used here to renormalize Eq. (3.1) are given in the Appendix. Furthermore, in the notation used in Eq. (3.1) for the momentum integrals, which are expressed in Euclidean spacetime, we have in the finitetemperature case that with divergent vacuum momentum integral terms regularized in the MS scheme, with d = 3 − 2 , γ E is the Euler-Mascheroni constant and M is the arbitrary mass regularization scale. The sum in Eq. (3.3) is performed over the Matsubara frequencies for bosons, ω n = 2πnT , n ∈ Z. The field propagators in Eq. (3.1) are such that P 2 + Ω 2 ≡ ω 2 n + E 2 (p), with dispersion relation Likewise, when the external magnetic field B is coupled to the system, we now have instead for the regularized momentum integrals that where the first sum in Eq. (3.5) is over the Landau levels, with k ∈ N and e is the elementary charge. Here, the dispersion relation for the bosons is Note also that in writing Eqs. (3.5) and (3.6), without loss of generality, a gauge was chosen such that the magnetic field is in the z direction.
A. The renormalized effective potential at O(δ) After performing the Matsubara sums, we obtain for the renormalized effective potential the result at O(δ), where we are using a simplified notation in Eq. (3.7) for the Y φ,ψ and X φ,ψ terms, which will depend on the three cases that we will be studying here: At finite temperature and in the absence of an external magnetic field (B = 0), we have that and where At finite temperature and in the presence of an external magnetic field (B = 0), we have that where ζ (s, a) is the s derivative of the Hurwitz zeta function [47], and while for X we now have that where Finally, at zero temperature (T = 0) but in the presence of the external magnetic field, we have that in the above expressions J(B, T ) = 0 and I(B, T ) = 0 in Eqs. (3.12) and (3.15), respectively. Thus, and

22) and
(3.23) We can define the effective square masses of the fields which will include all effects of T and B at O(δ) in the OPT-PMS. From the second derivative in the background fields of the effective potential Eq. (3.7), we obtain Finally, note that Eqs. (3.24) and (3.25), which give the curvature of the effective potential Eq. (3.7) in the φ and ψ directions, respectively, should not be confused with the physical masses for the fields. For instance, the real and imaginary components for the φ and ψ fields, φ 1 , φ 2 , for the complex scalar field φ, and ψ 1 , ψ 2 , for the complex scalar field ψ, will have effective physical masses at O(δ) given, respectively, by With our choice of shifting the fields around their real components, we have that m φ1,eff and m ψ1,eff correspond to the Higgs modes for the φ and ψ fields, respectively, while m φ2,eff and m ψ2,eff correspond to the Goldstone modes for the fields. Note that these masses remain positive definite always. In particular, from the previous equations we can easily demonstrate the validity of the Goldstone theorem in the OPT method (see, in particular, Ref. [44] for an explicit proof in the case of one complex scalar field. An analogous demonstration also follows here). Equations which are quite general results obtained from the application of the OPT-PMS at O(δ).
A. Some simple OPT-PMS results at O(δ) By considering the case (a), T = 0, B = 0, from Sec. III A 1, we have that whenever the system is at the critical point whereφ(T c,φ ) = 0 and m 2 φ,eff (T c,φ ) = 0, then Ω φ → 0. Likewise, in the case of symmetry restoration in the direction of ψ, whereψ(T c,ψ ) = 0 and m 2 ψ,eff (T c,ψ ) = 0, then Ω ψ → 0. In each of these cases, from Eq. (3.10) we find that (4.4) Note that due to the difference in sign for the function X in each of the above cases, we will have opposite effects on the symmetry behavior of the system, whenever thermal effects or the external magnetic field dominates. It is also useful to combine Eqs. (3.22) and (3.23) and also use Eqs. (3.24) and (3.25) to arrive at for m 2 ψ,eff > 0, (4.6) where in the effective masses, Eqs. (3.24) and (3.25), X φ and X ψ are evaluated at η φ =η φ and at η ψ =η ψ , respectively [i.e., with Ω φ and Ω ψ as given by Eqs. (4.1) and (4.2) at O(δ)].
From the above expressions, we can for example set the coupling between φ and ψ to zero, λ = 0. Let us consider for example the case (a), with T = 0, B = 0. Assuming symmetry restoration in the directions of the φ and ψ fields happening at T c,φ and at T c,ψ , whereφ(T c,φ ) = 0 andψ(T c,ψ ) = 0, respectively, we get from, e.g., Eq.
Similarly to the case above, where we set the coupling between φ and ψ to zero, λ = 0, we can make an analogous analysis for the case (c), with T = 0, B = 0. In this case, symmetry restoration in the φ and ψ directions happening at B c,φ and at B c,ψ , whereφ(T c,φ ) = 0 and ψ(T c,ψ ) = 0, respectively. We now get from Eq. (4.1) that Ω φ (B = B c,φ ) = m 2 φ,eff (B = B c,φ ) = 0. Then, from Eq. (4.4), we now get that eB c,φ = 24π 2 m 2 φ /(λ φ ln 2). Analogously, using Eq. (3.25), we have that eB c,ψ = 24π 2 m 2 ψ /(λ ψ ln 2). While in the case (a) symmetry is restored at some finite temperature when the system is initially in the broken phase (m 2 φ,ψ < 0) at T = 0, we have the opposite situation in an external magnetic field. Starting with the system in the symmetric phase (m 2 φ,ψ > 0) at B = 0, it will go to a symmetry broken phase at a sufficiently large external magnetic field. This is reminiscent of the magnetic catalysis effect [31].
When the intercoupling is nonvanishing, λ = 0, the situation gets more involved. The proper analysis in this case can only be studied numerically in the present study of the OPT-PMS. However, we can still get some useful results when analyzing the model in the usual perturbation theory and using either a high-temperature expansion, m φ,ψ /T 1, or a high magnetic field approximation, m 2 φ,ψ /(eB) 1. This is useful for later comparison with the OPT results.

B. Perturbation theory results
It is useful to compare the results to be obtained with the OPT for the model studied here with those obtained using perturbation theory in the high-temperature and/or high magnetic field case. In this case, the expressions for the effective masses Eqs. (3.24) and (3.25) still hold and we can also use Eqs. (4.5) and (4.6), where we set η φ = η ψ = 0 in those expressions whenever m 2 φ > 0 and m 2 ψ > 0, i.e., when we start in the symmetry restored phase in both φ and ψ directions. The case of starting in the symmetry broken phase for both φ and ψ, i.e., m 2 φ < 0 and m 2 ψ < 0, is worked out more conveniently within the one-loop approximation, instead of using perturbation theory, and it is given in the next subsection below. In the high-temperature approximation, m φ,ψ /T 1, the thermal integral contributing to the effective masses can be approximated at leading order as [48] Similarly, in the case of a finite magnetic field, using Eq. (3.18), we get in the high magnetic field approximation that When λ < 0 and choosing coupling constants such that, e.g., 2λ φ /3+λ < 0, because of the boundedness condition λ φ λ ψ > 9λ 2 , it necessarily implies that 2λ ψ /3 + λ > 0. Hence, when m 2 φ > 0, i.e., initially choosing the system to be in the symmetry restored phase in the φ direction, Eq. (4.9) will imply that there will be an ISB at the approximated critical temperature where we have considered thatψ(T = T c,φ ) = 0, i.e., the symmetry in the ψ direction is restored at T c,ψ < T c,φ , when m 2 ψ < 0, and remains restored, or when m 2 ψ > 0, in which case we always haveψ = 0. In the opposite case, when m 2 φ < 0, i.e., initially choosing the system to be in the symmetry broken phase in the φ direction, then it will remain in this state at arbitrarily large temperatures 3 . This is the case of SNR. In the case of considering a finite magnetic field, B = 0, but at T = 0, the situation is similar, although the roles of ISB and SNR get reversed because of the minus sign in Eq. (4.8). The situation is analogous when choosing that there is a symmetry restoration in the ψ direction, i.e., considering that m 2 φ , m 2 ψ > 0, with 2λ ψ /3+λ < 0, which now implies that 2λ φ /3 + λ > 0. Hence, the symmetry remains restored in the φ field direction at high temperatures, but in the ψ field direction there is ISB with the critical temperature for ψ given by As already discussed before, in the case of a finite magnetic field and at zero temperature, because of the minus sign in Eq. (4.8) the situation gets reversed with respect to what it is obtained in the finite-temperature case. In the case of perturbation theory at large magnetic fields, we can extend the results in Eqs. (4.9) and (4.10) by making the replacement T 2 /12 → −eB ln(2)/(4π) 2 . Hence, starting in the symmetry restored phase for both the φand ψ-field directions, i.e., considering m 2 φ > 0 and m 2 ψ > 0, and assuming 2λ φ /3 + λ < 0 (with λ < 0), we have that the symmetry will tend to remain restored in the φ direction, while in the ψ direction there is ISB at a critical magnetic field, In the case of starting in the symmetry broken phase for both field directions, m 2 φ < 0 and m 2 ψ < 0, and assuming again that 2λ φ /3 + λ < 0 the situation gets more involved since the results also depend on how the background fieldsφ andψ will behave at large B. As we are going to see numerically in Sec. V B, the tendency is both φ and ψ to remain in a symmetry broken state.

C. One-loop effective potential
In the case of starting in the symmetry broken phase for both φ and ψ, i.e., m 2 φ < 0 and m 2 ψ < 0, it is more convenient to work with the one-loop effective potential. The one-loop effective potential for the background fields φ 0 and ψ 0 given by the general expression where M 2 i (φ 0 , ψ 0 ) denotes the physical masses eigenstates of the model, i.e., expressing the mass eigenvalues in the Higgs and Goldstone basis. These are explicitly given by the diagonalization of the quadratic mass matrix for the fields, whose mass eigenvalues are When substituting the background fields by their treelevel vacuum expectation values, we can recognize that M 1 and M 2 are the two Higgs modes, associated with the φ 1 and ψ 1 real components for the two complex scalar fields, while M 3 (φ 0 , ψ 0 ) = M 4 (φ 0 , ψ 0 ) = 0 are the corresponding two Goldstone modes, associated with the φ 2 and ψ 2 imaginary components of φ and ψ. The renormalized one-loop effective potential is then given by In the next section we will study numerically the results obtained using the OPT-PMS and contrast these results with those obtained when using perturbation theory and/or using the one-loop approximation.

V. ISB AND SNR IN THE OPT-PMS NONPERTURBATIVE METHOD
Let us start the exploration of the results for ISB and SNR for the present model in the OPT-PMS nonperturbative method. For comparison purposes and for easiness of analysis, we will divide our presentation for the cases of T = 0 and B = 0, T = 0 and B = 0 and, finally include the combined effects of temperature and external magnetic field.
A. ISB and SNR in the OPT-PMS nonperturbative method at T = 0 and B = 0 In the case of only including finite-temperature effects, we need the thermal functions defined in Sec. III A 1. The respective momentum integrals are solved numerically directly. We start by considering the case with m 2 φ > 0 and m 2 ψ > 0. For illustration purposes, we choose the couplings such that it predicts ISB in the φ direction, while in the ψ direction the system remains in the symmetry restored phase.
In the example considered in Fig. 1(a) there is an ISB transition in the φ direction at a critical temperature T c,φ /M 29.3 in the OPT-PMS case, while for PT the critical temperature is smaller, T c,φ /M 26.1. In Fig.1(b) it is shown how the critical temperature for ISB in each case behaves as a function of the intercoupling λ. Note that ISB tends to disappear at less negative values of λ (which is a consequence of approaching the boundedness condition for the couplings). This is signaled by a diverging behavior of T c , which happens at a smaller (in modulus) value for λ in the PT case as compared to the OPT.
We consider the same conditions used in Fig.1 where we have seen ISB in the φ direction, while ψ remains in its symmetry restored phase, also in Fig. 2. In Fig. 2 we now show the effective masses Eqs. (3.24) and (3.25) and we contrast them with the physical ones, given by Eqs. (3.26), (3.27), (3.28), and (3.29) in the OPT-PMS scheme. We can see from Fig. 2(a) that at the transition point the effective mass square m 2 φ,eff changes from positive (symmetry restored) to negative (symmetry breaking), while m 2 ψ,eff remains positive throughout the temperature range shown, indicating that the symmetry remains restored in that direction. In Fig.2(b) we see that the Goldstone theorem applies, as it should, for each of the field directions. In this case, both Higgs and Goldstone modes agree with each other for all range of temperature in the ψ direction. For the φ direction, the Higgsand the Goldstone-like masses agree with each other in the φ direction for temperatures T ≤ T c,φ . Now, let us consider the case of SNR. Hence, we consider m 2 φ < 0 and m 2 ψ < 0, i.e., we start with the system in the symmetry broken state in both the field directions. In Fig. 3(a) we show both background fieldsφ andψ as a function of the temperature in the OPT case. The parameters chosen are such that m 2 φ = m 2 ψ = −M 2 < 0, with the same couplings as in Fig. 1, λ φ = 0.018, λ ψ = 0.6, and λ = −0.03. For these parameters the symmetry remains broken, i.e., there is SNR in the direction of φ, while there is a symmetry restoration at a critical temperature in the direction of ψ. For the parameters considered in this example, we have that T c,ψ /M 16.4 in the OPT case, while in the one-loop approximation we find T c,ψ /M 15.9. In this same example showing SNR, it is useful to give the effective masses Eqs. (3.24) and (3.25) and to contrast them with the physical ones, given by Eqs. (3.26), (3.27), (3.28) and (3.29) in the OPT scheme. These are shown in Figs. 3 (b) and 3(c), respectively. We can clearly see from Fig. 3(b) that at the transition point the effective mass square m 2 ψ,eff changes from negative (symmetry breaking) to positive (symme- try restored) at the critical temperature T c,ψ /M 16.4, while m 2 φ,eff remains negative throughout the temperature range shown, indicating SNR in the direction of φ.
In Fig. 3(c) we see that the Goldstone theorem applies, as it should, for each of the field directions. The physical masses remain positive definite as also expected. In particular, after symmetry restoration (in the direction of ψ), both the Higgs and Goldstone modes agree with each other, m ψ1,eff = m ψ2,eff for T ≥ T c,ψ . Turning now to the case of zero temperature, but finite external magnetic field, we start by looking at the results in the case when m 2 φ > 0 and m 2 ψ > 0, i.e., starting in the symmetry restored phase for both φ and ψ. The couplings are again chosen such that 2λ φ /3 + λ < 0 (with λ < 0). A representative example of this case is shown in Fig. 4. As expected, the situation gets reversed with respect to what has been shown in the corresponding case at finite temperature but zero external magnetic field and shown in Fig. 1. Here, we have ISB in the direction of the ψ field, while the symmetry remains restored in the φ direction, i.e.,φ = 0 throughout the range of magnetic field values considered. In the direction of the ψ field, the critical magnetic field, for the case of the parameters λ φ = 0.018, λ ψ = 0.6, and λ = −0.03 is given by eB c,ψ /M 2 = 615.1 in the case of the OPT, while in PT we obtain that eB c,ψ /M 2 = 623.1. Let us now consider the case of starting in the symmetry broken phase for both field directions, i.e., taking m 2 φ < 0 and m 2 ψ < 0. Again we assume that 2λ φ /3 + λ < 0. The results in this case are shown in Fig. 5. Note here that both φ and ψ remain in a symmetry broken state as a consequence of the magnetic field favoring symmetry breaking. This is akin to the magnetic catalysis effect seen in general due to a magnetic field, which tends to enhance the symmetry breaking [31]. In Fig. 5 we have considered both OPT and the one-loop approximation for comparison purposes.
C. ISB and SNR in the OPT-PMS nonperturbative method at T = 0 and B = 0 Let us now consider the effects from both temperature and external magnetic field. As already observed, this is a particularly interesting case, since there is a competition between the thermal effects and the magnetic field, which one acting in an opposite direction as far as symmetry breaking and restoration are concerned. Before entering in our numerical results for this case, let us recall that in the expressions involving the magnetic field, e.g., Eqs. (3.14) and (3.16), the thermal integrals also require a sum over the Landau levels. While at large magnetic field values, eB Ω 2 , one only requires to sum over a few Landau levels for proper convergence of the expressions, in the weak magnetic field regime, eB Ω 2 , one typically requires to consider a large number of levels in the sum. There have been a few, but very reliable approaches in the literature to handle this issue in the weak magnetic field regime. For example, Ref. [44] made use of the Euler-Maclaurin formula as a way to work out the sum over the large number of Landau's levels. In Ref. [49], the authors have proposed instead a weak magnetic field approximation for the bosonic propagator. In this case, the propagator in Euclidean spacetime for the charged scalar fields in the presence of an external magnetic field is approximated as [49] 1 where E 2 (p) = p 2 + Ω 2 , with p 2 = p 2 z + p 2 ⊥ and p 2 ⊥ = p 2 x + p 2 y . In this work, we have made use of both the Euler-Maclaurin formula considered in Ref. [44] as also the weak magnetic field approximation for the bosonic propagator given by Eq. (5.1). In the weak magnetic field regime eB Ω 2 both approaches are found to agree quite well. For the case of strong magnetic fields, we have considered a sufficient number of Landau's levels such to have convergence for the results. With these due cares taken into account, we now present our results.
In Fig. 6 we explore the effect of the external magnetic field in the ISB case for the φ field. Note that from the result shown in Fig. 6(a), PT exhibits a much stronger departure from the OPT result for the background field φ as the temperature increases beyond the critical temperature for ISB when in a strong magnetic field regime.  In Fig. 6(b) we show how the critical temperature for ISB in this example changes with an increasing magnetic field. One notices that for very strong fields, the critical point is obtained in the OPT and PT approximations as they approach each other, showing that nonperturbative effects brought by the OPT tend to be less important.
Finally, in Fig. 7, we study the SNR case when in the presence of strong magnetic fields. Here, we compare again the results from the OPT with the one-loop approximation and the parameters are chosen such that there is symmetry restoration in the direction of the ψ field, while the symmetry remains broken (SNR) in the direction of φ. In the ψ direction and for an external magnetic field of eB/M 2 = 10 3 , the critical temperature for symmetry restoration in the direction of ψ is found to be T c,ψ /M 23.9M in the OPT case, while in the oneloop approximation we find that T c,ψ /M 17.7. The way the critical temperature T c,ψ changes in both cases as a function of the external magnetic field is shown in Fig. 8. Here, opposite to the previous case of ISB shown in Fig. 6(b), the one-loop approximation tends to underestimate the critical temperature when compared to the OPT result. The difference between the results also increases the larger is the magnetic field.

VI. CONCLUSIONS
In this work, we have shown that the introduction of an external magnetic field can induce significant changes in the phase structure of a system composed by two complex scalar fields with both self-and inter-interactions. This is a consequence of the fact that thermal effects tend to act oppositely to those obtained for an external magnetic field. This situation gets particularly more involved whenever there is a possibility of ISB and/or SNR phenomena, which results by an appropriate choice for the parameters of the model, e.g., the coupling constants, in particular for the intercoupling between the fields.
In our study we have used the effective potential for the analysis of the phase structure of the model. The derivation of the effective potential and the results obtained from it were performed considering the nonperturbative method of the OPT. The results obtained in the OPT method were compared to those obtained in the perturbation theory as well as with those obtained when using the one-loop approximation for the effective potential.
Our results show that the effect of the magnetic field tends to always increase the critical temperature, making symmetry breaking more easily achieved whenever symmetry restoration is involved. In the case of ISB, the critical temperature also increases with the magnetic field, although here we can interpret that the external magnetic field makes it more difficult to achieve ISB, thus producing a higher T c .
In summary, our results still corroborate with the existence of ISB/SNR, even when including the effects of an external magnetic field and in the context of the nonperturbative OPT scheme. In the absence of a magnetic field, the critical temperature in the OPT case tends always to be larger than in the perturbation theory, or in the one-loop approximation, for both ISB and SNR. However, in the absence of thermal effects, but in a finite magnetic field, the critical magnetic field tends to be smaller in the OPT than in the perturbative and loop approximation cases. When both the external magnetic field and thermal effects are present, the differences between the OPT and the perturbative and loop approximation cases vary according to the magnitude of the external magnetic field. Nevertheless, the presence of an external magnetic field, due to the magnetic catalysis effect, pushes the energy scale for phase transition at finite temperature in both cases to be larger than in the absence of external fields.
As a possible future extension of this work, one could for example try to study the issue of inverse magnetic catalysis in the decoupled case, i.e., when λ = 0, in the OPT. Inverse magnetic catalysis seems only to be possible to be realized in a nonperturbative context [32]. In particular, by having effective coupling constants that include both thermal and magnetic field effects, it has been shown to be able to manifest inverse magnetic catalysis, as, e.g., in the model of a complex scalar field studied in Ref. [50]. A similar setting could possibly also be studied in the context of the OPT, but likely only when carrying out the calculations up to at least second order in the OPT, such that vertices (couplings) could also receive thermal and magnetic field effects. This is an interesting but much more technical and difficult problem that we hope to address in the future.

Appendix A: Renormalization
In the OPT scheme one needs only the standard renormalization terms for fully renormalizing the effective potential for the model, e.g., vacuum, masses, and coupling constant counterterms. These counterterms need then to be derived at the appropriate order in the OPT method. In the present case, we have derived the effective potential up to order δ in the OPT scheme, which implies requiring only a vacuum and appropriate masses counterterms [counterterms for the couplings are only necessary when carrying out the derivation at O(δ 2 ) and higher].
To fully renormalize the model at O(δ) we need the mass counterterms: These two mass counterterms, in particular, generate the last two terms shown in Eq. (3.1) that also contribute at O(δ 1 ). In particular, note that these terms are fundamental to remove the crosses divergence with temperature dependence that appears in the computation of the two-loop contributions in Eq. (3.1).
Finally, at the order δ in the OPT, the renormalization is completed with the inclusion of the vacuum counterterms,