Emergence of a bicritical end point in the random crystal field Blume-Capel model

TL;DR

This paper explores how disorder affects the phase diagram of the Blume-Capel model, revealing a transition from tricritical to bicritical end points as disorder strength increases.

Contribution

It introduces the phase diagram of the disordered Blume-Capel model, highlighting the emergence of bicritical end points with increasing disorder.

Findings

Tricritical point exists only at weak disorder.

Bicritical end point appears at higher disorder levels.

Phase diagram complexity increases with disorder strength.

Abstract

We obtain the phase diagram for the Blume-Capel model with bimodal distribution for random crystal fields, in the space of three fields: temperature, crystal field and magnetic field. We find that three critical lines meet at a tricritical point, but only for weak disorder. As disorder strength increases there is no tricritical point in the phase diagram. We instead find a bicritical end point, where only two of the critical lines meet on a first order surface in the H=0 plane. For intermediate strengths of disorder, the phase diagram has critical end points along with the bicritical end point. One needs to look at the phase diagram in the space of three fields to identify various such multicritical points.

Tables

Table 1. Table 1 : Width of the wing lines for different p 𝑝 p . T l ​ c subscript 𝑇 𝑙 𝑐 T_{lc} and Δ l ​ c subscript Δ 𝑙 𝑐 \Delta_{lc} represent the values of T 𝑇 T and Δ Δ \Delta for H = 0 𝐻 0 H=0 where the λ + subscript 𝜆 \lambda_{+} and λ − subscript 𝜆 \lambda_{-} lines meet and T u ​ c subscript 𝑇 𝑢 𝑐 T_{uc} is the value along the critical line as Δ → ∞ → Δ \Delta\rightarrow\infty and H → ∞ → 𝐻 H\rightarrow\infty .

$p$	$T_{lc}$	${\Delta }_{lc}$	$T_{uc}$	$\delta T$
0.0453	0.28043	0.501175	0.23866	0.0417665
0.05	0.276396	0.50468	0.237473	0.038923
0.07	0.26185	0.518896	0.23245	0.029399
0.1	0.2451	0.538417	0.224972	0.020128
0.2	0.2058	0.596376	0.2	0.0058
0.3	0.17643	0.6490843	0.174978	0.001452
0.4	0.15024	0.69968	0.1499	0.000248
0.5	0.125016	0.7499884	0.12498	0.000036

Table 2. Table 2 : Co-ordinates of the BEP and CEP’s for 0.022 < p < 0.107 0.022 𝑝 0.107 0.022<p<0.107 .

$p$	$T_{BEP}$	${△}_{BEP}$	$T_{CEP1}$	${△}_{CEP1}$	$T_{CEP2}$	${△}_{CEP2}$
0.03	0.2961208	0.489187	0.295197	0.489166	0.03	0.4977229
0.044	0.28153	0.500195	0.27585	0.500186	0.04401	0.521585
0.07	0.26185	0.518896	0.24036	0.519398	0.07099	0.533953
0.107	0.24166	0.542	0.15972	0.547514	0.13975	0.549068

Equations

$$H(C_N)=-\frac{1}{2N}(\sum _i{s}_i{)}^2-\sum _i{\Delta }_i{s}_i^2-H\sum _i{s}_i$$

$$P({\Delta }_i)=p\delta ({\Delta }_i-\Delta )+(1-p)\delta ({\Delta }_i+\Delta )$$

$$P(C_N:\sum _i{s}_i=x_{1}N;\sum _i{s}_i^2=x_{2}N)\sim \exp (-NI(x_1,x_2))$$

$$I(x_1,x_2)$$

$$\frac{x_2}{z}=\frac{p{e}^{\beta \Delta }}{1+z{e}^{\beta \Delta }}+\frac{(1-p)e^{-\beta \Delta }}{1+z{e}^{-\beta \Delta }}$$

$$\tanh (\beta (m+H))=\frac{m}{q}$$

$$z=\frac{2}{\sqrt{1-m^2/q^2}}$$

$$\frac{q}{z}=\frac{p{e}^{\beta \Delta }}{1+z{e}^{\beta \Delta }}+\frac{(1-p)e^{-\beta \Delta }}{1+z{e}^{-\beta \Delta }}$$

$$5-4\beta =2(\beta p-1)e^{\beta \Delta }+2(\beta -\beta p-1)e^{-\beta \Delta }$$

$$12\beta -17+(3\beta -10)\cosh (\beta \Delta )-3\beta (1-2p)\sinh (\beta \Delta )=0$$

$$cosh(\beta \Delta )=\frac{12\beta -19}{8}$$

$$\tilde{f}(m)$$

$$m$$

$$\frac{p(2x^2+xy)}{(1+2xy{)}^2}+\frac{(1-p)(2+xy)}{(x+2y{)}^2}=\frac{1}{2\beta }$$

$$\frac{p(x-8x^3-2x^{2}y)}{(1+2xy{)}^3}+\frac{(1-p)(x^2-8-2xy)}{(x+2y{)}^3}=0$$

$$\frac{2+xy}{[x+2y{]}^2}=\frac{1}{2\beta }$$

$$\frac{x^2-8-2xy}{[x+2y{]}^3}=0$$

$$y$$

$$x$$

$$m$$

$$H=\pm \frac{1}{\beta }\log (\frac{\beta -2+\sqrt{{\beta }^2-3\beta }}{\sqrt{4-\beta }})-m$$

$$m=\pm 2\sqrt{y^{2}-1}\Big{[}\frac{px}{1+2xy}+\frac{(1-p)}{y+2x}\Big{]}$$

$$\tilde{f}(m)=a_2{m}^2+a_4{m}^4+a_6{m}^6+a_8{m}^8$$

$$\displaystyle a_{2}=\frac{\beta}{2}\Bigg{(}1+\frac{2\beta(p-1)}{2+e^{\beta\Delta}}-\frac{2\beta pe^{\beta\Delta}}{1+2e^{\beta\Delta}}\Bigg{)}$$

$$\displaystyle a_{4}=\frac{\beta^{4}}{12}\Bigg{(}\frac{(-4+e^{\beta\Delta})(p-1)}{(2+e^{\beta\Delta})^{2}}+\frac{pe^{\beta\Delta}(-1+4e^{\beta\Delta})}{(1+2e^{\beta\Delta})^{2}}\Bigg{)}$$

$$\displaystyle a_{6}=\frac{\beta^{6}}{360}\Bigg{(}\frac{(64-26e^{\beta\Delta}+e^{2\beta\Delta})(p-1)}{(2+e^{\beta\Delta})^{3}}$$

$$\displaystyle\;\;\;\;\;-\frac{pe^{\beta\Delta}(1-26e^{\beta\Delta}+64e^{2\beta\Delta})}{(1+2e^{\beta\Delta})^{3}}\Bigg{)}$$

$$\displaystyle a_{8}=\frac{\beta^{8}}{20160}\Bigg{(}\frac{(1188e^{\beta\Delta}-2176-120e^{2\beta\Delta}+e^{3\beta\Delta})(p-1)}{(2+e^{\beta\Delta})^{4}}$$

$$\displaystyle\;\;\;\;\;+\frac{pe^{\beta\Delta}(-1+120e^{\beta\Delta}-1188e^{2\beta\Delta}+2176e^{3\beta\Delta})}{(1+2e^{\beta\Delta})^{4}}\Bigg{)}$$

$$1+\frac{2\beta (p-1)}{2+e^{\beta \Delta }}=\frac{2\beta p{e}^{\beta \Delta }}{1+2e^{\beta \Delta }}$$

$$\frac{a_6^2}{4a_4{a}_8}=1$$

Full text

I Introduction

Multicritical points typically occur in systems described by three or more thermodynamic fields. In these systems, there can be critical points that can be reached only by fixing three or more thermodynamic parameters. Hence the full phase diagram of such systems is multi-dimensional knobler ; fisher ; griffithsh . Such critical points are ubiquitous in nature, in systems like binary fluids widom ; binarypd , metamagnets metamagnet , alloys of magnetic and non magnetic materials alloys , $He^{3}-He^{4}$ mixtures helium , quantum metals kirkpatrick , polymer collapse degennes and quantum chromodynamics qcd . Among them the tricritical point (TCP) is one of the most widely studied and well understood multicritical points lawrie . Solvable models which display higher order critical points are useful in outlining the topology of the phase diagrams barbosa . In this context, the mean field Blume-Capel modelblume ; capel has been very useful and is one of the most well studied models. It is the simplest model to exhibit a TCP. TCP is an example of a multicritical point, which is a point of confluence of three critical lines in the space of three fields $(T,\Delta,H)$. Here $T$ and $H$ are the temperature and external field respectively and $\Delta$ is a non-ordering field, known as the crystal field griffiths ; lawrie . In the $(T,\Delta)$ plane (with $H=0$), TCP shows itself as a point where the critical line ends in a first order line.

Introducing randomness in bond strength or field strength is known to affect the phase diagram. For example, in two dimensions it was shown that even an infinitesimal amount of random field disorder can change a first order transition to a continuous transition or can destroy it altogether aizenman ; huiberker . In dimensions higher than two, similar behaviour is expected for strong disorder cardyj .

Blume-Capel model was introduced as an extension of Ising model where in addition to $\pm 1$, spin at each site can also be [math]. It has an extra term to take care of crystal field anisotropy blume ; capel . It was originally used to explain the first order magnetic transition in materials like $UO_{2}$ frazer . Since then the model has found application in explaining the behaviour of wide range of physical systems like $He_{3}-He_{4}$ mixtures helium , ternary fluidstfluids , semiconductor alloyssalloys , phenomena of inverse meltingimelting and so on. Study of these systems in random media is modelled by Blume Capel model with random crystal field buzano . The phase diagram of the model is known to change under the effect of disorder. In particular, mean-field random crystal field Blume-Capel model has been studied extensively using many different techniques like mean field approximation mfdilution ,effective field theoriesyuksel , renormalisation groupbranco ; snowman , Bethe latticebethe ,pair approximation method lara ,replica salmon ; santos and large deviations sumedha1 . All these work have focussed on the $(T,\Delta)$ plane. These different methods do not agree with each fully in the prediction of the phase-diagram, but they all report that the first order line and hence the TCP disappears for higher strengths of disorder. Some of them mfdilution ; salmon ; santos predicted different topologies of the phase-diagram depending on the strength of disorder, with multicritical points like critical end points, ordered critical point and double critical point.

Since TCP is a point of confluence of three critical lines, it is important to look at the effect of disorder on the other two lines meeting at the TCP. Hence, we revisit the problem and obtain the phase diagram in the space of three fields, on a fully connected graph by solving the model exactly, using large deviation theory touchette . For the Blume Capel model, the mean field solution is known to give the correct values of the critical exponents at the TCP in three dimensions lawrie and predicts the correct topology of the phase diagram for $d\geq 2$, in the pure casezirenberg ; silva . Similarly, we expect the change in topology of the phase diagram as a function of disorder strength derived in this paper to be robust and not restricted to mean-field solutions. Note that Landau approach cannot be used when the external field $H\neq 0$, as the value of the magnetisation is finite along the critical lines. We hence have to make use of the full free energy functional even to determine the critical lines.

We find that the TCP persists for only very weak disorder strengths. As the disorder strength increases, the TCP vanishes and a different multicritical point, bicritical end point(BEP) emerges where only two of the three critical lines end on a first order surface nelson ; helena-barbosa ; plascak . This point has been wrongly reported as an ordered critical point in earlier studies in $(T,\Delta)$ plane santos . BEP has been comparatively less observed and studied in the literature. Two well known examples where BEP has been observed are: anisoptropic continuous spin systems as an end point of the spin flop linenelson and in spin 3/2 systems with crystal fieldplascak ; butera . We also find that the model exhibits critical end points(CEP) for intermediate strengths of disorder as reported in earlier studies santos . Critical end point is a critical point where a line of second order transitions terminates at a line of first order transitions chaikin . Alternately, it can also be defined as a point where two phases become critical in the presence of one or more ordered phases, known as the spectator phases upton , in systems with multiple phases. We thus find three different phase diagrams depending on the strength of disorder. Recall that for the pure Blume Capel model, the phase diagram has three critical lines ($\lambda,\lambda_{+},\lambda_{-}$) which all meet at the tricritical point lawrie . In the mean field limit, all the critical points along these three lines fall in the Ising universality class. Along $\lambda$-line there is a spontaneous symmetry breaking transition from state with magnetisation $m=0$ to a state with $|m|\neq 0$ in $(T,\Delta)$ plane. Switching on the magnetic field introduces bias toward $m_{+}$ or $m_{-}$ state, depending on the sign of the magnetic field. This results in the coexistence of $m=0$ state with $m_{+}(m_{-})$ states for low $H_{+}(H_{-})$ respectively. These two coexistence surfaces meet along a triple line in the $(T,\Delta)$ plane. $\lambda_{+}/\lambda_{-}$ lines separate this coexistence surface from the ferromagnetically ordered phases with opposite magnetisation. Hence the pure Blume Capel phase diagram has only one ferromagnetic state and one paramagnetic state in the $(T,\Delta)$ plane.

We find that for weak disorder, the three critical lines meet at a tricritical point(see Fig1(a)). But a new ferromagnetic state appear now at very low temperature for all $p>0$, which is separated from usual ferromagnetic state via first order quadruple line, that ends in a CEP. The value of magnetisation in this phase depends on the strength of disorder and increases with increasing disorder strength.

For intermediate disorder strengths, the two critical lines with $H\neq 0$ meet at a BEP and the line of continuous transition in $(T,\Delta)$ plane (known as $\lambda$ line) meets a line of first order transition at a CEP. CEP and BEP are connected via a quadruple line, along which the four phases co-exist(see Fig 1(b)). The quadruple line separates two ferromagnetic phases. This happens due to abrupt change in the number of $\pm 1$ spins across this line. Unlike the new ferromagnetic state that occurs at low temperature, entropy is important for this ferromagnetic state as it occurs at a relatively higher temperature and hence entropically it becomes useful to have more spin particles. The system hence has three ferromagnetic and one paramagnetic phase in this range of disorder.

For strong disorder, the BEP persists but CEP vanishes and the $\lambda$ line continues to $\Delta\rightarrow\infty$ (see Fig. 1(c)). This phase has the same two ferromagnetic phases as the weak disorder case and a paramagnetic phase. We will study these three topologies in this paper.

The plan of the paper is as follows: In Section II we discuss the Blume-Capel model in the presence of external field and derive the equations for critical lines in $(T,\Delta,H)$ space. In Section III we study the phase diagram for strong disorder and intermediate disorder by using the full free energy functional. We also look at the magnetisation, density and magnetic susceptibility near the BEP, to understand the nature of BEP. In Section IV we briefly discuss the case of weak disorder and in Section V we show that a Landau expansion of the free energy functional cannot describe the BEPs and CEPs of this model. We conclude with a short discussion in Section VI.

II Model

We study the Blume-Capel model with random crystal field disorder in the presence of external field on a fully connected graph. The Hamiltonian can be written as

[TABLE]

where $\Delta_{i}$ represent quenched random crystal field at each site, $H$ is the external field and $s_{i}$ are spin$-1$ random variables which can take $\pm 1,0$ values. There are two order parameters: magnetisation, $m=\bar{s}$ and density of $\pm 1$ spins, $q=\bar{s^{2}}$. These are obtained by taking a quenched average of the random variables $s$ and $s^{2}$ respectively cardy . We draw random crystal fields from bimodal distribution of the kind:

[TABLE]

Since $p=0$ or $1$ will imply no disorder and $p=1/2$ would be the most random case, it is enough to look for $0\leq p\leq 0.5$.

It can be shown that the probability of a configuration $C_{N}$ satisfies large deviation principle(LDP) in the presence of random crystal field disordertouchette ; ellis ; sumedha1 ,i.e

[TABLE]

The rate function $I(x_{1},x_{2})$ for bimodal random crystal field disorder in the absence of external field was calculated using tilted LDP recentlysumedha1 . Using the same method, the rate function in the presence of external field is:

[TABLE]

where $z$ is the solution of the equation:

[TABLE]

In the limit, $N\rightarrow\infty$, for a given $\beta,\Delta$ and $H$, the value of $x_{1}$ and $x_{2}$ that minimise $I(x_{1},x_{2})$ will give the value of magnetisation ($m$) and density($q$). The minima of the rate function in $(x_{1},x_{2})$ plane gives the free energy for a given $\beta(=1/T),\Delta$ and $H$. Hence the values of $x_{1}$ and $x_{2}$ which minimise $I(x_{1},x_{2})$ are the value of $m$ and $q$ respectively for a given set of thermodynamic variables. Minimising $I(x_{1},x_{2})$ with respect to $x_{1}$ and $x_{2}$ results in the following equations for $m$ and $q$:

[TABLE]

[TABLE]

where $z$ is related to $q$ via the Eq. 5, i.e:

[TABLE]

In subsection II.1 we will recap the results in the absence of external field and then build the equations for phase diagram in $(T,\Delta,H)$ space in subsection II.2

II.1 Two field phase diagram in the $(T,\Delta)$ plane

For $H=0$, the phase diagram has been studied earlier sumedha1 ; santos . We will briefly recap those results here: Assuming $m$ to be small the fixed point equations, Eq. 6 and 7 can be linearized around $m=0$. This gives $q=1/\beta$ and $z=2$ at the critical point. Substituting these values in Eq. 8, gives the equation for a line of continuous transition in the $H=0$ plane. The line of continuous transition in $H=0$ plane is known as the $\lambda$-line and satisfies the following equation

[TABLE]

This is valid only when the higher order terms in the expansion can be ignored. Taking $q=(1+\epsilon)/\beta$ and expanding in powers of $\epsilon$ we find that the coefficient of linear term in $\epsilon$ becomes zero when

[TABLE]

Solving Eq. 9 and 10 together we get the condition for break down of linear approximation as

[TABLE]

Hence for a given $\Delta$, there will be either no transition or a first order transition, beyond the value of $\beta$ that satisfy Eq. 11

The value of $(\beta,\Delta)$(or equivalently $(T,\Delta)$) which satisfy Eqs. 9 and Eq. 11 simultaneously gives the location of TCP for a given $p$. It was found in sumedha1 that beyond $p_{c}=0.0454$ the two equations cannot be satisfied simultaneously and hence there is no TCP, and the $\lambda$ line in the $(T,\Delta)$ plane extends to $\Delta\rightarrow\infty$. This treatment is equivalent to Taylor expanding the rate function to get an equivalent Landau free energy functional, which we will discuss in Sec. V.

II.2 Three field phase diagram in $(T,\Delta,H)$ space

Let us now take $H\neq 0$ and look for the critical points in the full ($T,\Delta,H$) space. We know that at the TCP there are two other continuous lines with $H\neq 0$ which meet the $\lambda$ line. We call these, depending on the value of $H$, as $\lambda_{+}$ and $\lambda_{-}$.

We wish to understand the effect of disorder on the two critical lines $\lambda_{+}$ and $\lambda_{-}$. We will focus on the effect of disorder on these two critical lines in this paper. Along these lines, $m\neq 0$ and one cannot look for continuous transition by expanding the free energy functional like we did in the section II.1.

Note that at the fixed point the value of $m$ and $q$ are related via Eq. 6. Since we are only interested in the fixed points, at fixed points $q$ is completely determined by $m$, hence the rate function can be replaced by a one parameter functional $\tilde{f}(m)$, which comes out to be:

[TABLE]

From this we get the following self-consistent equation for $m$:

[TABLE]

Since $m\neq 0$ along the $\lambda_{+}$ and $\lambda_{-}$ lines, expanding $\tilde{f}(m)$ in powers of $m$ to get a Landau free energy functional will not give the correct critical behaviour. But, in general along a critical line, the first three derivatives of the free energy functional with respect to the order parameter should be zero. This is because between two successive minimas , there must exist two points of inflexion, i.e $f^{{}^{\prime\prime}}=0$ and hence also a point where $f^{{}^{\prime\prime\prime}}=0$. Hence at the continuous transition, all three derivatives should vanish simultaneously. Hence to study $\lambda_{+}$ and $\lambda_{-}$ critical lines we equate the first three derivatives of $\tilde{f}(m)$ w.r.t $m$ to zero lawrie (and fourth derivative should be greater than zero). This is true also for the $\lambda$ line, as for $H=0$ and $m=0$ the third derivative is trivially zero and second derivative gives the same condition as Eq. 9.

In general, equating second and third derivative of $\tilde{f}(m)$ to zero we get the following two conditions respectively:

[TABLE]

here $x=\exp(\beta\Delta)$ and $y=\cosh\beta(m+H)$. For $p\neq 0$, the two equations are quartic and hexic in $x$.

For $p=0$, they reduce to the following simpler equations :

[TABLE]

Solving these equations we get

[TABLE]

Hence, we reproduce the classic results of Blume,Emery and Griffiths beg : There is a line of critical points for $4\geq\beta\geq 3$ for $H>0$ and another for $H<0$. Both critical lines extend to $\Delta\rightarrow\infty$. These two lines enclose two first order surfaces which meet in the $H=0$ plane along a triple line(line with three phase co-existence). Above $\beta=4$ there is no value of $x$ and $y$ that can satisfy Eqs. 16 and 17 simultaneously. The magnetisation along these two critical lines is not zero and is equal to

[TABLE]

This can be used to get the value of $H$ along the critical lines, which comes out to be

[TABLE]

These two critical lines meet in the $H=0$ plane at a point with $T_{TCP}=1/3$ and $\bigtriangleup_{TCP}=0.462098$. This is the well known TCP in ($T,\Delta$) plane for $p=0$(can be obtained by solving Eq. 9 and 11 simultaneously for $p=0$).

For $p\neq 0$, we use Mathematica mathematica to solve Eq(14) and Eq(15) simultaneously to get the two critical lines numerically. To solve the equations for any arbitrary $p$, we scan different values of $\beta$ and $\Delta$ and hence $x$ and solve Eq(14) (corresponding to $\tilde{f}^{\prime\prime}(m)=0$) exactly to get the corresponding value of $y$. Then we substitute the value of $x$ and $y$ in Eq(15) to check if $(x,y)$ satisfy the condition, $\tilde{f}^{\prime\prime\prime}(m)=0$.

For each set of $(x,y)$ that satisfy Eq. 14 and Eq. 15 simultaneously, we can calculate $m$ using the equation:

[TABLE]

The above equation is derived from Eq. 13 by taking $\cosh\beta(m+H)=y$ and $\exp(\beta\Delta)=x$. The corresponding value of $H$ along the critical lines can then be calculated by inverting $y=\cosh\beta(m+H)$.

For a TCP to exist the two critical lines in the $H\neq 0$ plane should meet in $H=0$ plane at the point where second order line ends in a first order transition line in the ($T,\Delta$) plane. We can put $H=0$ and $m=0$ in Eqs. 14 and 15 to directly look for this point. Hence, we separately solve the two equations for $y=1$. Interestingly, we find that for $y=1$, the two equations can be solved simulatenously only for $p\leq p_{c}(=0.0454)$. This is also the value of $p$ beyond which linear stability analysis breaks down and Eq. 9 is not valid anymore. More interestingly even though the two equations can be solved for $H=0$ till $p\leq 0.0454$, we find that for $p>0.022$, one more solution shows up, with $m\neq 0$ and $H=0$. For $p>0.0454$, all possible solutions have $m\neq 0$.

We find that the two critical lines, $\lambda_{+}$ and $\lambda_{-}$ meet $\lambda$ line at a TCP for $p<0.022$. For $p>0.022$, the two critical lines, $\lambda_{+}$ and $\lambda_{-}$ meet inside the first order surface, i.e at a point where $m\neq 0$. This point hence is not a TCP, but a BEP. Furthermore, we find that for $p>0.022$ there are two different kinds of phase diagrams possible: For $0.022<p\leq 0.1078$ the phase diagram is as shown in Fig. 1(b): In $H=0$ plane there is a four phase coexistence line starting from the BEP which separates the two ordered phases. This line meets the $\lambda$-line defined via Eq. 9, giving rise to a CEP. From CEP there is a three phase coexistence line which ends in another CEP. For $0.1078<p\leq 0.5$ the phase diagram is as shown in Fig. 1(c): There is a four phase coexistence line from BEP which never crosses the $\lambda$ line defined via Eq 9 and goes all the way to $T=0$. Moreover, we find that $\lambda_{+}$ and $\lambda_{-}$ critical lines exist for all strengths of disorder(i.e for all values of $p$). We give more details of these multicritical points and phase topologies in the next few sections.

III Strong disorder and BEP

For $p>0.022$, the two critical lines for $H\neq 0$ do not meet at the potential TCP point as given by simultaneous solution of Eq 9 and 11. Instead they meet inside the ordered plane. We find that the two wings are separated by a first order line in $H=0$ plane, which behaves differently for $0.022<p\leq 0.1078$ and for $0.1078<p\leq 0.5$. Hence we will look at these two regimes separately.

III.1 $0.1078<p\leq 0.5$

For this range of $p$, along the first order line in $H=0$ plane there is a four phase coexistence, which ends in a bicritical end point(see Fig1(c)). This line is a first order transition line between two ordered states with different values of magnetisations. These two different ordered states are a result of disorder and are not present in the pure system. At low temperatures, the system prefers $\pm 1$ spin states when $\Delta$ is small. As $\Delta$ increases, due to disorder, states with finite fraction of zero spins compete with the states with only $\pm 1$ spins. This can be seen by looking at the order parameter $q$ as a function of $\Delta$, as shown in Fig. 2.

One can see all the transitions clearly by plotting $\tilde{f}(m)$ in different regions of the phase diagram as shown in Fig. 3 for $p=0.2$. From the plots we can see that the $H=0$ line separates the two ordered phases. Along $H\neq 0$ critical lines, two of these phases become critical and at BEP the two critical phases coexist.

To understand the nature of transition especially at BEP, we looked at the magnetisation$(m)$ and magnetic susceptibility $(\chi)=\frac{\partial m}{\partial H}|_{H\rightarrow 0}$. Let us first look at the magnetisation as a function of $T$ in the $H=0$ plane for different fixed values of $\Delta$ (see Fig. 4). We find that for $\Delta<\Delta_{BEP}$, the magnetisation changes its slope near $T=T_{BEP}$, the change becomes sharper as one approaches $\Delta=\Delta_{BEP}$. For $\Delta>\Delta_{BEP}$(but close to $\Delta_{BEP}$), the magnetisation undergoes a first order transition as it crosses the quadruple line and then changes slope near $T=T_{BEP}$. For $\Delta$ much larger than $\Delta_{BEP}$, as we increase $T$ there is no first order jump or change of slope around $T=T_{BEP}$. We also looked at $m$ as a function of $\Delta$ for three different values of $T$ (see fig. 5). First order jump as one crosses the quadruple line is clear for $T<T_{BEP}$. For $T>T_{BEP}$ there is no signature of any transition.

It is hard to deduce the nature of transition at BEP by looking at the magnetisation alone. Hence we studied the magnetic susceptibility near BEP. First we look at it for fixed value of $\Delta$. As we fix $\Delta=\Delta_{BEP}$ and vary $T$, we find that there is an infinite peak at the $T$ of $\lambda$ transition. There is another peak at $T=T_{BEP}$, but this peak is finite(see Fig. 6). This behaviour can be contrasted with the behaviour at $\Delta>\Delta_{BEP}$ as shown in Fig 7. We find a discontinuity where it crosses the first order line and a finite peak near $T=T_{BEP}$.

We also studied magnetic susceptibility as we vary $\Delta$ at $T_{BEP}$. As expected, we found that there is a finite peak at $\Delta=\Delta_{BEP}$ For $T<T_{BEP}$, there was instead a first order jump in magnetic susceptibility. We scanned a large region in $(T,\Delta)$ plane near BEP. We find that the effect of the presence of BEP is felt even far away from the point. But the magnetisation and susceptibility plots are smooth near BEP(though susceptibility shows a cusp). It was shown via scaling arguments helena-barbosa that if the two critical lines meeting at BEP are in the same universality class and are symmetric, then the singular behaviour contribution to the phase boundary cancels out helena-barbosa ; plascak . In our case the two critical lines $\lambda_{+}$ and $\lambda_{-}$ lie in the Ising universality class. Looking at the three dimensional phase diagram it is clear that there is only one phase in the system in the sense that there exist a path between any two non-singular points in the phase diagram which does not have to encounter a singularity. At BEP the first three derivatives of $\tilde{f}(m)$ w.r.t $m$ are zero and hence the free energy is not analytic at this point. Hence, we conclude that BEP is a point of two phase co-existence and there is no critical transition from one phase to another at BEP.

As $p$ increases we find that the critical lines enclosing the wings become flatter and the temperature at which they meet in $H=0$ plane decreases. We have tabulated the range of $T$ for different $p$ in Table 1.

III.2 $0.022<p\leq 0.1078$

In this region the wings meet at BEP as before, but the first order quadruple line now intersects the $\lambda$-line at a critical end point(we will call this critical end point as CEP1 to distinguish it from the other critical end point in the phase diagram at a lower temperature, which we will call as CEP2). After that it becomes a line of triple point(see Fig1(b)). In Fig.8, we plot the free energy functional along this line. Along the first order line there is a line of four phase coexistence between BEP and CEP1 and then there is a usual triple line between CEP1 and CEP2. As shown in Fig. 8(3), CEP1 itself is neither a quaduple or a triple point. It is instead a point where a critical state coexists with two ordinary stable phases. Between CEP2 and [math] temperature there is again a quadruple line as shown in Fig 8(5).

The CEP is a point where two phases become critical in the presence of one or more non critical spectator phase. At CEP, $\tilde{f}(m)$ for $m=0$ and for $m\neq 0$ should be equal (i.e $\tilde{f}(m=0)=\tilde{f}(m\neq 0)$) along with their derivative with respect to $m$ ($\tilde{f}^{\prime}(m=0)=\tilde{f}^{\prime}(m\neq 0)$). If this point lies on the $\lambda$ line, then we get the condition for CEP. Hence to find CEP, we explore the $\lambda$-line for a point where $\tilde{f}(m=0)=\tilde{f}(m\neq 0)$ along with $\tilde{f}^{\prime}(m=0)=\tilde{f}^{\prime}(m\neq 0)$. We find that for $p>0.1078$ the condition cannot be satisfied.

In santos , Santos et al also reported the presence of CEP1 and CEP2 for $0.022<p<0.074$, by looking at the point of intersection of the $\lambda$-line with the first order line. We find that this topology extends till $p=0.1078$. In order to understand the discrepancy, we have plotted $\lambda$-line given by Eq. 9 along with a line parallel to $T$-axis at $\Delta=(1+p)/2$ in Fig. 9. The line $\Delta=(1+p)/2$ is a good approximation to the first order line in the $(T,\Delta)$ plane as we found that the first order line is almost parallel to $T$-axis. As shown in Fig. 9, $\Delta=(1+p)/2$ line crosses the $\lambda$-line once till $p\approx 0.07$ and thrice for $0.07<p<0.11$. For $p>=0.11$ there is no intersection. More careful analysis using the full free energy functional, gives us the value to be around $p=0.1078$. This matches with the value obtained by equating the free energy functional and its first derivative along the $\lambda$-line, as described in the previous paragraph.

In Table 2 we tabulate the location of BEP, CEP1 and CEP2 for different values of $p$. The first order line between CEP2 and $\Delta$- axis is similar to the first order line reported in Sec III.1, which separates the states with almost all $\pm 1$ spins from a state with $p$ fraction of $\pm 1$ spins. The presence of CEP1 and a four phase co-existence line between BEP and CEP1 is due to the occurence of a new magnetic state. This state has more than $p$ fraction of $\pm 1$ spins, as it occurs at a higher temperature, very close to the $\lambda$ line.

We plot the magnetic susceptibility as a function of $\Delta$ for $T=T_{BEP}$ and for $T=T_{CEP}$ in $H=0$ plane. As expected first plot shows two peaks:a finite peak at $BEP$ and an infinite peak at intersection with the $\lambda$-line (see Fig. 10), while the second plot shows one peak only at CEP1 (see Fig. 11).

IV Weak disorder and TCP

Along the region $0\leq p\leq 0.022$ the wings meet the $\lambda$-line at the TCP and the phase diagram is similar to the pure case. Along the first order line there is three phase coexistence. As $p$ increases, the TCP shifts towards smaller $T$ and larger $\Delta$. At $p=0.022$ the TCP becomes a fourth order critical point. Again, at very low temperature there is a CEP, similar to the case discussed in the Section III.2 for all $p>0$, which is a state with $p$ fraction of $\pm 1$ spins and occurs at low temperatures on a complete graph. For $0.017<p<0.022$, there is re-entrance region in the phase diagram, as the TCP does not coincide with the maximum of $\lambda$-line given by Eq. 9.

V Landau theory

In the previous sections we studied the phase diagram by looking at the full free energy functional and its derivatives. Usually, Landau theory is a very useful tool to classify different kinds of transitions and even though it might not be accurate quantitatively, it helps in understanding different possible topologies of the phase diagram. But while very successful in explaining ordinary critical point, it is not always possible to find a Landau description for higher order critical point, i.e. it is perhaps possible to define a free functional always, but it might not always be Taylor expandable deepak . In this section we expand the free energy functional to check if we can explain the phase diagrams based on the coefficients of different powers of the order parameter. For example, the Ising universality class critical point can be determined easily by expanding upto fourth power in $m$, provided that the next higher order coefficient is positive. For TCP one needs to expand till sixth order. A sixth order Landau theory hence allows only for ordinary critical points and TCPs. We expect that we need to keep more terms in the expansion, if we expect to find higher order critical points like CEP and BEP japan . Hence we expanded the free energy functional till eighth power of $m$. We get

[TABLE]

where $a_{i}$’s are Landau coefficients, as follows:

[TABLE]

[TABLE]

The second order transition is given by $a_{2}=0$, provided $a_{4}>0$. Equating $a_{2}=0$ gives us:

[TABLE]

This equation is same as Eq. 9, obtained by linear expansion around $m=0$. According to the Landau theory, a new universality class, namely the TCP occurs when $a_{4}$ becomes equal to [math], provided $a_{6}>0$. We find that the condition for $a_{4}=0$ along the $\lambda$-line is the same as given by substituting Eq. 11 into Eq. 9. For $p>p_{c}=0.0454$, $a_{4}$ is never [math] and hence beyond $p_{c}$ the condition for occurence of TCP cannot be satisfied. For $p>0.022$, $a_{6}<0$ at the point where $a_{4}=0$. Hence sixth order Landau theory while sufficient for $p<0.022$, is not enough for $p>0.022$.

Hence for $a_{6}<0$, we consider the expansion till eighth order, since $a_{8}>0$ for all ranges of the parameters. CEP will be a point along the $\lambda$-line (given by Eq. 25) where the $\tilde{f}(T_{c},m_{c})=0$ and $\tilde{f}^{\prime}(T_{c},m_{c})=0$ and $m_{c}\neq 0$. Solving these, we get the condition for the existence of CEP to be

[TABLE]

We find that Eq. 26 can be satisfied only for $0.022<p\leq 0.0454$, and that too at a point very close to the point where $a_{4}=0$. For example, for $p=0.044$ from Eq. 26, we get $(T_{CEP1},\Delta_{CEP1})=(0.267,0.497)$ and for $p=0.03$ we get $(T_{CEP1},\Delta_{CEP1})=(0.294,0.489)$. Hence we find that the value obtained via Eq. 26 are different from the ones obtained by looking at the full free energy functional in Section III.2(see Table II). The difference increases with increasing $p$. More importantly, in Section III.2 we had found numerically that CEP is present for a much larger range of $p$: $0.022<p\leq 0.1078$.

To estimate BEP using truncated $\tilde{f}(m)$, we equate the first three derivatives of the truncated $\tilde{f}(m)$ in Eq. 23 w.r.t $m$ to [math]. For $m\neq 0$, this gives the condition for BEP to be: $a_{6}=-\sqrt{\frac{8a_{4}a_{8}}{3}}$. Again this condition gets satisfied only for $0.022<p\leq 0.0454$. This gives a BEP very close to CEP and the actual location does not match with the numerical estimates of Section III. Hence, a Landau description of this system predicts the phase diagram correctly for $p<0.022$(except for CEP present at very low temperatures for all $p>0$) and gives qualitatively similar diagram for $0.022<p<0.0454$, though the location of BEP and CEPs does not match the actual value. For $p>0.0454$ it is inadequate in predicting the phase diagram. We tried including more terms in the expansion of $\tilde{f}(m)$, but we could not locate BEP using a truncated $\tilde{f}(m)$, suggesting that full $\tilde{f}(m)$ is needed for locating the BEP.

VI Discussion

Blume-Capel model is a very useful model due to its simplicity and rich phase diagram. Its phase diagram in the presence of disorder in $(T,\Delta)$ plane has been studied extensively using many different techniques. In this paper we studied the three field phase diagram in the presence of disorder, which has not been studied earlier. This is useful especially to correctly predict the nature of multicritical points. We found as the disorder strength increases, the two wings meet at a BEP. We showed that this is actually a point of co-existence of two critical phases, where the magnetic susceptibility is finite. Hence in-spite of the three derivatives of the free energy being zero at BEP, it is not critical. This point was identified as an ordered critical point in earlier studies santos ; salmon . Also, we corrected the estimate of onset of topology III as a function of disorder strength compared to earlier estimate santos . In mfdilution a different bimodal distribution of the random crystal field ($P(\Delta_{i})=p\delta(\Delta_{i}-\Delta)+(1-p)\delta(\Delta)$) was studied using Landau theory. They observed that as the disorder anisotropy increases, the first order line meets the second order line such that there is a re-entrant part in the phase diagram. By analogy with binary fluids, they conjectured that when the second order line has a re-entrant part, it will end in a double critical point(or bicritical end point) followed by a CEP as in Fig 1(b). Since their work was based on Landau expansion, they could not identify the bicritical end point and critical end point precisely. Also they did not report a topology similar to Fig. 1(c) for strong disorder.

It would be interesting to see if similar phase diagrams are realized in finite dimensions using numerical simulations fytas ; deserno ; sethna . The origin of BEP in our model is different than that for the pure anisotropic continuous spin systems, where the BEP was seen as an end-point of spin flop transition line nelson . A two parameter Landau theory description exists for spin-flop chaikin . Our free energy functional is a one parameter function which shows BEP. It would be interesting to see if a one parameter Landau theory can be built, which has a BEP as seen in the strong disorder case in our work. A study of three fields diagram for random field Blume Capel model would also be useful to understand the nature of TCPs reported in rf in the $(T,\Delta)$ phase diagram wip .

VII Acknowledgement

We thank Deepak Dhar for valuable comments on the manuscript and Mustansir Barma for discussions. S thanks ICTP for hospitality during the completion of this work.