Numerical Solution for the Mean-Field Critical Curves of the $SO(\mathfrak{N})$ Polynomial $\lambda\varphi^{4}$ Models

Jorge L. deLyra

Department of Mathematical Physics, Physics Institute, University of São Paulo

Version 1, May 2006

One of the important quantities one can calculate using the $N=1$ mean-field method for the polynomial models is the position of the critical curve within the parameter plane of each model. This is done in section 3 of chapter 2 of the book [1], for the case of the $O(1)$ model. The result can be extended to the $SO(\mathfrak{N})$ models, and in this general case the $N=1$ mean-field result for the critical curve $\lambda(\alpha)$ of the $\lambda\varphi^{4}$ polynomial model is:

$$\frac{\mathfrak{N}}{2d}=\frac{\displaystyle \int_{0}^{\infty... ...left[(d+\alpha/2)\varphi^{2} +(\lambda/4)\varphi^{4}\right]}},$$ (1)

which determines $\lambda(\alpha)$ implicitly. Since it is not possible to solve this equation analytically in order to obtain $\lambda(\alpha)$ in explicit form, we will be concerned here with the numerical solution of this equation, so that we may plot the graph of this function. The integrals can be written in terms of the parabolic cylinder functions ${\bf D}_{\nu}$ [3],

$$\frac{\sqrt{2\lambda}}{2d}=\frac{{\bf D}_{-\left(\frac{\math... ...}}{2}\right)} \left(\frac{2d+\alpha}{\sqrt{2\lambda}}\right)}.$$

However, this does not help us to solve the equation numerically because these functions are not readily available in a computer implementation. Note that for even $\mathfrak{N}=2p$ the integrals can be written as integrals on a variable $\chi=\varphi^{2}$,

$$\frac{\mathfrak{N}}{2d}=\frac{\displaystyle \int_{0}^{\infty... ...1} \;e^{-\left[(d+\alpha/2)\chi +(\lambda/4)\chi^{2}\right]}}.$$

Hence it is not surprising that for even $\mathfrak{N}$ they can be written in terms of the error function, which can in fact be done through the ${\bf D}_{-n}$ functions, with integer $n$. The ${\bf D}_{-n}$ functions can be written in terms of the error function $\Phi(x)$, which is given by [4]

$$\Phi(x)=\frac{2}{\sqrt{\pi}}\int_{0}^{x}{\rm d}t\;e^{-t^{2}}.$$

For the first few functions ${\bf D}_{-n}$ we have [5]

$$\begin{eqnarray*} {\bf D}_{0}(x) & = & e^{-x^{2}/4} {\bf D}_{-1}(x) & = & \sq... ...)\left[1-\Phi\left(\frac{x}{\sqrt{2}}\right)\right] \right\} \\ \end{eqnarray*}$$

From the first two one can get all the functions in this sequence by repeated use of the recurrence relation [6]

$${\bf D}_{p+1}(x)-x{\bf D}_{p}(x)+p{\bf D}_{p-1}(x)=0,$$

which can also be written in the form

$${\bf D}_{-(n+2)}(x)=\frac{{\bf D}_{-(n)}(x)-x{\bf D}_{-(n+1)}(x)}{n+1}.$$

Using all this one can solve the equation numerically for the case of even $\mathfrak{N}$, one case at a time, using the error functions erf or derf, which are readily available for use with the g77 Fortran compiler. For example, for $\mathfrak{N}=2$, which is the simplest case, the equation is

$$\sqrt{\frac{2}{\pi}}\;e^{-\frac{(2d+\alpha)^{2}}{4\lambda}}=... ...t[1-\Phi\left(\frac{2d+\alpha}{2\sqrt{\lambda}}\right)\right].$$

Given a value of $\alpha$, one can solve this equation numerically for $\lambda$, and vice-versa. However, this solves only the even-$\mathfrak{N}$ cases, not the odd-$\mathfrak{N}$ cases, including $\mathfrak{N}=1$.

For all values of $\mathfrak{N}$ one can make a change of variables in order to simplify the original integrals in a way that is convenient for numerical purposes. Since the solution for $\lambda=0$ is obviously $\alpha=0$, we can assume that $\lambda>0$ and consider the variable $\chi$ given by

$$\chi=\left(\frac{\lambda}{4}\right)^{1/4}\varphi,\mbox{    }... ...\rm d}\varphi=\left(\frac{4}{\lambda}\right)^{1/4}{\rm d}\chi,$$

in terms of which the equation becomes

$$\frac{\mathfrak{N}}{2d}=\frac{2}{\sqrt{\lambda}}\;\frac{\dis... ...ft[\frac{2d+\alpha}{\sqrt{\lambda}}\chi^{2}+\chi^{4}\right]}}.$$

Defining now a new parameter $\xi$ given by

$$\xi=-\frac{2d+\alpha}{2\sqrt{\lambda}}\Longrightarrow \alpha=-2\left(\xi\sqrt{\lambda}+d\right),$$

we may write the equation in the form

$$\sqrt{\lambda}=\frac{4d}{\mathfrak{N}}\;\frac{\displaystyle ... ...\rm d}\chi\;\chi^{\mathfrak{N}-1}\;e^{2\xi\chi^{2}-\chi^{4}}},$$

or, multiplying above and below by $\exp(-\xi^{2})$,

$$\sqrt{\lambda}=\frac{4d}{\mathfrak{N}}\;\frac{\displaystyle ... ...}{\rm d}\chi\;\chi^{\mathfrak{N}-1}\;e^{-(\chi^{2}-\xi)^{2}}}.$$

Given a value of $\xi$, this allows one to calculate $\sqrt{\lambda}$ using a single pair of numerical integrations, and afterward to calculate the corresponding $\alpha$ from $\xi$, $d$ and $\sqrt{\lambda}$. The parameter $\xi$ is a new parameter for the critical curve. In order to determine its range of variation as we travel along the curve from end to end, consider first the asymptotic line of the critical curve, for large values of $-\alpha$ and $\lambda$, which can be obtained through the representation of the equation in terms of the parabolic cylinder functions [7], and is given by

$$\lambda=-\frac{2d+\alpha}{\beta_{c}},\mbox{    }\beta_{c}=\frac{\mathfrak{N}}{2d},$$

where $\beta_{c}$ is the critical point of the sigma model over the arc at infinity. In this case we have for $\xi$ the expression

$$\xi=-\frac{2d+\alpha}{2\sqrt{\lambda}} =\frac{\beta_{c}\lambda}{2\sqrt{\lambda}} =\frac{1}{2}\beta_{c}\sqrt{\lambda}.$$

Therefore $\xi\rightarrow\infty$ corresponds to $\lambda\rightarrow\infty$ and hence to the sigma-model critical point over the arc at infinity, and when

$$\lambda\rightarrow\infty\mbox{   we have   } \xi\rightarrow\infty\mbox{  as  }\sqrt{\lambda}.$$

If we now consider the tangent to the critical curve at the Gaussian point, which can be obtained by differentiating equation (1) implicitly with respect to $\alpha$ and $\lambda$, then expressing the resulting integrals in terms of the $\Gamma$ function [8], and which is given by

$$\lambda=-\frac{\alpha}{\beta_{c}},\mbox{    }\beta_{c}=\frac{\mathfrak{N}}{2d},$$

we see that in this case we get for $\xi$ the expression

$$\xi=-\frac{2d+\alpha}{2\sqrt{\lambda}} =-\frac{d}{\sqrt{\lam... ...=-\frac{d}{\sqrt{\lambda}}+\frac{1}{2}\beta_{c}\sqrt{\lambda}.$$

Therefore $\xi\rightarrow-\infty$ corresponds to $\lambda\rightarrow 0$ and hence to the Gaussian point, and when

$$\lambda\rightarrow 0\mbox{   we have   } \xi\rightarrow-\infty\mbox{  as  }\frac{1}{\sqrt{\lambda}}.$$

In short, we have the following limiting behaviors:

$$\begin{array}{cclcccr} -\alpha\mbox{ and }\lambda & \rightar... ...splaystyle \sqrt{\lambda}} & \rightarrow & -\infty. \end{array}$$

We have therefore the following general algorithm for the numerical solution of the equation:

1. Choose a value for $\xi$ in $(-\infty,\infty)$.

2. Calculate $\sqrt{\lambda}$ by numerical integration,

   
   

   $$\sqrt{\lambda}=\frac{4d}{\mathfrak{N}}\;\frac{\displaystyle ... ...}{\rm d}\chi\;\chi^{\mathfrak{N}-1}\;e^{-(\chi^{2}-\xi)^{2}}}.$$
   
   

3. Calculate the corresponding $\alpha$ by the formula

   
   

   $$\alpha=-2\left(\xi\sqrt{\lambda}+d\right).$$
   
   

The resulting pair $(\alpha,\lambda)$ is a point of the critical curve of the $SO(\mathfrak{N})$ model, in $d$ dimensions. Probably the best way to choose values of $\xi$ so that the resulting points $(\alpha,\lambda)$ are distributed fairly homogeneously along the critical curve is to use the equation of the tangent line at the Gaussian point, choosing some $\varepsilon$, integer values of $n=1,2,3,\ldots$, and using

$$\xi_{n}=\frac{\mathfrak{N}}{4d}\sqrt{n\varepsilon}-\frac{d}{\sqrt{n\varepsilon}}.$$

Since the functions being integrated are centered around a single maximum that moves to large values of the variable when $\xi$ is large, in order to integrate them efficiently we must know the location of the points of maximum of the functions. The basic idea for the integration of one of these functions is to start at the point of maximum of the function, and then to integrate to both sides, down towards zero and up towards infinity. Let us consider then the functions

$$\chi^{p}\;e^{-(\chi^{2}-\xi)^{2}},$$

where $p\geq 0$ is an integer. Except for the case $p=0$ this is zero for $\chi=0$, and in all cases for $\chi\rightarrow\infty$. Taking the derivative we get

$$\chi^{p-1}[p+4\xi\chi^{2}-4\chi^{4}]\;e^{-(\chi^{2}-\xi)^{2}}.$$

Besides vanishing at $0$ (for $p>1$) and at $\infty$, which are clearly points of minimum, this can be zero for

$$p+4\xi\chi^{2}-4\chi^{4}=0 \Longrightarrow \chi^{2}=\frac{\xi\pm\sqrt{\xi^{2}+p}}{2},$$

and since $\chi^{2}\geq 0$ only the $+$ sign can be kept, and we get

$$\chi_{\rm max}=\sqrt{\frac{\xi+\sqrt{\xi^{2}+p}}{2}}.$$

This is correct in all cases including $p=0$. For $p=0$ we have that, if $\xi>0$, then $\chi_{\rm max}=\sqrt{\xi}$, and if $\xi\leq 0$, then $\chi_{\rm max}=0$. Note that the only case in which the maximum is at zero is $p=0$ with $\xi\leq 0$, in all other cases the maximum is at a positive and non-vanishing value of $\chi$.

The functions examined above have the form of a pulse around their points of maximum, but this pulse has a variable width, which decreases as the parameter $\xi$ increases to positive values. This makes is harder to find an appropriate numerical integration interval for them. For the cases $p>0$ we can perform one more transformation of the integrals, in order to improve this situation. We make a transformation of integration variables to the variable $\gamma=\chi^{2}$ and hence write

$$\int_{0}^{\infty}{\rm d}\chi\;\chi^{p}\;e^{-(\chi^{2}-\xi)^{... ...{\rm d}\gamma \;\gamma^{\frac{p-1}{2}}\;e^{-(\gamma-\xi)^{2}},$$

where the integrand has now a fairly constant width around its maximum. In this case the location of the point of maximum is a bit different, that is, $\gamma_{\rm max}$ is not exactly $\chi_{\rm max}^{2}$. The integrand is now

$$\gamma^{\frac{p-1}{2}}\;e^{-(\gamma-\xi)^{2}},$$

and its derivative is given by

$$\left[\frac{p-1}{2}\gamma^{\frac{p-3}{2}} -2(\gamma-\xi)\gamma^{\frac{p-1}{2}}\right]e^{-(\gamma-\xi)^{2}}.$$

Note that, due to the possibly negative powers of $\gamma$, the cases $p=1$ and $p=2$ should be examined separately. For $p=1$ we have

$$-2(\gamma-\xi)e^{-(\gamma-\xi)^{2}},$$

which is zero at the point of minimum at infinity and at $\gamma=\xi$, which is therefore the point of maximum. Since $\gamma\geq 0$ this is valid only for $\xi\geq 0$, for $\xi<0$ the function is monotonically decreasing with its maximum at $\gamma=0$. For $p=2$ we have

$$\left[\frac{1}{2}\gamma^{-\frac{1}{2}} -2(\gamma-\xi)\gamma^{\frac{1}{2}}\right]e^{-(\gamma-\xi)^{2}}.$$

Note that in this case the function has an infinite derivative at the point $\gamma=0$, which makes it unsuitable for numerical integration. Its point of maximum can be obtained as one of the roots of a quadratic polynomial and is in fact included in the general case, which we proceed to analyze. If we multiply it by $\gamma^{3/2}$, the derivative of the general case becomes

$$\gamma^{\frac{p}{2}} \left[\frac{p-1}{2}-2(\gamma-\xi)\gamma\right]e^{-(\gamma-\xi)^{2}}.$$

Except for the cases $p=1$ and $p=2$ this vanishes at the points of minimum at $\gamma=0$ and at $\gamma\rightarrow\infty$, and for all cases it vanishes at the point of maximum given by

$$\frac{p-1}{4}+\xi\gamma-\gamma^{2}=0 \Longrightarrow \gamma=\frac{\xi\pm\sqrt{\xi^{2}+(p-1)}}{2},$$

where since $\gamma\geq 0$ only the $+$ sign can be kept, and we get

$$\gamma_{\rm max}=\frac{\xi+\sqrt{\xi^{2}+(p-1)}}{2}.$$

This is in fact correct in all cases including $p=1$ and $p=2$. For $p=1$ we see that this gives the correct answer: if $\xi>0$, then $\gamma_{\rm max}=\xi$, and if $\xi\leq 0$, then $\gamma_{\rm max}=0$. Note that the only case in which the maximum is at zero is $p=1$ with $\xi\leq 0$, in all other cases the maximum is at a positive and non-vanishing value of $\gamma$.