text:
# 1+∞ Dimensional Attractor Neural Networks

## 1 Introduction

Solvable models of recurrent neural networks are bound to be simplified representations of biological reality. The early statistical mechanical studies of such networks, e.g. , concerned mean-field models, whose statics and dynamics are by now well understood, and have obtained the status of textbook material . The focus in theoretical research has consequently turned to new areas, such as solving the dynamics of large recurrent networks close to saturation , the analysis of finite size phenomenology , solving biologically more realistic models , or networks with spatial structure . In this paper we analyse Ising spin models of recurrent networks with spatial structure, in which there are two coexistent classes of Hopfield-type  interactions: infinite-range ones (operating between any pair of neurons), and 1D short-range ones (operating between nearest neighbours only). The study of this type of structure is motivated by the interplay between long-range processing via (excitatory) pyramidal neurons and short-range processing via (inhibitory) inter-neurons, which is typically observed in cortical tissue. In the present model however, and in contrast to early papers on spatially structured networks, exact solutions based solely on simple mean-field approaches are ruled out. Due to short-range interactions analytical solutions require significantly more complicated methods and the present models can be solved exactly only by a combination of mean- and random-field techniques , whereas for the special case in which the system has stored a single pattern, a simple (Mattis) transformation allows us to derive the solution via a combination of mean-field methods and transfer-matrices.

Our paper is organised as follows. We first solve the one-pattern case in which pattern-variables can be transformed away (thus providing a convenient and exactly solvable benchmark case against which to test the general theory). This also hints at the interesting features induced by short- versus long-range competition in the more general model. In particular, already in the one-pattern model we find first-order phase transitions, regimes corresponding to multiple locally stable states and we find that sequential and parallel dynamics phase diagrams are related by simple transformations. We then proceed to the general case, with an arbitrary number of stored patterns, away from saturation regimes. For sequential dynamics (Boltzmann equilibrium) we adapt the 1D random-field techniques as originally developed for site-disordered Ising chains, and combine them with mean-field methods; for parallel dynamics (Peretto equilibrium) we adapt and combine with mean field methods the procedure in , based on $`4\times 4`$ random transfer matrices. The disorder-averaged free energy and the order parameters are found, in closed and exact form, as integrals over the distribution of a characteristic variable, which represents a specific ratio of conditioned partition functions. This distribution is calculated following .

In the region where the infinite-range versus short-range competition is prominent, our theory predicts a *series* of (continuous and discontinuous) dynamic phase transitions, and a free energy surface with multiple local minima. These features become more prominent, in both number and strength, when the number of stored patterns increases, in sharp contrast to Hopfield-type infinite-range networks , where ‘pure state’ solutions are independent of the number of patterns stored. The transition lines of sequential and parallel dynamics phase diagrams are found to be related by reflection symmetries and also parallel dynamics macroscopic equations can describe 2-cycles rather than fixed-points solutions. Finally we test our theory against several exactly solvable cases and against numerical simulations; the latter are found to exhibit interesting but extremely slow domain-induced dynamics, with plateau phases. Once equilibration has occurred, we obtain excellent agreement between theory and experiment.

## 2 Model Definitions

We study models with $`N`$ Ising spin neuron variables $`𝝈=(\sigma _1,\mathrm{},\sigma _N)\{1,1\}^N`$, which evolve in time stochastically on the basis of post-synaptic potentials $`h_i(𝝈)`$ (or local fields), following the Glauber-type rule

$$\mathrm{Prob}[\sigma _i(t+1)=\pm 1]=\frac{1}{2}\left[1\pm \mathrm{tanh}[\beta h_i(𝝈(t))]\right]h_i(𝝈)=\underset{ji}{}J_{ij}\sigma _j+\theta _i$$
(1)
The parameters $`J_{ij}`$ and $`\theta _i`$ represent synaptic interactions and firing thresholds. The (non-negative) parameter $`\beta `$ controls the amount of noise, with $`\beta =0`$ and $`\beta =\mathrm{}`$ corresponding to purely random and purely deterministic response, respectively. If the interaction matrix is symmetric, both a random sequential execution and a fully parallel execution of the stochastic dynamics (1) will evolve to a unique equilibrium state. The corresponding microscopic state probabilities can both formally be written in the Boltzmann form $`p_{\mathrm{}}(𝝈)\mathrm{exp}[\beta H(𝝈)]`$, with associated Hamiltonians  (since Peretto’s pseudo-Hamiltonian Hamiltonian $`H_{\mathrm{par}}`$ depends on $`\beta `$, the associated statistics are not of the Boltzmann form):

$`H_{\mathrm{seq}}(𝝈)`$ $`=`$ $`{\displaystyle \frac{1}{2}}{\displaystyle \underset{ij}{}}\sigma _iJ_{ij}\sigma _j{\displaystyle \underset{i}{}}\theta _i\sigma _i`$ (2)
$`H_{\mathrm{par}}(𝝈)`$ $`=`$ $`{\displaystyle \frac{1}{\beta }}{\displaystyle \underset{i}{}}\mathrm{log}2\mathrm{cosh}[\beta h_i(𝝈)]{\displaystyle \underset{i}{}}\theta _i\sigma _i`$ (3)
In both cases, expectation values of order parameters can be obtained by differentiation of the free energy per neuron $`f=lim_N\mathrm{}(\beta N)^1\mathrm{log}_𝝈\mathrm{exp}[\beta H(𝝈)]`$, which acts as a generating function. For the parameters $`J_{ij}`$ and $`\theta _i`$ we now make the following choice 111Competition between a different type of uniform infinite-range and random nearest-neighbour interactions has been studied recently in , for sequential dynamics only.:

$$J_{ij}=\frac{J_{ij}^{\mathrm{}}}{N}+J_{ij}^s(\delta _{i,j+1}+\delta _{i,j1})\theta _i=\theta $$
(4)
$$J_{ij}^{\mathrm{}}=J_{\mathrm{}}^{(1)}+J_{\mathrm{}}^{(2)}\underset{\mu =1}{\overset{p}{}}\xi _i^\mu \xi _j^\mu J_{ij}^s=J_s^{(1)}+J_s^{(2)}\underset{\mu =1}{\overset{p}{}}\xi _i^\mu \xi _j^\mu $$

This corresponds to the result of having stored a set of binary patterns $`\{𝝃^1,\mathrm{},𝝃^p\}`$ with $`\xi _i^\mu \{1,1\}`$ in a one-dimensional chain of neurons, through Hebbian-type learning, but with different (potentially conflicting) embedding strengths $`J_s^{(2)}`$ and $`J_{\mathrm{}}^{(2)}`$ associated with the short-range versus the infinite-range interactions. We will choose $`pN`$. The parameters $`J_s^{(1)}`$ and $`J_{\mathrm{}}^{(1)}`$ control uniform contributions to the interactions within their class. Taking derivatives of $`f`$ with respect to $`J_{\mathrm{}}^{(2)}`$ produces expressions involving the familiar ‘overlap’ order parameters:

$`\mathrm{sequential}:`$ $`𝒎^2=2{\displaystyle \frac{f}{J_{\mathrm{}}^{(2)}}}=\underset{N\mathrm{}}{lim}({\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}\sigma _i𝝃_i)^2_{\mathrm{eq}}`$
$`\mathrm{parallel}:`$ $`𝒎^2={\displaystyle \frac{f}{J_{\mathrm{}}^{(2)}}}=\underset{N\mathrm{}}{lim}({\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}\sigma _i𝝃_i)^2_{\mathrm{eq}}`$
where the brackets $`\mathrm{}_{\mathrm{eq}}`$ denote equilibrium averages, $`𝒎=(m_1,\mathrm{},m_p)`$, and $`𝒙_i𝒙_j=_\mu x_i^\mu \xi _j^\mu `$. Note that for $`J_s^{(2)}=J_{\mathrm{}}^{(2)}=0`$ we obtain the simpler model

$$J_{ij}=\frac{J_{\mathrm{}}}{N}+J_s(\delta _{i,j+1}+\delta _{i,j1})\theta _i=\theta $$
(5)
The Mattis transformation $`\sigma _i\sigma _i\xi _i`$ maps this model onto

$$J_{ij}=\frac{J_{\mathrm{}}}{N}\xi _i\xi _j+J_s(\delta _{i,j+1}+\delta _{i,j1})\xi _i\xi _j\theta _i=\theta \xi _i$$
(6)
which corresponds to the result of having stored just a single pattern $`𝝃=(\xi _1,\mathrm{},\xi _N)\{1,1\}^N`$. Taking derivatives of $`f`$ with respect to the parameters $`\theta `$ and $`J_s`$ in (5) produces our order parameters:

$`\mathrm{seq}:`$ $`m={\displaystyle \frac{f}{\theta }}=\underset{N\mathrm{}}{lim}{\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}\sigma _i_{\mathrm{eq}}a={\displaystyle \frac{f}{J_s}}=\underset{N\mathrm{}}{lim}{\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}\sigma _{i+1}\sigma _i_{\mathrm{eq}}`$
$`\mathrm{par}:`$ $`m={\displaystyle \frac{1}{2}}{\displaystyle \frac{f}{\theta }}=\underset{N\mathrm{}}{lim}{\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}\sigma _i_{\mathrm{eq}}a={\displaystyle \frac{1}{2}}{\displaystyle \frac{f}{J_s}}=\underset{N\mathrm{}}{lim}{\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}\sigma _{i+1}\mathrm{tanh}[\beta h_i(𝝈)]_{\mathrm{eq}}`$
where we have simplified the parallel dynamics observables with the identities

$$\sigma _{i+1}\mathrm{tanh}[\beta h_i(𝝈)]_{\mathrm{eq}}=\sigma _{i1}\mathrm{tanh}[\beta h_i(𝝈)]_{\mathrm{eq}}\mathrm{and}\mathrm{tanh}[\beta h_i(𝝈)]_{\mathrm{eq}}=\sigma _i_{\mathrm{eq}}$$

which follow from (1) and from invariance under the transformation $`iN+1i`$ (for all $`i`$). For model (5) $`m`$ is the average neuronal activity. For sequential dynamics $`a`$ describes the average equilibrium state covariances of neighbouring neurons, and for parallel dynamics it gives the average equilibrium state covariances of neurons at a given time $`t`$, and their neighbours at time $`t+1`$ (the difference between the two meanings of $`a`$ will be important in the presence of 2-cycles). For model (6) one similarly finds

$$m=\underset{N\mathrm{}}{lim}\frac{1}{N}\underset{i}{}\xi _i\sigma _i_{\mathrm{eq}}a=\underset{N\mathrm{}}{lim}\frac{1}{N}\underset{i}{}(\xi _i\sigma _i)(\xi _{i+1}\sigma _{i+1})_{\mathrm{eq}}$$

The observable $`m`$ is here the familiar overlap order parameter of associative memory models , which measures the quality of pattern recall in equilibrium. Note that $`m,a[1,1]`$.

## 3 Solution and Phase Diagrams for $`p=1`$

Before we proceed to the solution of the general model (4) we first solve the relatively simple situation, where a single pattern has been stored following (5) and (6). This has the following advantages: it being a simpler version of (4), it allows us to explore general features and build intuition, without as yet any serious technical subtleties, as brought up by the general model. At the same time it provides an excellent benchmark test of the general theory to which it should reduce for $`J_s^{(2)}=J_{\mathrm{}}^{(2)}=0`$. For the remainder of this section our analysis will refer to model (5).

### Solution via Transfer Matrices

In calculating the asymptotic free energy per neuron $`f`$ it is advantageous to separate terms induced by the long-range synapses from those induced by the short-range ones, via insertion of $`1=𝑑m\delta [m\frac{1}{N}_i\sigma _i]`$. Upon using the integral representation of the $`\delta `$-function, we then arrive at

$$f=\underset{N\mathrm{}}{lim}\frac{1}{\beta N}\mathrm{log}𝑑m𝑑\widehat{m}e^{\beta N\varphi (m,\widehat{m})}$$

with

$`\varphi _{\mathrm{seq}}(m,\widehat{m})`$ $`=`$ $`im\widehat{m}m\theta {\displaystyle \frac{1}{2}}J_{\mathrm{}}m^2{\displaystyle \frac{1}{\beta N}}\mathrm{log}R_{\mathrm{seq}}(\widehat{m})`$
$`\varphi _{\mathrm{par}}(m,\widehat{m})`$ $`=`$ $`im\widehat{m}m\theta {\displaystyle \frac{1}{\beta N}}\mathrm{log}R_{\mathrm{par}}(m,\widehat{m})`$
The quantities $`R`$ contain all complexities due to the short-range interactions in the model. They correspond to

$`R_{\mathrm{seq}}(\widehat{m})`$ $`=`$ $`{\displaystyle \underset{𝝈\{1,1\}^N}{}}e^{i\beta \widehat{m}_i\sigma _i}e^{\beta J_s_i\sigma _i\sigma _{i+1}}`$
$`R_{\mathrm{par}}(m,\widehat{m})`$ $`=`$ $`{\displaystyle \underset{𝝈\{1,1\}^N}{}}e^{i\beta \widehat{m}_i\sigma _i}{\displaystyle \underset{i}{}}\mathrm{log}[2\mathrm{cosh}[\beta J_{\mathrm{}}m+\beta \theta +\beta J_s(\sigma _{i+1}+\sigma _{i1})]]`$
They can be calculated using the transfer-matrix method giving

$$R_{\mathrm{seq}}(\widehat{m})=\mathrm{Tr}[\text{T}_{\mathrm{seq}}^N]\text{T}_{\mathrm{seq}}=\left(\begin{array}{cc}e^{\beta J_si\beta \widehat{m}}& e^{\beta J_s}\\ e^{\beta J_s}& e^{\beta J_s+i\beta \widehat{m}}\end{array}\right)$$

$$R_{\mathrm{par}}(m,\widehat{m})=\mathrm{Tr}[\text{T}_{\mathrm{par}}^N]\text{T}_{\mathrm{par}}=\left(\begin{array}{cc}2\mathrm{cosh}[\beta w_+]e^{i\beta \widehat{m}}& 2\mathrm{cosh}[\beta w_0]\\ 2\mathrm{cosh}[\beta w_0]& 2\mathrm{cosh}[\beta w_{}]e^{i\beta \widehat{m}}\end{array}\right)$$

where $`w_0=J_{\mathrm{}}m+\theta `$ and $`w_\pm =w_0\pm 2J_s`$. The identity $`\mathrm{Tr}[\text{T}^N]=\lambda _+^N+\lambda _{}^N`$, in which $`\lambda _\pm `$ are the eigenvalues of the $`2\times 2`$ matrix T enables us to take the limit $`N\mathrm{}`$ in our equations. The integral over $`(m,\widehat{m})`$ is then for $`N\mathrm{}`$ evaluated by steepest descent, and is dominated by the saddle points of the exponent $`\varphi `$. We thus arrive at the transparent result

$$f=\mathrm{extr}\varphi (m,\widehat{m})\{\begin{array}{c}\varphi _{\mathrm{seq}}(m,\widehat{m})=im\widehat{m}m\theta \frac{1}{2}J_{\mathrm{}}m^2\frac{1}{\beta }\mathrm{log}\lambda _+^{\mathrm{seq}}\hfill \\ \varphi _{\mathrm{par}}(m,\widehat{m})=im\widehat{m}m\theta \frac{1}{\beta }\mathrm{log}\lambda _+^{\mathrm{par}}\hfill \end{array}$$
(7)
where $`\lambda _+^{\mathrm{seq}}`$ and $`\lambda _+^{\mathrm{par}}`$ are the largest eigenvalues of $`\text{T}_{\mathrm{seq}}`$ and $`\text{T}_{\mathrm{par}}`$:

$`\lambda _+^{\mathrm{seq}}`$ $`=`$ $`e^{\beta J_s}\mathrm{cosh}[i\beta \widehat{m}]+\left[e^{2\beta J_s}\mathrm{cosh}^2[i\beta \widehat{m}]2\mathrm{sinh}[2\beta J_s]\right]^{\frac{1}{2}}`$
$`\lambda _+^{\mathrm{par}}`$ $`=`$ $`\mathrm{cosh}[\beta w_+]e^{i\beta \widehat{m}}+\mathrm{cosh}[\beta w_{}]e^{i\beta \widehat{m}}+[\mathrm{cosh}^2[\beta w_+]e^{2\beta i\widehat{m}}+`$
$`+\mathrm{cosh}^2[\beta w_{}]e^{2\beta i\widehat{m}}2\mathrm{cosh}[\beta w_+]\mathrm{cosh}[\beta w_{}]+4\mathrm{cosh}^2[\beta w_0]]^{\frac{1}{2}}`$
For simplicity, we will restrict ourselves to the case where $`\theta =0`$; generalisation of what follows to the case of arbitrary $`\theta `$, by using the full form of (7), is not significantly more difficult. The expressions defining the value(s) of the order parameter $`m`$ can now be obtained from the saddle point equations $`_m\varphi (m,\widehat{m})=_{\widehat{m}}\varphi (m,\widehat{m})=0`$. This is a straightforward differentiation task for the sequential case. For parallel, one obtains a set of coupled non-linear equations, namely

$$m=(m,\stackrel{~}{m})\stackrel{~}{m}=(\stackrel{~}{m},m)$$
(8)
where $`\stackrel{~}{m}=i\widehat{m}/J_{\mathrm{}}`$ and $`(,)`$ corresponds to

$$(p,q)=\mathrm{\Delta }(p,q)^1\left(e^{2\beta J_s}\mathrm{sinh}[\beta J_{\mathrm{}}(p+q)]e^{2\beta J_s}\mathrm{sinh}[\beta J_{\mathrm{}}(pq)]\right)$$

$$\mathrm{\Delta }(p,q)=\left[e^{2\beta J_s}\mathrm{sinh}^2[\beta J_{\mathrm{}}(p+q)]+e^{2\beta J_s}\mathrm{sinh}^2[\beta J_{\mathrm{}}(pq)]+2\mathrm{cosh}^2[\beta J_{\mathrm{}}p]+2\mathrm{cosh}^2[\beta J_{\mathrm{}}q]\right]^{\frac{1}{2}}$$

with $`\mathrm{\Delta }(p,q)=\mathrm{\Delta }(q,p)`$. We will now show that the parallel dynamics fixed point problem (8) admits the unique solution $`m=\mathrm{sgn}[J_{\mathrm{}}]\stackrel{~}{m}`$, by the following argument:

$$\begin{array}{ccc}\mathrm{For}J_{\mathrm{}}0:m=\stackrel{~}{m}\hfill & \mathrm{since}& 0(m\stackrel{~}{m})^2=\mathrm{\Omega }(m,\stackrel{~}{m})(m\stackrel{~}{m})\mathrm{sinh}[\beta J_{\mathrm{}}(\stackrel{~}{m}m)]0\\ \mathrm{For}J_{\mathrm{}}<0:m=\stackrel{~}{m}\hfill & \mathrm{since}& 0(m+\stackrel{~}{m})^2=\mathrm{\Omega }(m,\stackrel{~}{m})(m+\stackrel{~}{m})\mathrm{sinh}[\beta J_{\mathrm{}}(\stackrel{~}{m}+m)]0\end{array}$$

$$\mathrm{where}\mathrm{\Omega }(m,\stackrel{~}{m})2e^{2\beta J_s}\mathrm{\Delta }(m,\stackrel{~}{m})^1>0$$

Insertion of these solutions to the original function $`(,)`$ allows us to reduce (8) to a simple 1D fixed-point problem, similar in structure with what follows from the sequential case, namely

$$\begin{array}{ccccccc}\mathrm{sequential}:\hfill & & \widehat{m}=imJ_{\mathrm{}},\hfill & & m=G(m;J_{\mathrm{}},J_s)\hfill & & \\ \mathrm{parallel}:\hfill & & \widehat{m}=imJ_{\mathrm{}},\hfill & & m=G(m;J_{\mathrm{}},J_s)\hfill & & \mathrm{for}J_{\mathrm{}}0\hfill \\ & & \widehat{m}=imJ_{\mathrm{}},\hfill & & m=G(m;J_{\mathrm{}},J_s)\hfill & & \mathrm{for}J_{\mathrm{}}<0\hfill \end{array}$$
(9)
with

$$G(m;J_{\mathrm{}},J_s)=\frac{\mathrm{sinh}[\beta J_{\mathrm{}}m]}{\sqrt{\mathrm{sinh}^2[\beta J_{\mathrm{}}m]+e^{4\beta J_s}}}$$
(10)
The macroscopic observable $`a`$ is generated by differentiating the reduced free energy per neuron (7):

$$\begin{array}{ccccccc}\mathrm{sequential}:\hfill & & \widehat{m}=imJ_{\mathrm{}},\hfill & & a=F(m;J_{\mathrm{}},J_s)\hfill & & \\ \mathrm{parallel}:\hfill & & \widehat{m}=imJ_{\mathrm{}},\hfill & & a=F(m;J_{\mathrm{}},J_s)\hfill & & \mathrm{for}J_{\mathrm{}}0\hfill \\ & & \widehat{m}=imJ_{\mathrm{}},\hfill & & a=F(m;J_{\mathrm{}},J_s)\hfill & & \mathrm{for}J_{\mathrm{}}<0\hfill \end{array}$$
(11)
with

$$F(m;J_{\mathrm{}},J_s)=\frac{\mathrm{cosh}[\beta J_{\mathrm{}}m]\sqrt{\mathrm{sinh}^2[\beta J_{\mathrm{}}m]+e^{4\beta J_s}}+\mathrm{sinh}^2[\beta J_{\mathrm{}}m]e^{4\beta J_s}}{\mathrm{cosh}[\beta J_{\mathrm{}}m]\sqrt{\mathrm{sinh}^2[\beta J_{\mathrm{}}m]+e^{4\beta J_s}}+\mathrm{sinh}^2[\beta J_{\mathrm{}}m]+e^{4\beta J_s}}$$
(12)
Note that in the absence of short-range interactions we recover the familiar Curie-Weiss law $`m=\mathrm{tanh}[\beta J_{\mathrm{}}m]`$ whereas in the absence of long-range interactions (9) and (11) reduce to $`m=0`$ and $`a=\mathrm{tanh}[\beta J_s]`$, as they should. Finally, it is worth noting that equations (9,10,11,12) allow us to derive the physical properties of the parallel dynamics model from those of the sequential dynamics model via simple parameter transformations.

### Phase Transitions & Phase Diagrams

Our main order parameter $`m`$ is to be solved from an equation of the form $`m=G(m)`$, in which $`G(m)=G(m;J_{\mathrm{}},J_s)`$ for both sequential and parallel dynamics with $`J_{\mathrm{}}0`$, whereas $`G(m)=G(m;J_{\mathrm{}},J_s)`$ for parallel dynamics with $`J_{\mathrm{}}<0`$. Note that, due to $`G(0;J_{\mathrm{}},J_s)=0`$, the trivial solution $`m=0`$ always exists. In order to obtain a phase diagram we have to perform a bifurcation analysis of the equations (9,10), and determine the combinations of parameter values for which specific non-zero solutions are created or annihilated (the transition lines). Bifurcations of non-zero solutions occur when simultaneously $`m=G(m)`$ (saddle-point requirement) and $`1=_mG(m)`$ ($`m`$ is in the process of being created/annihilated). Analytical expressions for the lines in the $`(\beta J_s,\beta J_{\mathrm{}})`$ plane where second-order transitions occur between recall states (where $`m0`$) and non-recall states (where $`m=0`$) are obtained by solving the coupled equations $`m=G(m)`$ and $`1=_mG(m)`$ for $`m=0`$. This gives:

$$cont.trans.:\begin{array}{cc}\mathrm{sequential}:\hfill & \beta J_{\mathrm{}}=e^{2\beta J_s}\hfill \\ \mathrm{parallel}:\hfill & \beta J_{\mathrm{}}=e^{2\beta J_s}\mathrm{and}\beta J_{\mathrm{}}=e^{2\beta J_s}\hfill \end{array}$$
(13)
whereas for the macroscopic observable $`a`$ we obtain

$$\mathrm{sequential}/\mathrm{parallel}:a=\mathrm{tanh}[\beta J_s]$$

If along the lines (13) we inspect the behaviour of $`G(m)`$ close to $`m=0`$ we can anticipate the possible existence of first-order transitions, using the properties of $`G(m)`$ for $`m\pm \mathrm{}`$, in combination with $`G(m)=G(m)`$. Precisely at the lines (13) we have $`G(m)=m+\frac{1}{6}G^{^{\prime \prime \prime }}(0).m^3+𝒪(m^5)`$. Since $`lim_m\mathrm{}G(m)=1`$ one knows that, when $`G^{\prime \prime \prime }(0)>0`$, a discontinuous transition must have already taken place earlier, and that away from the lines (13) there will consequently be regions where one finds five solutions of $`m=G(m)`$ (two positive ones, two negative ones and $`m=0`$). Along the lines (13) the condition $`G^{\prime \prime \prime }(0)>0`$ translates into

$$\begin{array}{cc}\mathrm{sequential}:\hfill & \beta J_{\mathrm{}}>\sqrt{3}\mathrm{and}\beta J_s<\frac{1}{4}\mathrm{log}3\hfill \\ \mathrm{parallel}:\hfill & |\beta J_{\mathrm{}}|>\sqrt{3}\mathrm{and}|\beta J_s|<\frac{1}{4}\mathrm{log}3\hfill \end{array}$$
(14)
In the present models it turns out that one can also find an analytical expression for the discontinuous transition lines in the $`(\beta J_s,\beta J_{\mathrm{}})`$ plane, in the form of a parametrisation. For sequential dynamics one finds a single line, parametrised by $`x=\beta J_{\mathrm{}}m[0,\mathrm{})`$:

$$discont.trans.:\beta J_{\mathrm{}}(x)=\sqrt{\frac{x^3}{x\mathrm{tanh}(x)}},\beta J_s(x)=\frac{1}{4}\mathrm{log}\left[\frac{\mathrm{tanh}(x)\mathrm{sinh}^2(x)}{x\mathrm{tanh}(x)}\right]$$
(15)
This can be verified by explicit substitution into (9). Since this parametrisation (15) obeys $`\beta J_s(0)=\frac{1}{4}\mathrm{log}3`$ and $`\beta J_{\mathrm{}}(0)=\sqrt{3}`$, the discontinuous transition indeed starts precisely at the point predicted by the convexity of $`G(m)`$ at $`m=0`$, see (14). In the limit $`x\mathrm{}`$ the slope of (15) approaches $`\beta J_{\mathrm{}}(\mathrm{})/\beta J_s(\mathrm{})2`$. For sequential dynamics the line (15) gives all non-zero solutions of the bifurcation requirements $`m=G(m)`$ and $`1=_mG(m)`$. For parallel dynamics one finds, in addition to (15), a second ‘mirror’ transition line, generated by the transformation $`\{\beta J_{\mathrm{}},\beta J_s\}\{\beta J_{\mathrm{}},\beta J_s\}`$.

Having determined the transition lines in parameter space, we can turn to the phase diagrams. Figure 1 shows the phase diagram for the two types of dynamics, in the $`(\beta J_s,\beta J_{\mathrm{}})`$ plane (note: of the three parameters $`\{\beta ,J_s,J_{\mathrm{}}\}`$ one is redundant). For sequential dynamics we find ($`i`$) a region where $`m=0`$ only ($`ii`$) a region where the trivial as well as two non-trivial states (a positive and a negative one) are all locally stable (selection will thus be based on initial conditions) and ($`iii`$) a region where only $`m0`$ states are localy stable. The transitions ($`i)(iii`$) and ($`ii)(iii`$) are second-order ones (solid lines of figure 1) whereas the transition ($`i)(ii`$) is first-order (dashed line). In region ($`i`$), where $`m=0`$ only, one finds $`a=\mathrm{tanh}[\beta J_s]`$ whereas in regions $`(ii)`$ and $`(iii)`$ $`a`$ is given by the full expression of (11). For parallel dynamics we find, in addition to the sequential phase transitions, a second set of ‘mirror’ transition lines generated by $`\{\beta J_{\mathrm{}},\beta J_s\}\{\beta J_{\mathrm{}},\beta J_s\}`$. In contrast to the sequential case however, here in the region $`\beta J_{\mathrm{}}<0`$ and $`\beta J_s<0`$ where $`m0`$ can be a physical state (lower left corner of the phase diagram of figure 1) one finds 2-cycles between the two $`m0`$ (positive and negative) recall states. This can be inferred from the exact dynamical solution that is available along the line $`J_s=0`$ (see e.g. ), given by the deterministic map $`m(t+1)=\mathrm{tanh}[\beta J_{\mathrm{}}m(t)]`$. This map gives a stable period-2 oscillation for $`\beta J_{\mathrm{}}<1`$, of the form $`m(t)=(1)^tm^{}`$, where $`m^{}=\mathrm{tanh}[\beta |J_{\mathrm{}}|m^{}]`$. In the 2-cycle region one has $`a=lim_N\mathrm{}\frac{1}{N}_i\sigma _{i+1}\mathrm{tanh}[\beta h_i(𝝈)]<0`$. This can be understood on the basis of the (parallel dynamics) identity $`\sigma _{i+1}\mathrm{tanh}[\beta h_i(𝝈)]=\sigma _{i+1}(t)\sigma _i(t+1)`$.

We find that in contrast to models with nearest neighbour interactions only ($`J_{\mathrm{}}=0`$, where no pattern recall will occur), and to models with mean-field interactions only ($`J_s=0`$, where pattern recall can occur), the combination of the two interaction types leads to qualitatively new modes of operation, especially in the competition region, where $`J_{\mathrm{}}>0`$ and $`J_s<0`$ (Hebbian long-range synapses, combined with anti-Hebbian short range ones). The novel features of the diagram can play a useful role: the existence of multiple locally stable states ensures that only sufficiently strong recall cues will evoke pattern recognition; the discontinuity of the transition subsequently ensures that in the latter case the recall will be of a substantial quality. In the case of parallel, similar statements can be made in the opposite region of synaptic competition, but now involving 2-cycles. Since figure 1 cannot show the zero noise region $`(\beta =T^1=\mathrm{})`$, we have also drawn the interesting competition region of the sequential dynamics phase diagram in the $`(J_{\mathrm{}},T)`$ plane, for $`J_s=1`$ (see figure 2, left picture). At $`T=0`$ one finds coexistence of recall states $`(m0)`$ and non-recall states $`(m=0)`$ for any $`J_{\mathrm{}}>0`$, as soon as $`J_s<0`$. In the same figure (right picture) we show the magnitude of the discontinuity in the order parameter $`m`$ along the discontinuous transition line, as a function of $`\beta J_{\mathrm{}}`$.

Finally we show, by way of further illustration of the coexistence mechanism, the value of reduced exponent $`\varphi _{\mathrm{seq}}(m)`$ given in (7), evaluated upon elimination of the auxiliary order parameter $`\widehat{m}`$: $`\varphi (m)\varphi _{\mathrm{seq}}(m,imJ_{\mathrm{}})`$. The result, for the parameter choice $`(\beta ,J_{\mathrm{}})=(2,3)`$ and for three different short-range coupling stengths (corresponding to the three phase regimes: non-zero recall, coexistence and zero recall) is given in figure 3. In the same figure we also give the sequential dynamics bifurcation diagram displaying the value(s) of the overlap $`m`$ as a function of $`\beta J_{\mathrm{}}`$ and for $`\beta J_s=0.6`$ (a line crossing all three phase regimes in figure 1).

## 4 Solution and Phase Diagrams for Arbitrary $`p`$

In the general case (4), where $`p>1`$, the pattern variables $`\{\xi _i^\mu \}`$ can no longer be transformed away as in model (6); here in order to arrive at expressions for the order parameters one has to perform the disorder average of the free energy over the realisation of $`\{\xi _i^\mu \}`$. Due to the addition of nearest neigbour interactions, however, the disorder average has become significantly more complicated than that in infinite-range models, even in the loading regime away from saturation (i.e. $`lim_N\mathrm{}p/N=0`$). We perform the disorder average via a suitable adaptation of 1D random-field techniques, see e.g. . These are based on the derivation of a stochastic process in which observables of systems of size $`N`$ are mapped to observables of systems of size $`N+1`$. We will assume that in the thermodynamic limit this process leads to a unique stationary state. The free energy is then given as an integral over the distribution of a characteristic ratio of conditioned partition functions, and the order parameters follow via differentiation.

### Adaptation of RFIM Techniques

We introduce the notation $`𝒎^{}=(m_0,𝒎)`$, $`\widehat{𝒎}^{}=(\widehat{m}_0,\widehat{𝒎})`$, $`𝝃^{}=(\xi ^0,𝝃)`$ where $`\xi ^0=1`$. We separate the overlap order parameter via insertion of $`1=𝑑𝒎^{}\delta [𝒎^{}1/N_i\sigma _i𝝃_i^{}]`$, and we replace the $`\delta `$-functions by their integral representation giving

$$f=\underset{N\mathrm{}}{lim}\frac{1}{\beta N}\mathrm{log}𝑑𝒎^{}𝑑\widehat{𝒎}^{}e^{\beta N\varphi _N(𝒎^{},\widehat{𝒎}^{})}$$

where

$`\varphi _N^{\mathrm{seq}}(𝒎^{},\widehat{𝒎}^{})`$ $`=`$ $`i𝒎^{}\widehat{𝒎}^{}m_0\theta {\displaystyle \frac{1}{2}}J_{\mathrm{}}^{(1)}m_0^2{\displaystyle \frac{1}{2}}J_{\mathrm{}}^{(2)}𝒎^2{\displaystyle \frac{1}{\beta N}}\mathrm{log}R_N^{\mathrm{seq}}(\widehat{𝒎}^{})`$ (16)
$`\varphi _N^{\mathrm{par}}(𝒎^{},\widehat{𝒎}^{})`$ $`=`$ $`i𝒎^{}\widehat{𝒎}^{}m_0\theta {\displaystyle \frac{1}{\beta N}}\mathrm{log}R_N^{\mathrm{par}}(𝒎^{},\widehat{𝒎}^{})`$ (17)
The quantities $`R_N(𝒎^{},\widehat{𝒎}^{})`$ contain, as in model (5) in section 3, the summation over $`𝝈=(\sigma _1,\mathrm{},\sigma _N)`$ and the short-range neuron interactions

$$R_N^{\mathrm{seq}}(\widehat{𝒎})=\underset{𝝈\{1,1\}^N}{}_{\mathrm{seq}}(𝝈)R_N^{\mathrm{par}}(𝒎^{},\widehat{𝒎}^{})=\underset{𝝈\{1,1\}^N}{}_{\mathrm{par}}(𝝈),$$
(18)
$`_{\mathrm{seq}}(𝝈)`$ $`=`$ $`e^{i\beta _i\sigma _i\widehat{𝒎}^{}𝝃_i^{}}e^{\beta _i\sigma _iJ_{i,i+1}^s\sigma _{i+1}}`$
$`_{\mathrm{par}}(𝝈)`$ $`=`$ $`e^{i\beta _i\sigma _i\widehat{𝒎}^{}𝝃_i^{}}{\displaystyle \underset{i}{}}\mathrm{log}[2\mathrm{cosh}[\beta (J_i^{\mathrm{}}+\sigma _{i1}J_{i1,i}^s+\sigma _{i+1}J_{i,i+1}^s)]]`$
where the short-hands $`J_{i,i+1}^s`$ and $`J_i^{\mathrm{}}`$ correspond to:

$$J_{i,i+1}^s=J_s^{(1)}+J_s^{(2)}𝝃_i𝝃_{i+1}J_i^{\mathrm{}}=\theta +J_{\mathrm{}}^{(1)}m_0+J_{\mathrm{}}^{(2)}𝝃_i𝒎$$

Each of the quantities in (18) is now separated into different constituent parts, defined by conditioning on either the state of the last neuron in the chain (sequential dynamics), or on the states of the last two neurons in the chain (parallel dynamics):

$$\left[\begin{array}{c}R_{N,}^{\mathrm{seq}}\\ R_{N,}^{\mathrm{seq}}\end{array}\right]=\underset{𝝈}{}_{\mathrm{seq}}(𝝈)\left[\begin{array}{c}\delta _{\sigma _N,1}\hfill \\ \delta _{\sigma _N,1}\hfill \end{array}\right]$$
(19)
$$\left[\begin{array}{c}R_{N,}^{\mathrm{par}}\\ R_{N,}^{\mathrm{par}}\\ R_{N,}^{\mathrm{par}}\\ R_{N,}^{\mathrm{par}}\end{array}\right]=\underset{𝝈}{}_{\mathrm{par}}(𝝈)\left[\begin{array}{c}\delta _{\sigma _{N1},1}\delta _{\sigma _N,1}\hfill \\ \delta _{\sigma _{N1},1}\delta _{\sigma _N,1}\hfill \\ \delta _{\sigma _{N1},1}\delta _{\sigma _N,1}\hfill \\ \delta _{\sigma _{N1},1}\delta _{\sigma _N,1}\hfill \end{array}\right]$$
(20)
We now add an extra neuron to the chain. After some simple bookkeeping and assuming non-periodic boundary conditions we derive the following recurrence relations:

$$\mathrm{seq}:\left(\begin{array}{c}R_{N+1,}\hfill \\ R_{N+1,}\hfill \end{array}\right)=\left(\begin{array}{cc}e^{\beta (J_{N,N+1}^si\widehat{𝒎}^{}𝝃_{N+1}^{})}\hfill & e^{\beta (J_{N,N+1}^s+i\widehat{𝒎}^{}𝝃_{N+1}^{})}\hfill \\ e^{\beta (J_{N,N+1}^si\widehat{𝒎}^{}𝝃_{N+1}^{})}\hfill & e^{\beta (J_{N,N+1}^s+i\widehat{𝒎}^{}𝝃_{N+1}^{})}\hfill \end{array}\right)\left(\begin{array}{c}R_{N,}\hfill \\ R_{N,}\hfill \end{array}\right)$$
(21)
$`\mathrm{par}:`$ $`(\begin{array}{c}R_{N+1,}\hfill \\ R_{N+1,}\hfill \end{array})=2Q_{N,N+1}^+(\begin{array}{cc}\frac{e^{i\beta \widehat{𝒎}^{}𝝃_{N+1}^{}}}{Q_{N1,N}^+}A_{(+,+)}^{(N)}& \frac{e^{i\beta \widehat{𝒎}^{}𝝃_{N+1}^{}}}{Q_{N1,N}^{}}A_{(,+)}^{(N)}\\ \frac{e^{i\beta \widehat{𝒎}^{}𝝃_{N+1}^{}}}{Q_{N1,N}^+}A_{(+,)}^{(N)}& \frac{e^{i\beta \widehat{𝒎}^{}𝝃_{N+1}^{}}}{Q_{N1,N}^{}}A_{(,)}^{(N)}\end{array})(\begin{array}{c}R_{N,}\hfill \\ R_{N,}\hfill \end{array})`$ (35)
$`(\begin{array}{c}R_{N+1,}\hfill \\ R_{N+1,}\hfill \end{array})=2Q_{N,N+1}^{}(\begin{array}{cc}\frac{e^{i\beta \widehat{𝒎}^{}𝝃_{N+1}^{}}}{Q_{N1,N}^+}A_{(+,+)}^{(N)}& \frac{e^{i\beta \widehat{𝒎}^{}𝝃_{N+1}^{}}}{Q_{N1,N}^{}}A_{(,+)}^{(N)}\\ \frac{e^{i\beta \widehat{𝒎}^{}𝝃_{N+1}^{}}}{Q_{N1,N}^+}A_{(+,)}^{(N)}& \frac{e^{i\beta \widehat{𝒎}^{}𝝃_{N+1}^{}}}{Q_{N1,N}^{}}A_{(,)}^{(N)}\end{array})(\begin{array}{c}R_{N,}\hfill \\ R_{N,}\hfill \end{array})`$
with the short-hand notation:

$$A_{(\pm ,\pm )}^{(n)}=A_{(\pm ,\pm )}^{(n)}(𝝃_{n1}𝝃_n,𝝃_n𝝃_{n+1},𝝃_n𝒎)=\mathrm{cosh}[\beta (J_n^{\mathrm{}}\pm J_{n1,n}^s\pm J_{n,n+1}^s)]$$

$$Q_{n,n+1}^\pm =Q^\pm (𝝃_n𝝃_{n+1},𝝃_n𝒎)=\mathrm{cosh}[\beta (J_{n+1}^{\mathrm{}}\pm J_{n,n+1}^s)]$$
(36)
The fact that parallel dynamics leads to a recurrence process given by $`4\times 4`$ random matrices, instead of the simpler $`2\times 2`$ random matrices of the sequential case, is due to the appearance of short-range couplings of the form $`cosh[\beta (J_i^{\mathrm{}}+\sigma _{i1}J_{i1,i}^s+\sigma _{i+1}J_{i,i+1}^s)]`$ in (18), rather than the more familiar exponentials. This, in turn, is an immediate consequence of the form of the pseudo-Hamiltonian (3). We also note that in the parallel case the two reduced $`2\times 2`$ matrices differ only by the prefactors $`Q_{N,N+1}^\pm `$. The recurrence matrices (21-35) will form the basis for evaluating the free energy per neuron, which follows from

$$R_N^{\mathrm{seq}}=R_{N,}^{\mathrm{seq}}+R_{N,}^{\mathrm{seq}}$$
(37)
$$R_N^{\mathrm{par}}=R_{N,}^{\mathrm{par}}+R_{N,}^{\mathrm{par}}+R_{N,}^{\mathrm{par}}+R_{N,}^{\mathrm{par}}$$
(38)
### The Stochastic Process for Conditioned Partition Functions

In the spirit of  we next define the ratios between the conditioned quantities (19) and (20), and study their stochastic ‘evolution’ generated by adding new neurons successively to the chain:

$$\mathrm{sequential}:k_{n+1}=e^{2\beta i\widehat{𝒎}^{}𝝃_{n+1}^{}}\frac{R_{n+1,}}{R_{n+1,}}$$
(39)
$$\mathrm{parallel}:k_{n+1}^{(1)}=e^{2\beta i\widehat{𝒎}^{}𝝃_{n+1}^{}}\frac{R_{n+1,}}{R_{n+1,}}k_{n+1}^{(2)}=e^{2\beta i\widehat{𝒎}^{}𝝃_{n+1}^{}}\frac{R_{n+1,}}{R_{n+1,}}$$

$$k_{n+1}^{(3)}=\frac{Q_{n,n+1}^+}{Q_{n,n+1}^{}}\frac{R_{n+1,}}{R_{n+1,}}$$
(40)
Due to their dependence on the random variables $`\{𝝃\}`$, each of these quantities is a stochastic variable. Insertion of (21-35) into the above relations gives

$`\mathrm{sequential}:`$ $`k_{n+1}={\displaystyle \frac{e^{\beta J_{n,n+1}^s}+k_ne^{2\beta i\widehat{𝒎}^{}𝝃_n^{}}e^{\beta J_{n,n+1}^s}}{e^{\beta J_{n,n+1}^s}+k_ne^{2\beta i\widehat{𝒎}^{}𝝃_n^{}}e^{\beta J_{n,n+1}^s}}}`$ (41)
$`\mathrm{parallel}:`$ $`k_{n+1}^{(1)}={\displaystyle \frac{A_{(+,)}^{(n)}k_n^{(2)}+A_{(,)}^{(n)}k_n^{(1)}k_n^{(3)}}{A_{(+,+)}^{(n)}k_n^{(2)}+A_{(,+)}^{(n)}k_n^{(1)}k_n^{(3)}}}k_{n+1}^{(2)}={\displaystyle \frac{A_{(+,)}^{(n)}+A_{(,)}^{(n)}k_n^{(3)}}{A_{(+,+)}^{(n)}+A_{(,+)}^{(n)}k_n^{(3)}}}`$ (42)
$`k_{n+1}^{(3)}={\displaystyle \frac{A_{(+,)}^{(n)}k_n^{(1)}k_n^{(2)}+A_{(,)}^{(n)}k_n^{(1)}k_n^{(2)}k_n^{(3)}}{A_{(+,)}^{(n)}k_n^{(2)}+A_{(,)}^{(n)}k_n^{(1)}k_n^{(3)}}}e^{2\beta i\widehat{𝒎}^{}𝝃_n^{}}`$
For parallel dynamics we observe that, if the identity $`k_n^{(1)}=k_n^{(2)}`$ is true, then also $`k_{n+1}^{(1)}=k_{n+1}^{(2)}`$. Furthermore, it can be easily checked that this is the case for $`n=2`$. We are thus guaranteed that $`k_n^{(1)}=k_n^{(2)}`$ for all $`n2`$, and now all three quantities $`\{k_i^{(1)},k_i^{(2)},k_i^{(3)}\}`$ can (for any $`i=2,\mathrm{},N`$) be expressed in terms of a single stochastic variable which we will take to be $`k_i^{(1)}k_i`$. This simplifies considerably the description of the parallel dynamics stochastic process:

$$k_{n+1}^{(1)}=k_{n+1}k_{n+1}^{(2)}=\frac{A_{(+,)}^{(n)}+A_{(,)}^{(n)}k_{n1}e^{2\beta i\widehat{𝒎}^{}𝝃_{n1}^{}}}{A_{(+,+)}^{(n)}+A_{(,+)}^{(n)}k_{n1}e^{2\beta i\widehat{𝒎}^{}𝝃_{n1}^{}}}k_{n+1}^{(3)}=k_ne^{2\beta i\widehat{𝒎}^{}𝝃_n^{}}$$

We have thus obtained 1D stochastic maps for both sequential and parallel dynamics:

$`\mathrm{seq}:k_{i+1}`$ $`=`$ $`\psi _{\mathrm{seq}}\left(k_j;𝝃_j𝝃_{j+1},\widehat{𝒎}^{}𝝃_j^{}|ji\right)`$
$`=`$ $`{\displaystyle \frac{e^{\beta J_{n,n+1}^s}+k_ne^{2\beta i\widehat{𝒎}^{}𝝃_n^{}}e^{\beta J_{n,n+1}^s}}{e^{\beta J_{n,n+1}^s}+k_ne^{2\beta i\widehat{𝒎}^{}𝝃_n^{}}e^{\beta J_{n,n+1}^s}}}`$
$`\mathrm{par}:k_{i+1}`$ $`=`$ $`\psi _{\mathrm{par}}\left(k_{j1};𝝃_{j1}𝝃_j,𝝃_j𝝃_{j+1},𝒎𝝃_j,\widehat{𝒎}^{}𝝃_{j1}^{}|ji\right)`$
$`=`$ $`{\displaystyle \frac{\mathrm{cosh}[\beta (J_i^{\mathrm{}}+J_{i1,i}^sJ_{i,i+1}^s)]+e^{2\beta i\widehat{𝒎}^{}𝝃_{i1}^{}}k_{i1}\mathrm{cosh}[\beta (J_i^{\mathrm{}}J_{i1,i}^sJ_{i,i+1}^s)]}{\mathrm{cosh}[\beta (J_i^{\mathrm{}}+J_{i1,i}^s+J_{i,i+1}^s)]+e^{2\beta i\widehat{𝒎}^{}𝝃_{i1}^{}}k_{i1}\mathrm{cosh}[\beta (J_i^{\mathrm{}}J_{i1,i}^s+J_{i,i+1}^s)]}}`$
We observe that, in contrast to the sequential dynamics mapping, parallel dynamics distinguishes between even and odd sites, and produces two independent sets of stochastic variables $`\{k_j\}`$.

### Disorder Averaging and the Free Energy

One can now express the non-trivial part of the free energy in terms of the above ratio’s $`k_i`$, using (37-38):

$$\frac{1}{\beta N}\mathrm{log}R_N^{\mathrm{seq}}=\frac{1}{\beta N}\mathrm{log}R_{N,}^{\mathrm{seq}}+𝒪(\frac{1}{N})\frac{1}{\beta N}\mathrm{log}R_N^{\mathrm{par}}=\frac{1}{\beta N}\mathrm{log}R_{N,}^{\mathrm{par}}+𝒪(\frac{1}{N})$$

Upon also using (39-40) and (41-42) iteratively, in order to map the quantities $`\{R_{n,}^{\mathrm{seq}},R_{n,}^{\mathrm{par}}\}`$ onto $`\{R_{n1,}^{\mathrm{seq}},R_{n1,}^{\mathrm{par}}\}`$ (for all $`nN`$), one can write the above expressions in the form:

$`{\displaystyle \frac{1}{\beta N}}\mathrm{log}R_{N,}^{\mathrm{seq}}(𝒎^{})`$ $`=`$ $`{\displaystyle \frac{1}{\beta N}}{\displaystyle \underset{i}{}}\mathrm{log}\left[e^{\beta J_{i,i+1}^s}+e^{\beta J_{i,i+1}^s}e^{2\beta i\widehat{𝒎}^{}𝝃_i^{}}k_i\right]`$
$`+{\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}i\widehat{𝒎}^{}𝝃_{i+1}^{}+𝒪({\displaystyle \frac{1}{N}})`$
$`{\displaystyle \frac{1}{\beta N}}\mathrm{log}R_{N,}^{\mathrm{par}}(𝒎^{},\widehat{𝒎}^{})`$ $`=`$ $`{\displaystyle \frac{1}{\beta N}}{\displaystyle \underset{i}{}}\mathrm{log}[2\mathrm{cosh}[\beta (J_{i+1}^{\mathrm{}}+J_{i,i+1}^s+J_{i+1,i+2}^s)]+`$
$`+2\mathrm{cosh}[\beta (J_{i+1}^{\mathrm{}}J_{i,i+1}^s+J_{i+1,i+2}^s)]e^{2\beta i\widehat{𝒎}^{}𝝃_i^{}}k_i]`$
$`+{\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}i\widehat{𝒎}^{}𝝃_{i+1}^{}{\displaystyle \frac{1}{\beta N}}{\displaystyle \underset{i}{}}\mathrm{log}\left[Q_{i+1,i+2}^+/Q_{i+1,i+2}^{}\right]+𝒪({\displaystyle \frac{1}{N}})`$
where the terms $`Q_{j,k}^\pm `$ are defined at (36). For $`N\mathrm{}`$ the above expressions are self-averaging and the boundary terms vanish, giving

$$\mathrm{sequential}:\underset{N\mathrm{}}{lim}\frac{1}{\beta N}\mathrm{log}R_N^{\mathrm{seq}}(\widehat{𝒎}^{})=$$

$$=i\widehat{𝒎}^{}𝝃^{}_𝝃\frac{1}{\beta }\underset{𝝃^{}𝝃}{}𝑑k^{}P^\mathrm{s}(k^{},𝝃^{},𝝃)\mathrm{log}\left[e^{\beta (J_s^{(1)}+J_s^{(2)}𝝃^{}𝝃)}+e^{\beta (J_s^{(1)}+J_s^{(2)}𝝃^{}𝝃)}e^{2\beta i\widehat{𝒎}^{}𝝃_{}^{}{}_{}{}^{}}k^{}\right]$$
(45)
$$\mathrm{parallel}:\underset{N\mathrm{}}{lim}\frac{1}{\beta N}\mathrm{log}R_N^{\mathrm{par}}(𝒎^{},\widehat{𝒎}^{})=$$

$$=i\widehat{𝒎}^{}𝝃_{}^{}{}_{}{}^{}_𝝃^{}\frac{1}{\beta }\mathrm{log}2\frac{1}{\beta }\mathrm{log}\left[\frac{Q^+(𝝃^{}𝝃,𝝃𝒎)}{Q^{}(𝝃^{}𝝃,𝝃𝒎)}\right]_{\{𝝃^{},𝝃\}}$$

$$\frac{1}{\beta }\underset{𝝃^{\prime \prime }𝝃^{}𝝃}{}𝑑k^{\prime \prime }P^\mathrm{p}(k^{\prime \prime },𝝃^{\prime \prime },𝝃^{},𝝃)\mathrm{log}[A_{(+,+)}(𝝃^{\prime \prime }𝝃^{},𝝃^{}𝝃,𝝃^{}𝒎)+A_{(,+)}(𝝃^{\prime \prime }𝝃^{},𝝃^{}𝝃,𝝃^{}𝒎)k^{\prime \prime }e^{2\beta i\widehat{𝒎}^{}𝝃_{}^{}{}_{}{}^{\prime \prime }}]$$
(46)
in which the probability distributions $`P^\mathrm{s}(\mathrm{})`$ and $`P^\mathrm{p}(\mathrm{})`$ are defined as follows:

$`\mathrm{sequential}:P^\mathrm{s}(k^{},𝝃^{},𝝃)`$ $`=`$ $`\underset{N\mathrm{}}{lim}{\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}\delta [k^{}k_i]\delta [𝝃^{}𝝃_i]\delta [𝝃𝝃_{i+1}]_{\mathrm{eq}}`$ (47)
$`=`$ $`{\displaystyle \frac{1}{2^p}}P^\mathrm{s}(k^{},𝝃^{})`$
$`\mathrm{parallel}:P^\mathrm{p}(𝒌^{\prime \prime },𝝃^{\prime \prime },𝝃^{},𝝃)`$ $`=`$ $`\underset{N\mathrm{}}{lim}{\displaystyle \frac{1}{N}}{\displaystyle \underset{i}{}}\delta [𝒌^{\prime \prime }𝒌_{i1}]\delta [𝝃^{\prime \prime }𝝃_{i1}]\delta [𝝃^{}𝝃_i]\delta [𝝃𝝃_{i+1}]_{\mathrm{eq}}`$ (48)
$`=`$ $`{\displaystyle \frac{1}{2^{2p}}}P^\mathrm{p}(𝒌^{\prime \prime },𝝃^{\prime \prime })`$
In order to arrive at factorisation of these joint probability distributions we have used the fact that, due to the form of the stochastic maps (4) and (4), the quantities $`\{k_n\}`$ are (for both types of dynamics) independent of $`\{𝝃_{n+1}\}`$, i.e. of the pattern variables at the next site. The dependence on $`\{𝝃_n\}`$, however, does not allow us to write down immediately (4-4) in terms of a Markov process. We thus introduce the auxiliary variables $`\{𝝀_n\}`$, so that $`\{𝝃_n\}`$ itself becomes one of the stochastically evolving variables, and the extended process becomes indeed Markovian:

$$\mathrm{seq}:(\begin{array}{c}k\hfill \\ 𝝀\hfill \end{array})_{i+1}=\mathrm{\Psi }_{\mathrm{seq}}[(\begin{array}{c}k\hfill \\ 𝝀\hfill \end{array})_i;𝝃_{i+1}]=(\begin{array}{c}\psi _{\mathrm{seq}}(k_i;𝝀_i𝝃_{i+1},\widehat{𝒎}^{}𝝀_i^{})\\ 𝝃_{i+1}\end{array})$$

$$\mathrm{par}:(\begin{array}{c}k\hfill \\ 𝝀\hfill \end{array})_{i+1}=\mathrm{\Psi }_{\mathrm{par}}[(\begin{array}{c}k\hfill \\ 𝝀\hfill \end{array})_{i1};𝝃_{i+1},𝝃_i]=(\begin{array}{c}\psi _{\mathrm{par}}(k_{i1};𝝀_{i1}𝝃_i,𝝃_i𝝃_{i+1},𝒎𝝃_i,\widehat{𝒎}^{}𝝀_{i1}^{})\\ 𝝃_{i+1}\end{array})$$

Equivalently:

$$\mathrm{sequential}:P_{i+1}^\mathrm{s}(k,𝝀)=\frac{1}{2^p}\underset{𝝀^{}}{}dk^{}W_\mathrm{s}[\left(\begin{array}{c}k^{}\hfill \\ 𝝀^{}\hfill \end{array}\right)\left(\begin{array}{c}k\hfill \\ 𝝀\hfill \end{array}\right)]P_i^\mathrm{s}(k^{},𝝀^{})$$
(49)
$$\mathrm{parallel}:P_{i+1}^\mathrm{p}(k,𝝀)=\frac{1}{2^p}\underset{𝝀^{\prime \prime }}{}dk^{\prime \prime }W_\mathrm{p}[\left(\begin{array}{c}k^{\prime \prime }\hfill \\ 𝝀^{\prime \prime }\hfill \end{array}\right)\left(\begin{array}{c}k\hfill \\ 𝝀\hfill \end{array}\right)]P_i^\mathrm{p}(k^{\prime \prime },𝝀^{\prime \prime })$$
(50)
where for the transition probabilities $`W_\mathrm{s}[\mathrm{}]`$ and $`W_\mathrm{p}[\mathrm{}]`$ we define:

$$\mathrm{sequential}:W_\mathrm{s}[\left(\begin{array}{c}k^{}\hfill \\ 𝝀^{}\hfill \end{array}\right)\left(\begin{array}{c}k\hfill \\ 𝝀\hfill \end{array}\right)]=\delta _{𝝃,𝝀}\delta [k\psi _{\mathrm{seq}}(k^{},𝝀^{}𝝃,\widehat{𝒎}^{}𝝀_{}^{}{}_{}{}^{})]_𝝃$$

$$\mathrm{parallel}:W_\mathrm{p}[\left(\begin{array}{c}k^{\prime \prime }\hfill \\ 𝝀^{\prime \prime }\hfill \end{array}\right)\left(\begin{array}{c}k\hfill \\ 𝝀\hfill \end{array}\right)]=\delta _{𝝃,𝝀}\delta [k\psi _{\mathrm{par}}(k^{\prime \prime },𝝀^{\prime \prime }𝝃^{},𝝃^{}𝝃,𝒎𝝃^{},\widehat{𝒎}^{}𝝀_{}^{}{}_{}{}^{\prime \prime })]_{𝝃𝝃^{}}$$

The evaluation of the non-trivial part of the free energies (45,46) has now been reduced to determining the (stationary) distributions of (49-50), which will be denoted as $`P_{\mathrm{}}^\mathrm{s}(k)`$ and $`P_{\mathrm{}}^\mathrm{p}(k)`$. To achieve this final objective we follow  and introduce the integrated densities $`\widehat{P}(k)=_0^k𝑑zP(z)`$. After some algebra to eliminate the $`\delta `$-functions of the transition probabilities via the identity $`𝑑x\delta [g(x)]f(x)=f\left(g^{inv}(0)\right)/|_xg(0)|`$, we then obtain the following recursive relation:

$`\mathrm{sequential}:\widehat{P}_{i+1}^\mathrm{s}(k,𝝀)`$ $`=`$ $`{\displaystyle \frac{1}{2^p}}{\displaystyle \underset{𝝀^{}}{}}\widehat{P}_i^\mathrm{s}(B_\mathrm{s}(k;𝝀^{}𝝀,\widehat{𝒎}^{}𝝀_{}^{}{}_{}{}^{}),𝝀^{})`$ (51)
$`\mathrm{parallel}:\widehat{P}_{i+1}^\mathrm{p}(k,𝝀)`$ $`=`$ $`{\displaystyle \frac{1}{2^{2p}}}{\displaystyle \underset{𝝀^{\prime \prime }}{}}{\displaystyle \underset{𝝃^{}}{}}\widehat{P}_i^\mathrm{p}(B_\mathrm{p}(k;𝝀^{\prime \prime }𝝃^{},𝝃^{}𝝀,𝒎𝝃^{},\widehat{𝒎}^{}𝝀_{}^{}{}_{}{}^{}),𝝀^{\prime \prime })`$ (52)
where the functions $`B_\mathrm{s}(k)`$ and $`B_\mathrm{p}(k)`$ correspond to

$$B_\mathrm{s}(k)=\frac{ke^{\beta J_s}e^{\beta J_s}}{e^{\beta J_s}ke^{\beta J_s}}e^{2\beta i\widehat{𝒎}^{}𝝀_{}^{}{}_{}{}^{}}B_p(k)=\frac{kA_{(+,+)}A_{(+,)}}{A_{(,)}kA_{(,+)}}e^{2\beta i\widehat{𝒎}^{}𝝀_{}^{}{}_{}{}^{\prime \prime }}$$

For notation simplicity we have removed the arguments of $`J_s`$ and $`A_{(\pm ,\pm )}`$. These correspond to

$$\begin{array}{ccc}J_s\hfill & =& J_s(𝝀^{}𝝀)\hfill \\ & =& J_s^{(1)}+J_s^{(2)}𝝀^{}𝝀\hfill \end{array}\begin{array}{ccc}A_{(\pm ,\pm )}\hfill & =& A_{(\pm ,\pm )}(𝝀^{\prime \prime }𝝃^{},𝝃^{}𝝀,𝒎𝝃^{})\hfill \\ & =& \mathrm{cosh}[\beta (J^{\mathrm{}}(𝒎𝝃^{},m_0,\theta )\pm J_s(𝝀^{\prime \prime }𝝃^{})\pm J_s(𝝃^{}𝝀))]\hfill \end{array}$$

We have thus derived as yet fully exact expressions for the disorder-averaged free energies (45,46), as integrals over the distribution of the stochastic quantities $`\{k_{\mathrm{seq}},k_{\mathrm{par}}\}`$, which can be evaluated by numerical iteration of the recursive relation (51-52).

### The ‘Pure State’ Ansatz and Symmetries of the Model

We will now look for the so-called ‘pure-state’ solutions: we will assume that for $`N\mathrm{}`$ there is only one component $`m_\mu `$ of order $`𝒪(1)`$ (for simplicity we will take $`\mu =1`$) whereas for all $`\mu >1`$ $`m_\mu `$ is of order $`𝒪(N^{\frac{1}{2}})`$: $`𝒎^{}=(m_0,m,0,\mathrm{},0)`$ and $`\widehat{𝒎}^{}=(\widehat{m}_0,\widehat{m},0,\mathrm{},0)`$. This is the standard ansatz made in (infinite-range) associative memory models, which gives the dominant states of the system. The fact that this is also true for the present type of models is supported by numerical simulations, see section 6. Using this simplification we will now prove that the stationary integrated densities $`\widehat{P}_{\mathrm{}}^\mathrm{s}(\mathrm{})`$ and $`\widehat{P}_{\mathrm{}}^\mathrm{p}(\mathrm{})`$ are independent of all non-condensed pattern components. This will be shown by induction. We will show that if $`\widehat{P}_i^\mathrm{s}(k,𝝀)`$ and $`\widehat{P}_i^\mathrm{p}(k,𝝀)`$ are independent of $`\lambda _\mu `$ for some $`\mu >1`$, then $`\widehat{P}_{i+1}^\mathrm{s}(k,𝝀)`$ and $`\widehat{P}_{i+1}^\mathrm{p}(k,𝝀)`$ will also be independent of $`\lambda _\mu `$. This will then immediately imply that if the densities are independent of $`\{\lambda _\mu \}`$ for all $`\mu >1`$ at step $`i`$ of the process (51-52), this will remain so for all $`\mu >1`$ at any step $`j>i`$. Usage of the spin-flip operator $`F_\mu 𝒙=(x_1,\mathrm{},x_\mu ,\mathrm{},x_p)`$ allows us to write

$$\begin{array}{ccc}\mathrm{seq}:\mathrm{\Delta }_{i+1}^\mathrm{s}\hfill & =& \widehat{P}_{i+1}^\mathrm{s}(k,F_\mu 𝝀)\widehat{P}_{i+1}^\mathrm{s}(k,𝝀)\hfill \\ & =& \frac{1}{2^p}_𝝀^{}\left[\widehat{P}_i^\mathrm{s}(B_\mathrm{s}(k;𝝀^{}𝝀,\widehat{m}_0+\widehat{m}\lambda ^{}),F_\mu 𝝀^{})\widehat{P}_i^\mathrm{s}(B_\mathrm{s}(k;𝝀^{}𝝀,\widehat{m}_0+\widehat{m}\lambda ^{}),𝝀^{})\right]\hfill \\ \mathrm{par}:\mathrm{\Delta }_{i+1}^\mathrm{p}\hfill & =& \widehat{P}_{i+1}^\mathrm{p}(k,F_\mu 𝝀)\widehat{P}_{i+1}^\mathrm{p}(k,𝝀)\hfill \\ & =& \frac{1}{2^{2p}}_{𝝀^{\prime \prime },𝝃^{}}[\widehat{P}_i^\mathrm{p}(B_\mathrm{p}(k;𝝀^{\prime \prime }𝝃^{},𝝃^{}𝝀,m_0+m\xi ^{},\widehat{m}_0+\widehat{m}\lambda ^{\prime \prime }),F_\mu 𝝀^{})\hfill \\ & & \widehat{P}_i^\mathrm{p}(B_\mathrm{p}(k;𝝀^{\prime \prime }𝝃^{},𝝃^{}𝝀,m_0+m\xi ^{},\widehat{m}_0+\widehat{m}\lambda ^{\prime \prime }),𝝀^{})]\hfill \end{array}$$

If at step $`i`$ the identities $`\mathrm{\Delta }_i^\mathrm{s}=\mathrm{\Delta }_i^\mathrm{p}=0`$ are true, then $`\widehat{P}_i(B(k),F_\mu 𝝀)=\widehat{P}_i(B(k),𝝀)`$ (for both types of dynamics) and thus $`\mathrm{\Delta }_{i+1}^\mathrm{s}=\mathrm{\Delta }_{i+1}^\mathrm{p}=0`$. Upon choosing suitable initial conditions we can thus construct equilibrium integrated densities with the stated property; combination with the assumed ergodicity of the process then implies that the unique solution must have the property. This completes the proof. We can consequently write the stationary integrated densities in the form: $`\widehat{P}_{\mathrm{}}^\mathrm{s}(k,\lambda )`$ and $`\widehat{P}_{\mathrm{}}^\mathrm{p}(k,\lambda )`$, where $`\lambda \{1,1\}`$ corresponds to the condensed pattern component:

$`\mathrm{sequential}:\widehat{P}_{i+1}^\mathrm{s}(k,\lambda )`$ $`=`$ $`{\displaystyle \frac{1}{2}}{\displaystyle \underset{\lambda ^{}=\pm 1}{}}\widehat{P}_i^\mathrm{s}(B_\mathrm{s}(k;\lambda ^{}\lambda ,\widehat{m}_0+\widehat{m}\lambda ^{}),\lambda ^{})`$ (53)
$`\mathrm{parallel}:\widehat{P}_{i+1}^\mathrm{p}(k,\lambda )`$ $`=`$ $`{\displaystyle \frac{1}{4}}{\displaystyle \underset{\lambda ^{\prime \prime },\xi ^{}=\pm 1}{}}\widehat{P}_i^\mathrm{p}(B_\mathrm{p}(k;\lambda ^{\prime \prime }\xi ^{},\xi ^{}\lambda ,m_0+m\xi ^{},\widehat{m}_0+\widehat{m}\lambda ^{\prime \prime }),\lambda ^{\prime \prime })`$ (54)
One can also exploit the fact that the pattern variables $`\{𝝀\}_{\mathrm{seq}}`$ and $`\{𝝀\}_{\mathrm{par}}`$ appear in inner-products only, to simplify the non-trivial integrated expressions of equation (45,46). Since the argumentation will be qualitative we will forget about the details of these expressions and we will denote the integrated logarithmic expressions simply by

$$\mathrm{\Phi }_\mathrm{s}(𝝀^{}𝝀)=\mathrm{log}\left[e^{\beta J_s(𝝀^{}𝝀)}+e^{\beta J_s(𝝀^{}𝝀)}+k^{}e^{2\beta i(\widehat{m}_0+\widehat{m}\lambda ^{})}\right]$$

$$\mathrm{\Phi }_\mathrm{p}(𝝀^{\prime \prime }𝝃^{},𝝃^{}𝝀)=\mathrm{log}\left[A_{(+,+)}(𝝀^{\prime \prime }𝝃^{},𝝃^{}𝝀,\xi ^{}m)+A_{(,+)}(𝝀^{\prime \prime }𝝃^{},𝝃^{}𝝀,\xi ^{}m)k^{\prime \prime }e^{2\beta i(\widehat{m}_0+\widehat{m}\lambda ^{\prime \prime })}\right]$$

Upon using the gauge transformations: $`\lambda _\mu =\tau _\mu \lambda _\mu ^{}`$ for sequential dynamics and $`\lambda _\mu ^{\prime \prime }=\tau _\mu \xi _\mu ^{}`$, $`\lambda _\mu =\eta _\mu \xi _\mu ^{}`$ for parallel dynamics, where $`\tau _\mu `$ and $`v_\mu `$ are auxiliary Ising variables, the non-trivial parts of expressions (45) and (46) –corresponding to the integrals over the stochastic variables $`\{k_j\}`$ after the distribution factorisation (47) and (48) and the ‘pure state’ ansatz– take the form:

$`\mathrm{sequential}:`$ $`{\displaystyle \frac{1}{2^p}}{\displaystyle \underset{𝝀^{}}{}}{\displaystyle \underset{𝝀}{}}{\displaystyle 𝑑k^{}P_{\mathrm{}}^\mathrm{s}(k^{},\lambda ^{})\mathrm{\Phi }_\mathrm{s}(𝝀^{}𝝀;k^{},\widehat{m}_0+\widehat{m}\lambda ^{})}`$
$`={\displaystyle \frac{p1}{2^p}}{\displaystyle \underset{\lambda ^{}}{}}{\displaystyle \underset{𝝉}{}}{\displaystyle 𝑑k^{}P_{\mathrm{}}^\mathrm{s}(k^{},\lambda ^{})\mathrm{\Phi }_\mathrm{s}(\underset{\mu =1}{\overset{p}{}}\tau _\mu ;k^{},\widehat{m}_0+\widehat{m}\lambda ^{})}`$
$`={\displaystyle \frac{p1}{2^p}}{\displaystyle \underset{\lambda ^{}=\pm 1}{}}{\displaystyle \underset{𝒩_\tau =1}{\overset{p}{}}}{\displaystyle 𝑑k^{}(\begin{array}{c}p\hfill \\ 𝒩_\tau \hfill \end{array})P_{\mathrm{}}^\mathrm{s}(k^{},\lambda ^{})\mathrm{\Phi }_\mathrm{s}(2𝒩_\tau p;k^{},\widehat{m}_0+\widehat{m}\lambda ^{})}`$
$`\mathrm{parallel}:`$ $`{\displaystyle \frac{1}{2^{2p}}}{\displaystyle \underset{𝝀^{\prime \prime }}{}}{\displaystyle \underset{𝝃^{}}{}}{\displaystyle \underset{𝝀}{}}{\displaystyle 𝑑k^{\prime \prime }P_{\mathrm{}}^\mathrm{p}(k^{\prime \prime },\lambda ^{\prime \prime })\mathrm{\Phi }_\mathrm{p}(𝝀^{\prime \prime }𝝃^{},𝝃^{}𝝀,m\xi ^{};k^{\prime \prime },\widehat{m}_0+\widehat{m}\lambda ^{\prime \prime })}`$
$`={\displaystyle \frac{p1}{2^{2p}}}{\displaystyle \underset{𝝉}{}}{\displaystyle \underset{\xi ^{}}{}}{\displaystyle \underset{𝜼}{}}{\displaystyle 𝑑k^{\prime \prime }P_{\mathrm{}}^\mathrm{p}(k^{\prime \prime },\tau _1\xi ^{})\mathrm{\Phi }_\mathrm{p}(\tau _1+\underset{\mu =2}{\overset{p}{}}\tau _\mu ,\underset{\mu =1}{\overset{p}{}}\eta _\mu ,m\xi ^{};k^{\prime \prime },\widehat{m}_0+\widehat{m}\tau _1\xi ^{})}`$
$`={\displaystyle \frac{(p1)^2}{2^{2p}}}{\displaystyle \underset{\tau _1=\pm 1}{}}{\displaystyle \underset{\xi ^{}=\pm 1}{}}{\displaystyle \underset{𝒩_\tau =2}{\overset{p}{}}}{\displaystyle \underset{𝒩_\eta =1}{\overset{p}{}}}{\displaystyle }dk^{\prime \prime }(\begin{array}{c}p1\hfill \\ 𝒩_\tau \hfill \end{array})(\begin{array}{c}p\hfill \\ 𝒩_\eta \hfill \end{array})P_{\mathrm{}}^\mathrm{p}(k^{\prime \prime },\tau _1\xi ^{})\times `$
$`\times \mathrm{\Phi }_\mathrm{p}(\tau _1+2𝒩_\tau p+1,2𝒩_\eta p,m\xi ^{};k^{\prime \prime },\widehat{m}_0+\widehat{m}\tau _1\xi ^{})`$
where $`𝒩_\tau `$ and $`𝒩_\eta `$ represent the number of neurons with states equal to $`+1`$ in the configurations of $`𝝉`$, $`𝜼`$ respectively. We have thus replaced all summations over the $`2^p`$ configurations of the vectors $`\{𝝀,𝝀^{}\}_{\mathrm{seq}}`$ and $`\{𝝀^{\prime \prime },𝝃^{},𝝀\}_{\mathrm{par}}`$ by summations over binary- and $`p`$-state variables.

### Phase Diagrams

In order to calculate phase transitions and draw phase diagrams we will first calculate the free energy surfaces (16-17), which at this stage are still functions of the order parameters $`m`$ and $`\widehat{m}`$. For simplicity we will now set $`\theta =0`$. The distributions $`P_{\mathrm{}}^\mathrm{s}(\mathrm{})`$ and $`P_{\mathrm{}}^\mathrm{p}(\mathrm{})`$ can be calculated numerically via iteration of (53-54), and bifurcations of the non-trivial values for the pure state overlap from the trivial solution (if they exist) will then be given as the solutions of the following fixed point problems:

$$\begin{array}{cccccc}\hfill \mathrm{sequential}:& m=_{i\widehat{m}}F_{\mathrm{seq}}(i\widehat{m})\hfill & & \mathrm{at}& & i\widehat{m}=J_{\mathrm{}}^{(2)}m\hfill \\ \hfill \mathrm{parallel}:& m=_{i\widehat{m}}F_{\mathrm{par}}(i\widehat{m},m)\hfill & & \mathrm{at}& & i\widehat{m}=_mF_{\mathrm{par}}(i\widehat{m},m)\hfill \end{array}$$
(57)
with

$`F_{\mathrm{seq}}(i\widehat{m})`$ $`=`$ $`{\displaystyle \frac{p1}{2^p\beta }}{\displaystyle \underset{\lambda ^{}\pm 1}{}}{\displaystyle \underset{𝒩_\tau =1}{\overset{p}{}}}{\displaystyle }dk^{}(\begin{array}{c}p\\ 𝒩_\tau \end{array})P_{\mathrm{}}^\mathrm{s}(k^{},\lambda ^{})\times `$ (61)
$`\times \mathrm{log}\left[e^{\beta [J_s^{(1)}+J_s^{(2)}(2𝒩_\tau p)]}+e^{\beta [J_s^{(1)}+J_s^{(2)}(2𝒩_\tau p)]}k^{}e^{2\beta i\widehat{m}\lambda ^{}}\right]`$
$`F_{\mathrm{par}}(i\widehat{m},m)`$ $`=`$ $`{\displaystyle \frac{(p1)^2}{2^{2p}\beta }}{\displaystyle \underset{\tau _1=\pm 1}{}}{\displaystyle \underset{\xi ^{}=\pm 1}{}}{\displaystyle \underset{𝒩_\tau =2}{\overset{p}{}}}{\displaystyle \underset{𝒩_\eta =1}{\overset{p}{}}}{\displaystyle }dk^{\prime \prime }(\begin{array}{c}p1\\ 𝒩_\tau \end{array})(\begin{array}{c}p\\ 𝒩_\eta \end{array})P_{\mathrm{}}^\mathrm{p}(k^{\prime \prime },\tau _1\xi ^{})\times `$ (69)
$`\times \mathrm{log}[2\mathrm{cosh}[\beta J_{\mathrm{}}^{(2)}m\xi ^{}+2\beta J_s^{(1)}+\beta J_s^{(2)}(\tau _12p+2𝒩_\tau +2𝒩_\eta +1)]+`$
$`+2\mathrm{cosh}[\beta J_{\mathrm{}}^{(2)}m\xi ^{}+\beta J_s^{(2)}[2𝒩_\eta 2𝒩_\tau +\tau _1+1]]k^{\prime \prime }e^{2\beta i\widehat{m}\tau _1\xi ^{}}]`$
$`{\displaystyle \frac{p1}{2^{2p}\beta }}{\displaystyle \underset{\xi =\pm 1}{}}{\displaystyle \underset{𝒩=1}{\overset{p}{}}}\left(\begin{array}{c}p\\ 𝒩\end{array}\right)\mathrm{log}\left[{\displaystyle \frac{\mathrm{cosh}[\beta J_{\mathrm{}}^{(2)}m\xi +J_s^{(1)}+J_s^{(2)}(2𝒩p)]}{\mathrm{cosh}[\beta J_{\mathrm{}}^{(2)}m\xi J_s^{(1)}J_s^{(2)}(2𝒩p)]}}\right]`$
As in the solution of model (5) in section 3, the parallel dynamics fixed-point problem (57,69) takes the form of a set of coupled equations, which makes the evaluation of bifurcation points essentially more laborious.

The solutions of the above equations, for $`J_s^{(1)}=J_{\mathrm{}}^{(1)}=0`$ and for $`p=2`$, are shown in the phase diagrams of figure 4. One distinguishes between four different regions, dependent on the number of locally stable ‘pure state’ solutions: region A with $`m=0`$ only, region B with one locally stable $`m>0`$ state (and one $`m<0`$), regions C with two locally stable $`m>0`$ states (and two $`m<0`$ ones) and region F with three $`m>0`$ states (and three $`m<0`$ ones). Note that region F is created at the point where regions C start overlapping. The transition A$``$B is second-order, whereas B$``$C and F$``$C are first-order. The two qualitatively different types of bifurcations are also shown in figure 5 (left picture), where we draw the solution(s) of the overlap $`m`$ as a function of $`\beta J_{\mathrm{}}`$ along the line $`\beta J_s=1.8`$ (a line crossing regions A,B and C). The zero noise region ($`T=\beta ^1=0`$) for the phase diagram of figure 4 is shown in figure 5 (right picture) where we draw the transition lines separating recall regimes A ,B, C and F in the $`(J_{\mathrm{}}^{(2)},T)`$ plane for $`J_s^{(2)}=4`$.

In equations (61) and (69) we observe that, due to the explicit appearance of the variable $`p`$ in the solution of the overlap order parameter (which is due to short-range interactions, originating from expressions of the form $`\mathrm{exp}[\beta (J_s^{(1)}+J_s^{(2)}_{\mu =1}^p\xi _i^\mu \xi _{i+1}^\mu )]`$, see e.g. (18)), it will no longer be true that the pure state ansatz leads to solutions which are independent of the number of stored patterns, as is the case for standard mean-field Hopfield networks. This is also shown in the phase diagrams of figure 6, which have been constructed from (61,69) with $`J_s^{(1)}=J_{\mathrm{}}^{(1)}=0`$ and for $`p=15`$. We observe a significant increase in the number of transition lines, as well as additional transition lines appearing in the quadrant $`J_{\mathrm{}}^{(2)},J_s^{(2)}>0`$ (and also in $`J_{\mathrm{}}^{(2)},J_s^{(2)}<0`$ for the parallel case). Such effects become more and more prominent as the number of patterns increases. It is also worth noting that figures 4 and 6 imply that all ‘pairs’ of first-order transition lines point to the origin of the $`\{J_s^{(2)},J_{\mathrm{}}^{(2)}\}`$ plane.

Finally, due to the occurrence of imaginary saddle-points in (57) and our strategy to eliminate the variable $`\widehat{m}`$ by using the equation $`_m\varphi (m,\widehat{m})=0`$, it need not be true that the saddle-point with the lowest value of $`\varphi (m,\widehat{m})`$ is the minimum of $`\varphi `$ (complex conjugation can induce curvature sign changes, and in addition the minimum could occur at boundaries or as special limits). To remove this uncertainty we have evaluated the sequential dynamics free energy $`\varphi _{\mathrm{seq}}(m,\widehat{m})`$ (16) after elimination of the conjugate variable $`\widehat{m}`$ by using $`_m\varphi (m,\widehat{m})=0`$ (middle row of figure 7) and $`_{\widehat{m}}\varphi (m,\widehat{m})=0`$ (lower row of figure 7) which shows that, although the convexity of the free energy graphs is indeed affected, the location of the minima is not.

## 5 Benchmark Tests

We now compare our results with simple benchmark cases. First, the solution of model (4) should reduce to the model of Amit et al , for regimes where $`pN`$, upon removing short-range connectivity, i.e. for $`J_s^{(1)}=J_s^{(2)}=0`$. Indeed we find that in this limit the probability distributions $`P_{\mathrm{}}^{\mathrm{seq}}(k)`$ and $`P_{\mathrm{}}^{\mathrm{par}}(k)`$ of the stochastic variables (4) and (4) both reduce to the delta peak: $`\delta [k1]`$. This simplifies the solution of our problem and allows us to write for the free energies:

$`\varphi _{\mathrm{seq}}(m,\widehat{m})`$ $`=`$ $`im_0\widehat{m}_0im\widehat{m}m_0\theta {\displaystyle \frac{1}{2}}J_{\mathrm{}}^{(1)}m_0^2+{\displaystyle \frac{1}{2}}J_{\mathrm{}}^{(2)}m^2{\displaystyle \frac{1}{\beta }}\mathrm{log}2\mathrm{cosh}[\beta i(\widehat{m}_0+\widehat{m}\xi )]_\xi `$
$`\varphi _{\mathrm{par}}(m,\widehat{m})`$ $`=`$ $`im_0\widehat{m}_0im\widehat{m}m_0\theta `$
$`{\displaystyle \frac{1}{\beta }}\mathrm{log}2\mathrm{cosh}[\beta i(\widehat{m}_0+\widehat{m}\xi )]_\xi {\displaystyle \frac{1}{\beta }}\mathrm{log}2\mathrm{cosh}[\beta (\theta +J_{\mathrm{}}^{(1)}m_0+J_{\mathrm{}}^{(2)}m\xi )]_\xi `$
Simple differentiation with respect to $`\{m_0,\widehat{m}_0,m,\widehat{m}\}`$ verifies that in the mean-field limit the pure order parameter solutions reduce to $`m=\xi \mathrm{tanh}[\beta (\theta +J_{\mathrm{}}^{(1)}m_0+\xi J_{\mathrm{}}^{(2)}m)]_\xi `$ and $`m_0=\mathrm{tanh}[\beta (\theta +J_{\mathrm{}}^{(1)}m_0+\xi J_{\mathrm{}}^{(2)}m)]_\xi `$ as they should.

Our second benchmark test is provided by the exact solution of model (5), section 3. We can immediately map model (4) to model (5) by setting $`J_{\mathrm{}}^{(2)}=J_s^{(2)}=0`$. We then find that the key variables $`\{k_j\}`$ of equations (41-42) are given by a simple deterministic map. In fact, for both sequential and parallel dynamics we find that $`k_{\mathrm{seq}}`$ and $`k_{\mathrm{par}}`$ evolve towards the same fixed point:

$$k=e^{2\beta J_s}e^{\beta i\widehat{m}}\left[\mathrm{sinh}[\beta i\widehat{m}]+\sqrt{\mathrm{sinh}^2[\beta i\widehat{m}]+e^{4\beta J_s}}\right]$$

which (at the relevant saddle points) can be verified to lead to the fixed-point equation (9) of model (5), namely

$$m=G(m;J_{\mathrm{}},J_s)\mathrm{with}G(m;J_{\mathrm{}},J_s)=\frac{\mathrm{sinh}[\beta J_{\mathrm{}}m]}{\sqrt{\mathrm{sinh}^2[\beta J_{\mathrm{}}m]+e^{4\beta J_s}}}$$

Thirdly, we have also compared the free energies of model (4), as given by our present solution, to that which one finds when using the alternative random-field technique of . The latter relies on performing the spin summations in $`R=_𝝈(𝝈)`$ (18) and deriving appropriate functions $`A(𝝃𝝃^{})`$ and $`B(𝝃𝝃^{})`$ such that the identity $`\mathrm{cosh}[\beta (J_s(𝝃𝝃^{})\sigma ^{}i\widehat{m}\xi )]=\mathrm{exp}[\beta (A(𝝃𝝃^{})\sigma ^{}+B(𝝃𝝃^{}))]`$ is true for $`\sigma ^{}\{1,1\}`$. For instance, for the expression (18) of sequential dynamics this leads to:

$$\underset{N\mathrm{}}{lim}\frac{1}{\beta N}\mathrm{log}R_{\mathrm{seq}}=\underset{N\mathrm{}}{lim}\frac{1}{2\beta N}\underset{i=1}{\overset{N}{}}\mathrm{log}\left[4\mathrm{cosh}[\beta (J_{i,i+1}^s+h_i)]\mathrm{cosh}[\beta (J_{i,i+1}^sh_i)]\right]$$

$$\mathrm{where}h_{i+1}=i\widehat{m}\xi _{i+1}\frac{1}{2}\mathrm{log}\left[\frac{\mathrm{cosh}[\beta (J_{i,i+1}^sh_i)]}{\mathrm{cosh}[\beta (J_{i,i+1}^s+h_i)]}\right]J_{lk}^s=J_s^{(1)}+J_s^{(2)}𝝃_l𝝃_k$$

(with $`h_1=i\widehat{m}\xi _1`$) and is in complete agreement with the free energy as found from (16,45,51).

Finally, for the special case $`J_s^{(1)}=J_{\mathrm{}}^{(1)}=J_{\mathrm{}}^{(2)}=0`$ and $`p=1`$ (short-range bond disorder and absence of long-range interactions) our model reduces to the classical short-range random-bond Ising model , in which we expect the integrated density $`\widehat{P}_{\mathrm{}}(k)`$ (53) to acquire, at least in certain parameter regions, the form of the highly non-analytic Devil’s Staircase . In this special case the density (53) reduces (at saddle points $`\{i\widehat{m}_0=\theta ,i\widehat{m}=0\}`$) to

$$\widehat{P}_{i+1}(k,\lambda )=\frac{1}{2}\left\{\widehat{P}_i(e^{2\beta \theta }\frac{ke^{\beta J_s\lambda }e^{\beta J_s\lambda }}{e^{\beta J_s\lambda }ke^{\beta J_s\lambda }},1)+\widehat{P}_i(e^{2\beta \theta }\frac{ke^{\beta J_s\lambda }e^{\beta J_s\lambda }}{e^{\beta J_s\lambda }ke^{\beta J_s\lambda }},1)\right\}$$

where $`\lambda =\pm 1`$ represents bond-disorder. One can now prove by induction that if the identity $`\widehat{P}_i(k,1)=\widehat{P}_i(k,1)`$ is true then also $`\widehat{P}_{i+1}(k,1)=\widehat{P}_{i+1}(k,1)`$, so that (assuming ergodicity and uniqueness of the stationary density) the above expression reduces to a single recursive equation:

$$\widehat{P}_{i+1}(k)=\frac{1}{2}\left\{\widehat{P}_i\left(e^{2\beta \theta }\frac{ke^{\beta J_s}e^{\beta J_s}}{e^{\beta J_s}ke^{\beta J_s}}\right)+\widehat{P}_i\left(e^{2\beta \theta }\frac{ke^{\beta J_s}e^{\beta J_s}}{e^{\beta J_s}ke^{\beta J_s}}\right)\right\}$$

which is recognised as equation 10 of Ref. , upon a simple re-definition of our stochastic variables: $`\{k_n=e^{2\beta i\widehat{m}_0}R_{n,}/R_{n,}\}\{e^{2\beta i\widehat{m}_0}R_{n,}/R_{n,}\}`$ in (39). In the present benchmark case we have also verified the identity found in  relating short-range random-field models between sequential and parallel dynamics, namely

$$\psi _{\mathrm{seq}}(\psi _{\mathrm{seq}}(k;\xi ^{\prime \prime }\xi ^{},\theta ),\xi ^{}\xi ,\theta )=\psi _{\mathrm{par}}(k;\xi ^{\prime \prime }\xi ^{},\xi ^{}\xi ,\theta )$$

(which is the key identity to prove that in the thermodynamic limit sequential and parallel random-field models lead to the same physical states). Here, $`\psi _{\mathrm{seq}}(\mathrm{})`$ and $`\psi _{\mathrm{par}}(\mathrm{})`$ correspond to the functions defined in (4) and (4).

## 6 Theory Vs Simulations

In order to test our results further, we have performed extensive simulation experiments of the process (1), for model (4). In all cases the initial state is prepared randomly, with non-zero correlation only with pattern $`\{\xi _i^1\}`$. Our simulation results for the model which gives the phase diagram of figure 4 ($`p=2`$ and $`J_{\mathrm{}}^{(1)}=J_s^{(1)}=0`$) are shown in figure 7 (upper row) where we draw the equilibrium value of the recall overlap $`m_1(t\mathrm{})`$ as a function of the initial state $`m_1(t=0)`$. We have performed our experiments for three different regions of the phase diagram: region B (one $`m>0`$ stable state), region C (two $`m>0`$ stable states) and region F (three $`m>0`$ stable states). In regions B and C the simulation experiments verify the appearance and location of multiple ergodic sectors; to compare with the theoretical results, see also the free energy graphs in the middle and lower pictures of figure 7. In region F the simulation experiments show that the system can enter only two possible domains of attraction (excluding thus the theoretically predicted state $`m0.32`$). This is due to ($`i`$) the system’s finite size ($`N=1,000`$), in combination with ($`ii`$) the (relatively) small energy barrier separating the two physical states $`m_a0.32`$ and $`m_b0.68`$ (this allows the system to move from state $`m_a`$ to $`m_b`$ with a non-negligible probability). Our restriction to system size $`N=1,000`$ is prompted by the extremely long equilibration times (of the order of $`𝒪(10^6)`$ flips/spin). This, in turn, is due to domain formation: large clusters of neurons tend to freeze in specific configurations. As a consequence, in order for neurons to flip, the entire domain has to flip. In figure 8 (left graph) we show the value of the condensed overlap as a function of time in region F. We see that, starting from an initial state $`m_1(0)0.09`$, the system gradually approaches the theoretically predicted locally stable state, where it indeed stays for a period of $`410^5`$ flips/spin. Due to finite size effects, however, this state is thermodynamically unstable. A sudden transition to a new meta-stable state is then observed, generated by the flipping of entire domains. Equilibrium is reached in these simulations at about $`210^6`$ flips/spin, where a second and final jump transition takes place. In the right graph of figure 8 we show a simulation experiment carried out in region C, starting from initial conditions $`m_1(0)0.08`$. Here equilibrium is reached after about $`10^6`$ flips/spin, and due to a (relatively) high energy barrier separating the two $`m>0`$ physical states (see energy graphs in regions C, left column of figure 7) there is no domain-related transition. In all our experiments the value of the non-selected pattern overlap $`m_2(t)`$ is found to remain zero (open diamond points), which justifies a posteriori our pure state ansatz.

## 7 Discussion

In this paper we have presented an exact equilibrium solution for a specific class of spatially structured Ising spin (attractor) neural network models, in which there is competition induced by the presence of two qualitatively different types of synaptic interactions: those operating only between nearest neighbours in a 1D chain (short-range), and those operating between any pair of neurons (long-range). The values taken by the interactions present are given by Hebbian-type rules, as in the more familiar mean-field attractor networks. We have solved these models by using a combination of mean- and random-field techniques for both sequential and parallel dynamics. As in the standard 1D RFIM-type models our expressions for the disorder-averaged free energy per neuron take the form of integrals over the distribution of a random variable, which represents the ratio of conditioned partition functions. This distribution can be

evaluated numerically without much effort, and the key macroscopic observables then follow via simple differentiation.

We found that there are regions in parameter space where information processing between the two types of synaptic interactions can induce phenomena which are quite novel in the arena of associative memory models, such as the appearance of multiple locally stable states, and of first-order transitions between them, even for finite $`p`$ and upon making the ‘pure state’ ansatz. These peculiarities come to life particularly in regions where short- and long-range synapses compete most strongly, for instance, where one has Hebbian long-range interactions in combination with anti-Hebbian short-range ones, and they become more evident when increasing the number of stored patterns. Particularly in the upper left quadrant of parameter space $`\{J_{\mathrm{}}>0,J_s<0\}`$ one observes the appearance of an increasing number of dynamic transition lines (first- and second-order ones). This feature is in sharp contrast with the conventional (infinite-range) Hopfield-like networks , where for finite $`p`$ the ‘pure state’ ansatz automatically renders the remaining order parameter independent of the number of patterns stored. Phenomena such as simultaneous existence of multiple locally stable states (in which the quality of pattern recall depends crucially on initial conditions) can play a potentially useful role: poor cue signals can no longer evoke pattern recall. We also found that parallel dynamics transition lines in parameter space are exact reflection in the origin of those in sequential dynamics and that the relevant macroscopic observables can be obtained from those of sequential dynamics via simple transformations. Simulation experiments also show that the dynamics of the model are highly non-trivial, with plateaus and jump discontinuities, caused by complex domain formation and domain interaction, which would justify a study in itself.

In a similar fashion one could now also study more complicated systems, where (in addition to the long-range synapses) the short-range synapses reach beyond nearest neighbours. Such models can still be solved using the techniques employed here. A different type of generalisation would be to allow for a type of competition between synapses which would correspond to having stored patterns with different (pattern dependent) embedding strengths, as in . All these will be subjects of a future work.

### Acknowledgments

It is our pleasure to thank Jort van Mourik for useful suggestions on the simulations.


format: text

text:
# Cluster Coagulation and Growth Limited by Surface Interactions with Polymers

## Abstract

The physical and chemical properties of metal nanoparticles differ significantly from those of free metal atoms as well as from the properties of bulk metals, and therefore, they may be viewed as a transition regime between the two physical states. Within this nanosize regime, there is a wide fluctuation of properties, particularly chemical reactivity, as a function of the size, geometry, and electronic state of the metal nanoparticles. In recent years, great advancements have been made in the attempts to control and manipulate the growth of metal particles to pre-specified dimensions. One of the main synthetic methods utilized in this endeavor, is the capping of the growing clusters with a variety of molecules, e.g. polymers. In this paper we attempt to model such a process and show the relationship between the concentration of the polymer present in the system and the final metal particle size obtained. The theoretical behavior which we obtained is compared with experimental results for the cobalt-polystyrene system.

The physical and chemical properties of small clusters differ significantly from those of free atoms or molecules as well as from the properties of the bulk surfaces \[kn:mas1\]-\[kn:chered1\], and hence clusters may be viewed as a transition regime between the two physical states. Within the cluster regime, there is a wide fluctuation of properties, particularly chemical reactivity, depending on the size, geometry, electronic state and packing of the clusters \[kn:gatguc1\]-\[kn:parnie1\].

The motivating physical system for this work is the chemical synthesis of nanoscale cobalt clusters, via the thermal decomposition of cobalt carbonyl complexes in solutions containing either pure solvent or, in addition, polymers, which will result in a more viscous reaction environment. The chemical reactivity of metal clusters is largely determined by two main factors: (a) The size of the particles; and (b) The oxidation state of the metal clusters. The ability to control particle size would imply an ability to manipulate the reactivity and interfacial behavior of the clusters in our systems. The main factors governing the nucleation and growth of the metal particles during the decomposition process, are the metal cluster size; i.e., the number of reactive sites available for chemical reactions and the mobility of the reactive metal clusters and their ability to diffuse through the solution and collide with each other. Clearly, as particle size grows, the mobility of the particles decreases, and a preferred average pseudo-equilibrium size distribution is reached. If the medium in which the thermal decompositions take place is more viscous than a pure solvent, then it is expected that the upper bound on particle mobility will be reached faster, and hence the final average particle size will be smaller.

In this work, we develop a mathematical model which aims at elucidating the effect of the initial concentration of the polymer present in the reaction solution on the final cobalt particle size obtained via the thermal decomposition mechanism, and correlate the model to experimental results.

The general chemical synthesis of metal clusters has been designed to accomodate the manipulation of several process variables, and is best described in three stages (a) The preparation of homogeneous solutions of organometallic complexes in a carefully selected solvent; (b) The mixing of these solutions with a polymer solution resulting in a viscous medium which contains the homogeneously dispersed organometallic complexes; (c) The thermal decomposition of the organometallic complexes to form uniform dispersions of small particle size. The particular system studied in this work consists of cobalt carbonyl complexes thermally decomposed in the presence of polystyrene (PS) . It is important to note that experimentally, the solvent chosen for the decomposition reactions was toluene, which is a good solvent for polystyrene, and therefore it promoted the solvation of the polymer over a large compositional range.

The decomposition of metal carbonyls in various media has been reported in the literature for the preparation of fine metallic particles. The thermolysis of transition metal carbonyls in solution under an inert atmospheres is a well known technique for the preparation of pure metal powders \[kn:emeehr1\]-\[kn:tanfle1\]. Typically, the clusters formed following such thermal decomposition are $`50250`$ $`\AA `$ in diameter, and hence they belong to the mid-size cluster group, whose diameter lies within the $`0500`$ $`\AA `$ range \[kn:alihar1\]-\[kn:bru1\]. Some of the characteristic properties and advantages of these clusters prepared by the above method are: (a) The growth of the particles, and hence their final size, is determined mainly by diffusion distances and the driving forces between the growing clusters; (b) The chemical and physical control during synthesis of the clusters allows for the maintenance of a desired oxidation state and surface cleanliness; (c) The particles may be polycrystalline \[kn:reigol1, kn:tangol1\] and hence contain grain boundaries within the cluster. As a consequence, a large number of the atoms reside at the interface between these grain boundaries; (d) Another large fraction of the cluster’s atoms are on or close to the cluster’s surface, and hence these atoms do not order themselves identically to atoms within a ”bulk” material; and (e) There are various possibilities for the in situ reaction of these clusters with a wide range of particles or molecules, and moreover they are small enough to allow for a homogeneous dispersion within a polymer matrix.

Since during the decomposition reaction of metal carbonyls highly reactive intermediate species are formed \[kn:klatan1, kn:imikla1\], there are two major pathways for reaction available for these species: (a) They can aggregate or fragment to form larger or smaller clusters, according to the following general mechanisms respectively:

$`Co_j^{(0)}+Co_k^{(0)}Co_{j+k}^{(0)},`$ (1)
$`Co_{j+k}^{(0)}Co_j^{(0)}+Co_k^{(0)},`$
where $`j,k=1,2,\mathrm{}`$, and (b) They can interact with the polymer, according to the following general mechanism:

$$Co_j^{(0)}+\text{P}SCo\text{-}PS,$$

where $`j=1,2,\mathrm{}`$, resulting in metal attachment processes which may, in some cases, lead to polymer degradation and/or crosslinking. All the reactions described here are considered at $`90^oC`$ and in presence of $`N_2`$. The overall chemistry of the system will therefore be a combination of both pathways, and will be determined by the type of the reactive groups of the polymer in solution and the nature of the metallic species formed.

The overall prototype chemical reaction considered in this work is given by:

$$n\left[Co_x(CO)_y\right]+\text{PS}Co_m^{(0)}+\text{Co-PS}^{},$$
(2)
where  stands for complex compounds and $`x=1,2,4`$ or $`6`$ and $`y=4,8,12`$ or $`16`$, respectively, $`y/x4`$, and $`m=j+k`$ in (1). In the main reaction, the clusters are formed from smaller ones. The $`CoPS^{}`$ complexes formed during this process are the products of a parallel reaction which competes with the main coagulation-fragmentation process. The extent of the secondary reaction depends on the reactivity of the cobalt precursor and on the interaction coefficient between the reactive cobalt species and the polystyrene molecules. Therefore, the rate of the main coagulation reaction will depend not only on the concentration of the cobalt precursor and the rate coefficients of coagulation and fragmentation but also on the concentration of the available reactive sites on the polystyrene in the system which are capable of interacting with the cobalt precursor.

In order to translate the physical system considered here into mathematical terms, we first rescale the time variable multiplying it by the coefficient of reaction between cobalt fragments and polyestyrene, $`\lambda `$, a positive constant, and the initial cluster mass $`C`$, a positive constant. Then we define $`c_j`$ as the concentration of $`Co`$ clusters of size $`j`$ divided by $`C`$ and similarly, $`p`$ as the concentration of the reactive sites on the polystyrene, whose initial value is assumed to be a linear function of the polymer concentration, divided by $`C`$. For the model developed in this paper, we assumed that the only significant irreversible interaction between cobalt fragments and polystyrene will occur for $`c_1`$ with a coefficient of reaction $`\lambda `$, while for other clusters, $`j2`$, the interactions with the polystyrene will be reversible and weaker \[kn:cotwil1\]. Therefore, each $`c_1`$ particle which is involved in this interaction will be excluded from further participation in the main coagulation reaction. Moreover, the chemical bonding between the $`c_1`$ particles and the polystyrene reduces the number of reactive sites on the polymer and hence, limits the extent and rate of the cobalt-polystyrene interaction. It is important to stress that the overall chemistry of the system and the final particle size distribution of the cobalt clusters is dependent on the relative importance and contribution of the two competing reactions. The process can therefore be modelled by the following system of nonlinear differential equations \[kn:rotnov1\]:

$`\dot{c}_j={\displaystyle \frac{1}{2}}{\displaystyle \underset{k=1}{\overset{j1}{}}}[R_{jk,k}c_{jk}c_kQ_{jk,k}c_k]`$
$`{\displaystyle \underset{k=1}{\overset{Nj}{}}}[R_{j,k}c_jc_kQ_{j,k}c_{j+k}]\delta _{1,j}c_1p,`$ (3)
$`\dot{p}=c_1p,`$
for $`j=1,2,\mathrm{}`$, together with initial conditions satisfying $`_{j=1}^{\mathrm{}}jc_j(0)=1`$ and $`p(0)=p_0`$. In (3) $`p_0`$ is the initial dimensional concentration of reactive sites divided by $`C`$, $`R_{jk}`$ is the dimensional coefficient of coagulation divided by $`\lambda `$ and $`Q_{jk}`$ is the dimensional coefficient of fragmentation divided by $`\lambda C`$. All are non-negative and dimensionless. If we allow the formation of clusters of all sizes without any restriction, then $`N=\mathrm{}`$ in (3). Since, as we pointed out before, a pseudo-equilibrium state is reached as a consequence of the fragmentation process and the decrease in the mobility of large clusters, we assume that the concentration of clusters whose size is larger than some finite $`N`$ is negligible.

The process of metal coagulation and cluster formation (independent of the existence of polymer interactions) is a complex series of reactions which can be grouped into two categories: (1) The growth of clusters either by the formation of atom-atom bonds between two metal fragments and between a small cluster and a metal fragment, or by the development of surface-surface interactions between two large clusters. (2) The decomposition of large metal clusters into smaller fragments. In order to properly model these processes, it is important to construct the appropriate reaction rate expressions both for the cluster formation step and the cluster dissociation step (with or without polymer interaction).

The simplest coagulation expressions found in the literature which are pertinent for our model are: $`R_{jk}=j^\alpha k^\alpha `$ (predominantly chemical bond) \[kn:zif1\]-\[kn:ley2\], and $`R_{jk}=j^\alpha +k^\alpha `$ (predominantly surface interactions) \[kn:zif1, kn:zif5\], where in both cases $`0\alpha 1`$. In this model, we recognize the fact that there is a decline in particle mobility through the medium as a function of increased particle size. Moreover, the presence of the polystyrene molecules in the reaction solution contributes to the formation of a viscous medium which will impose an upper bound on the ability of larger clusters to interact via surface-surface interactions, and hence these types of interactions will have a negligible contribution to the overall coagulation process. These assumptions are consistent with the fact that in large clusters the fraction of atoms on the surface is small compared to the total number of atoms in the cluster, and therefore only this small fraction is essentially coordinatively unsaturated. In effect, we can view large clusters of being thermodynamically stable in real time. This is not the case for small clusters whose surface atoms represent a considerable fraction of the total metal mass. In this case, the model has to promote the growth of particles formed from smaller fragments, because this growth will minimize their free energy by reducing the number of coordinatively unsaturated surface atoms.

Therefore, the expression that we propose for the coefficients of coagulation of the cobalt fragments in our model (with or without polymer interactions) is

$$R_{j,k}=\frac{j+k}{\sigma }j^\alpha k^\alpha +\frac{\sigma }{j+k}\frac{\mu }{j^\beta k^\beta },$$
(4)
for $`j,k=1,2,\mathrm{}`$, where $`0\alpha ,\beta 1`$, $`\mu >0`$ is a weighting parameter that denotes the fraction of cobalt clusters whose coagulation contribution is negligible, and $`\sigma >0`$ denotes the critical size of the cobalt clusters beyond which the coagulation rates will drop as a function of size increase.

The dissociation of cobalt clusters into smaller particles has a direct impact on the final equilibrium cluster size distribution of the system. We assume that large clusters dissociate preferentially into smaller fragments of similar size rather than into fragments with a large variation in size. Therefore, the expression that we propose to describe the dissociation of large cobalt clusters into smaller fragments is

$$Q_{j,k}=\eta \frac{j^\gamma k^\gamma }{j^\gamma +k^\gamma }$$
(5)
for $`j,k=1,2,\mathrm{}`$, where: $`0\gamma 1`$, and $`\eta 0`$ is a weighting parameter that denotes the fraction of cobalt clusters which undergo dissociation.

In order to qualitatively test our model, system (3) has been solved numerically using a Runge-Kutta method of order four with 0.005 as time-step and for $`N=500`$, $`\alpha =0.1`$, $`\beta =0.1`$, $`\gamma =1`$, $`\eta =0.0001`$, $`\sigma =250`$ and $`\mu =1.0`$. These values of the parameters have been chosen so as the solutions reflect a pseudo-equilibrium situation in which there is a preferred average size or, in some cases, preferred average sizes. In Figure 1 we can see a graph of these preferred average sizes as a function of $`p_o`$.

Experimental results support the qualitative relationship between initial concentration of polymer in the reaction solution and the final particle average size. Experiments were performed by conducting the cobalt carbonyl decompositions in solutions with different polystyrene concentrations. A TEM microgrograph picture of cobalt particles obtained from a decomposition reaction in the absence of PS (solvent only) is shown in Figure 2. The average particle size obtained in this case is $`179.5\AA `$. Another TEM microgrograph picture of cobalt particles obtained from a decomposition reaction in a PS matrix ( 85 wt.% PS in solution) is shown in Figure 3. The average particle size obtained in this case is $`13.8\AA `$. Clearly, the presence of PS in the reaction medium at these extreme concentrations has a drastic effect on final average particle size, which decreases by more that an order of magnitude when the medium consists mainly of PS.

Figure 4 summarizes the experimental data obtained for media containing a variety of concentrations of PS. Not only does the average final particle size decrease with increasing initial concentration of PS, but the average cobalt size distribution decreases as well. This indicates that the presence of the polymer in the reaction medium prevents the formation of larger clusters thus inhibiting the coagulation process.