## 3 Two Coupled FHN Neural Systems with Time Delay  
In this section, in order to consider the effects of time delay in signal transmission between two neurons, we modify Equation (2) into:  
\[ \begin{array}{l}\dot{u}_{1}=-u_{1}(u_{1}-1)(u_{1}-a)-u_{2}+c(u_{1}-u_{3}(t- \tau_{1})),\\ \dot{u}_{2}=b(u_{1}-\gamma u_{2}),\\ \dot{u}_{3}=-u_{3}(u_{3}-1)(u_{3}-a)-u_{4}+c(u_{3}-u_{1}(t-\tau_{2})),\\ \dot{u}_{4}=b(u_{3}-\gamma u_{4}),\end{array}  \tag{12}\]  
where \(\tau_{1}(\tau_{2})\) represents the time delay in signal transmission from \(u_{3}(u_{1})\) to \(u_{1}(u_{3})\) and other parameters remain the same as that in Equation (2). In this paper, we place emphasis on presenting the influence of time delay in signal transmission between two neurons on synchronization conditions. We take \(\tau_{1}=\tau_{2}=\tau\) in Equation (12) to simplify the synchronization analysis. Similar to that of Equation (2), synchronization in Equation (12) is said to be achieved if:  
\[\lim_{t\rightarrow\infty}(u_{1}-u_{3})=0\ \,\ \ \lim_{t\rightarrow\infty}(u_{2}-u_{ 4})=0\ \.\]  
Substituting:  
\[e_{1}=\frac{u_{1}-u_{3}}{2},e_{2}=\frac{u_{1}+u_{3}}{2},e_{3}=\frac{u_{2}-u_{4 }}{2},e_{4}=\frac{u_{2}+u_{4}}{2},\]into Equation (12) leads to:  
\[ \begin{array}{l}\dot{e}_{1}=-\varepsilon_{1}^{3}+[-3\varepsilon_{2}^{2}+2(a+1)e _{2}+c-a]e_{1}+ce_{1}(t-\tau)-e_{3},\\ \dot{e}_{3}=b(e_{1}-\gamma e_{3}),\\ \dot{e}_{2}=-\varepsilon_{2}^{3}+(1+a)(\varepsilon_{1}^{2}+\varepsilon_{2}^{2} )-(3\varepsilon_{1}^{2}+a-c)e_{2}-ce_{2}(t-\tau)-e_{4},\\ \dot{e}_{4}=b(e_{2}-\gamma e_{4}).\end{array}  \tag{13}\]  
By eliminating \(e_{3}\) and \(e_{4}\), Equation (13) becomes:  
\[ \begin{array}{l}\ddot{e}_{1}+[3(\varepsilon_{1}^{2}+\varepsilon_{2}^{2})+ \alpha_{1}e_{2}+\alpha_{2}]\dot{e}_{1}+(6e_{2}+\alpha_{1})e_{1}\dot{e}_{2}+ \alpha_{3}e_{1}^{3}\\ \qquad+(3\alpha_{3}\varepsilon_{2}^{2}+2\alpha_{4}e_{2}+\alpha_{5})e_{1}- \alpha_{6}e_{1}(t-\tau)-\dot{c}\dot{e}_{1}(t-\tau)=0,\\ \ddot{e}_{2}+[3(\varepsilon_{1}^{2}+\varepsilon_{2}^{2})+\alpha_{1}e_{2}+ \alpha_{2}]\dot{e}_{2}+(6e_{2}+\alpha_{1})e_{1}\dot{e}_{2}+\alpha_{3}e_{2}^{3 }\\ \qquad+(3\alpha_{3}e_{1}^{2}+\alpha_{5})e_{2}+\alpha_{4}(\varepsilon_{1}^{2}+ \varepsilon_{2}^{2})+\alpha_{6}e_{2}(t-\tau)+\dot{c}\dot{e}_{2}(t-\tau)=0, \end{array}  \tag{14}\]  
where \(\alpha_{1,2,3,4,5,6}\), have been defined in Equation (3). We successfully applied the energy method to calculate periodic solutions of a delay differential equation [39], which involved the computations of the mechanical energy of oscillators. This illustrated that although the solution space of a delay differential system has infinite dimensions, approximate calculations of mechanical energy of oscillators are still feasible based on the energy method. Thus, we can use the energy method to analyze the mechanical energy variation of \(e_{1}\) in Equation (14). To this end, we redefine the following functions:  
\[ \begin{array}{l}g_{1}(e_{1})=\alpha_{3}e_{1}^{3},\ g_{2}(e_{2})\\ f_{1}\big{(}e,\dot{e},e(t-\tau),\dot{e}(t-\tau)\big{)}\\ \qquad+\big{(}3\alpha_{3}e_{2}^{2}+2\alpha_{4}e_{2}+\alpha_{5}\big{)}e_{1}- \alpha_{6}e_{1}(t-\tau)-\dot{c}\dot{c}_{1}(t-\tau),\\ f_{2}\big{(}e,\dot{e},e(t-\tau),\dot{e}(t-\tau)\big{)}\\ \qquad+\big{(}3\alpha_{3}e_{1}^{2}+\alpha_{5}\big{)}e_{1}+\alpha_{4}\big{(} \varepsilon_{1}^{2}+\varepsilon_{2}^{2}\big{)}+\alpha_{6}e_{2}(t-\tau)+\dot{c }\dot{e}_{2}(t-\tau).\end{array} \]  
The potential and mechanical energy of Equation (14) can be expressed by:  
\[ \begin{array}{l}V_{1}(e_{1})=\int_{0}^{e_{1}}{(\alpha_{3})e_{1}^{3}de_{1}}= \frac{1}{4}\alpha_{3}e_{1}^{4},\\ V_{2}(e_{2})=\int_{0}^{e_{2}}{(\alpha_{3})e_{2}^{3}de_{2}}=\frac{1}{4}\alpha_{3} e_{2}^{4},\\ E(e,\dot{e})=V(e)+\frac{1}{2}\dot{e}.\end{array} \]  
We define the transformation between Equation (14) and its energy coordinate \((E,\theta)\) as:  
\[