text:
# Magneto-Transport Properties of Doped RuSr2GdCu2O8

## Introduction

Soon after the first successful synthesis of RuSr2GdCu2O8 Bauernfeind1996  the material was found to display not only superconductivity ($`T_c`$45K) but coexisting magnetic order with $`T_{Curie}`$135KTallon1999 ; Bernhard1999 . Evidence accumulated from static magnetization, muon spin rotation ($`\mu `$SR)Bernhard1999  and from Gd-ESRFainstein1999  studies showing that the magnetism is a spatially uniform bulk property. Specific heat measurementsTallon1999  and the diamagnetic shielding fraction at low temperaturesBernhard1999 ; Bernhard2000 ; McCrone2000\_2  indicate that the superconductivity is also a bulk property, and that the two phases therefore coexist on a truly microscopic scale. An initial neutron diffraction study eliminated the possibility of ferromagnetic order with the Ru moments lying in the RuO2 plane, but did not rule out ferromagnetic alignment with the moments parallel to the $`c`$-axis, canted ferromagnetism or itinerant ferromagnetismChmaissem2000 . Subsequent polarized neutron diffraction dataLynn2000  have thrown the debate on RuSr2GdCu2O8 wide open by appearing to show that the underlying ordering of the Ru moments below the magnetic transition is in fact G-type antiferromagnetic (antiparallel nearest-neighbour ordering in all three crystallographic directions). Finally, the most recent neutron measurements on RuSr2YCu2O8 confirmed that there is indeed a ferromagnetic component of about 0.28 $`\mu _B`$ which is about 1/5th of the AF component of 1.2 $`\mu _B`$ Takigawa . The magnetic order shows a rather strong and unusual response to an applied magnetic field, with the FM component growing rapidly in strength and dominating over the AF already at 2T. Whatever the nature of its magnetism, the discovery of this material is an exciting development which poses new and important questions to our understanding of the interplay between magnetic and SC order.

Magnetoresistance (MR), Hall effect and thermopower (TEP) measurements on undoped sintered ceramic RuSr2GdCu2O8 were presented previouslyMcCrone1999 . Above $`T_{Magnetic}`$ the MR is negative and proportional to the square of the Ru magnetization and was ascribed to spin-scattering of the current carriers. A model for dilute magnetic alloys was used to extract from these data a value ($``$25 meV) for the exchange interaction between the Ru moments and the carriers. Below $`T_{Magnetic}`$ the Hall effect and TEP both fall anomalously. It will be shown that these data may be explained within a simple two-band model by a transition from localized to more itinerant behavior in the RuO2 layer at $`T_{Magnetic}`$. Evidence for delocalized carriers within the RuO layers has also been obtained from other transport and microwave absorption studiesWilliams2002  as well as from Ru-NMR measurements where clear anomalies in the Ru-NMR relaxation rate occur near Tc Tokunaga . This observation suggests that the Ru nuclear moments experience a sizeable hyperfine coupling to the charge carriers that enter the SC state.

The magneto-thermopower reveals an extremely unusual variation of $`T_c`$ with applied fieldMcCrone2000\_2 : $`T_c`$ actually increases by $`4K`$ as the applied field is increased to 2T. The increase saturates along with the Ru magnetization, suggesting that the onset of Ru magnetic order reduces a magnetic pair-breaking effect in the CuO2 layer.

The carrier concentration in RuSr2GdCu2O8 and its magnetic and SC properties, structural deformations, and so forth may be altered by cation substitution. Examining the transport properties of such samples should lead to a better understanding of the parent material. In this paper we present magneto-transport measurements on substituted RuSr2GdCu2O8. It will be shown that the data strongly support a simple two-band model in which the Hall effect and TEP of each sample are determined by the properties of the CuO2 and RuO2 layers, weighted appropriately by their conductivities. The model indicates that the RuO2 layer in the undoped material is very poorly conducting at room temperature, with $`\sigma _{Ru}0.1\sigma _{Cu}`$, increasing to $`0.3\sigma _{Cu}`$ or higher at low temperature. While in most of the samples studied the CuO2 layer remains the better conductor at all temperatures, we find that the RuO2 layer dominates the conductivity below $`T_{Magnetic}`$ in a sample with 10% Ce4+ substituted for Gd3+.

## Experimental Methods

Phase-pure sintered pellets of RuSr2GdCu2O8 were synthesized as described previously via solid-state reaction of a stoichiometric mixture of high-purity metal oxides and SrCO3Tallon1999\_2 ; Bernhard1999 . The doped samples, listed in Table I, were produced similarly; the compositions given are nominal. A final extended anneal at 1060C in flowing high-purity O2 produces a marked improvement in the crystallinity of the undoped material, resulting in a higher resistive $`T_c`$ (as defined by $`\rho (T)=0`$) but no significant change in the thermodynamic TcTallon1999 .

Bars of approximate dimensions 4$`\times `$1$`\times `$0.7mm3 were cut from the sintered pellets using a diamond wheel, then polished down to a thickness of $`150\mu `$m in order to increase the measured Hall voltage. They were mounted on quartz substrates in a standard six-contact configuration allowing both resistance and Hall voltage to be measured simultaneously. The contacts were made using 25$`\mu `$m gold wire and Dupont 6838 conducting epoxy, cured in air at 450C for six minutes, giving contact resistances $`<1\mathrm{\Omega }`$.

Resistivity and Hall effect measurements were made using an ac current source, low-noise transformers and lock-in amplifiers. A frequency of $`77`$Hz was used to avoid mains pick-up, with current densities of around 0.25Acm-2. The Hall coefficient, $`R_H`$, was usually measured by stabilizing the temperature and field (10T unless stated otherwise), then measuring the Hall voltage with the sample rotated by 0 and 180 with respect to the field. The Hall coefficient is then given by $`R_H=\frac{\left(V_0V_{180}\right)t}{2IB}`$ where B is the magnetic field, t the sample thickness and I the current. This method eliminates the MR of the sample, and the offset voltage from $`\rho _{xx}`$ due to contact misalignment. Where $`R_H`$ was measured as a function of field, this was swept to both positive and negative values and $`R_H(B)`$ determined from $`V_BV_B`$.

TEP measurements were made by the ‘toggled’ heating methodCooper1971 ; Resel1996 . Two 25$`\mu `$m chromel-alumel thermocouples, attached to the sample with small blobs of silver paint, measure both the thermal emf and temperature gradient, ensuring that these are measured between the same two points. The sample is first stabilized at the measurement temperature, a small thermal gradient is applied, and the resulting thermal emf measured. The thermal gradient is then reversed, allowing slowly changing thermal emfs in the cryostat wires to be nulled out. A ‘rest state’ was added whereby both ends of the sample were heated at half-power, providing two extra measurement points. Adding this state keeps the total power dissipation into the stage constant, avoiding fluctuation of its temperature when the heater currents are changed.

## Results

### Transport Measurements on Pure RuSr2GdCu2O8

Hall effect, thermopower and resistivity data for undoped RuSr2GdCu2O8 are shown in Fig. 1. The room-temperature value of the TEP implies a hole concentration, $`p_{Cu}`$, of 0.06-0.07 holes/CuObertelli1992 , while its temperature dependence is typical of other high-$`T_c`$ materials, with the exception of the unusual linear temperature dependence below $`T_{Magnetic}`$. The overall magnitude and temperature dependence of the Hall coefficient is consistent with a doping level, $`p_{Cu}`$, of $`0.07`$ holes/Cu, as inferred from the room-temperature TEP. $`R_H`$ displays a high-$`T_c`$-like temperature dependence well above $`T_{Magnetic}`$. However, below about 170K there is an anomalous downturn in $`R_H`$ which is not seen in typical high-$`T_c`$ data. The so-called ‘anomalous’ Hall effect observed in magnetic materials has been measured and discounted as the cause of this downturnMcCronethesis . Alternative possibilities are that it is due to charge delocalization in the RuO2 plane occurring near the magnetic transition, or to charge transfer into the CuO2 layers. It will be shown that a two-band model, with a localized to itinerant transition occurring at $`T_{Magnetic}`$ in the RuO2 layer, can explain both these and the TEP data.

#### The conductivity of the RuO2 layer

We now introduce a simple two-band model which successfully describes the features of the above data. The bands in this model are those formed by carriers in the Cu and Ru orbitals; the overall TEP and Hall effect are given by the sum of the CuO2 and RuO2 layer values, weighted by the layer conductivitiesHurd1972 ; Barnard :

$$R_H=\frac{R_H^{Ru}\left(\sigma _{xx}^{Ru}\right)^2+R_H^{Cu}\left(\sigma _{xx}^{Cu}\right)^2}{\left(\sigma _{xx}^{Ru}+\sigma _{xx}^{Cu}\right)^2}$$
(1)
$$S=\frac{S^{Ru}\sigma _{xx}^{Ru}+S^{Cu}\sigma _{xx}^{Cu}}{\sigma _{xx}^{Ru}+\sigma _{xx}^{Cu}}$$
(2)
With some reasonable estimates of the RuO2 and CuO2 layer properties, it is possible to use this model and the measured room-temperature Hall effect and TEP to put a limit on the conductivity of the RuO2 layer. To do this we assume that the Hall coefficient of the RuO2 layer is approximately zero (the maximum value observed in other two-dimensional Ru oxides studied to date is $`5\times 10^{10}`$m3C-1)Yang1999 ; Mackenzie1996 ; Perry2000 . With this assumption in Eqn. 1, the conductivity of the RuO2 layer may be estimated from

$$\begin{array}{c}\frac{\sigma _{xx}^{Ru}}{\sigma _{xx}^{Cu}}=\sqrt{\frac{R_H^{Cu}}{R_H}}1\hfill \\ \\ R_H^{Ru}R_H^{Cu};\sigma _{xx}^{Ru}<\sigma _{xx}^{Cu}\hfill \end{array}$$
(3)
The ratio of the Hall coefficient of the CuO2 layers $`R_H^{Cu}`$ to the measured value $`R_H`$ caused by the presence of the RuO2 layer, is hard to estimate due to the uncertain doping state in RuSr2GdCu2O8 and the spread of values of $`R_H^{Cu}`$, for a given doping level, in the literatureWuyts1996 ; Carrington1993 ; Fisher1994 . Given these uncertainties, a reasonable range of values of $`R_H^{Cu}/R_H`$ is 1-1.4, giving $`\sigma _{xx}^{Ru}`$ in the range $`00.18\sigma _{xx}^{Cu}`$. The summary of $`R_H^{Cu}`$ values in the review by Cooper and LoramCooper1996  would favor the low end of this range.

For this range of conductivity in the RuO2 layer, Equation 2 predicts that the measured net TEP lies some 0-8$`\mu `$VK-1 below the intrinsic CuO2 layer value, i.e. that $`75S_{290}^{Cu}83`$. It is very unlikely that $`S_{290}^{Cu}`$ lies in the upper half of this range: a value of $`S=83\mu `$VK-1 would imply an extremely small hole concentration for which a $`T_c`$ as high as 46K would be extraordinary.

Having placed a limit on the conductivity one can use a 2D model to determine $`k_Fl`$, the product of the Fermi wave-vector with the mean free path for the RuO2 layers. This quantity gives an indication as to whether the carriers are localized or itinerant and for a cylindrical Fermi surface may be written as

$$k_Fl=\sigma \frac{2\pi \mathrm{}c}{e^2}$$
(4)
where $`c`$ is the separation of the planes. Data in the literature for the ab-plane resistivity of under-doped YBa2Cu3O7-δ (YBCO) films and single crystals, with $`p0.07`$, give a consistent value of 1.2m$`\mathrm{\Omega }`$cm at room temperature Carrington1993 ; Sun1994 ; Jin1998 ; Carrington1992  giving $`k_Fl_{Cu}=1.3`$, near the limit of localization. In fact, in only slightly more under-doped samples one sees a semiconducting upturn at low temperatures. Given the range of ratios of $`\sigma _{xx}^{Ru}`$ to $`\sigma _{xx}^{Cu}`$ derived from the Hall effect, $`k_Fl_{Ru}=00.45`$ at room temperature. The TEP data suggest that the true value is at the low end of this range, indicating that the carriers in the RuO2 layers are at best very poorly metallic.

#### Temperature dependence of $`\sigma _{Ru}`$

Having established that the room-temperature conductivity of the RuO2 layer is close to zero, typical $`S(T)`$ and $`R_H(T)`$ data for high-$`T_c`$ superconductors may be scaled so that the room temperature values match those of RuSr2GdCu2O8. The differences below $`T_{Magnetic}`$ may then be used to follow $`\sigma _{Ru}`$ as a function of temperature.

Typical Hall effect data for the CuO2 layer have been taken from measurements on sintered Ca-doped YBa2Cu3O7-δ, while $`R_H^{Ru}`$ will be set to zero, its value in other RuO2 layer compounds being much lower than $`R_H^{Cu}`$Yang1999 ; Mackenzie1996 ; Perry2000 . Typical $`S^{Cu}`$ data are approximated by measurements on sintered YBa2Cu3O7-δ with $`\delta =0.53`$Cooper1996 , multiplied by 1.12 to match the high-temperature RuSr2GdCu2O8 data. Finally, $`S^{Ru}`$ is approximated by data measured on a sintered sample of SrRuO3, which displays a magnitude and temperature dependence similar to that of Sr2RuO4Yoshino1996 . All these data are shown in Figure 1, together with the resulting $`\beta (T)=\frac{\sigma _{Ru}}{\sigma _{Cu}}`$ calculated from Eqns. 1 and 2.

Given the uncertainties in the approximated RuO2 and CuO2 layer properties the two $`\beta (T)`$ curves calculated independently from the drops in $`S(T)`$ and $`R_H(T)`$ agree well qualitatively. If TEP data for a sample of 20% Sn-doped RuSr2GdCu2O8, in which we shall argue that $`\sigma _{Ru}`$ is strongly suppressed below $`T_{Magnetic}`$, is used to approximate $`S^{Cu}`$, the agreement is also quantitative. Because the TEP is a less sensitive function of $`\beta `$ than the Hall effect, the difference between $`S^{1212}`$ and $`S^{Cu}`$ is quite small compared with that between $`R_H^{1212}`$ and $`R_H^{Cu}`$. Thus the value of $`\beta `$ calculated from the TEP data is more sensitive to inaccuracy in the assumed $`S^{Cu}`$ data. This explains why using the (only slightly different) 20% Sn-doped data to approximate $`S^{Cu}(T)`$ results in a better match to $`\beta (T)`$ calculated from the Hall effect. Whichever data are used, the results show a rapid rise in the relative conductivity of the RuO2 layer below 150K, to $`0.3\sigma _{Cu}`$ or higher.

### Transport Measurements on Substituted RuSr2GdCu2O8

#### Sn-Doped RuSr2GdCu2O8

The diamagnetic Sn4+ ion substitutes for Ru in solid solution, and is slightly larger in size than Ru4+/5+. The effects of doping the Ru site are of extreme interest given the current debate regarding the spin and charge configuration of the Ru ionsWilliams2000 ; Butera2000 ; Liu2000 .

We note that the Sn-doped samples studied here were from two sources prepared with slightly different annealing strategies. Comparison of their sample resistivities is therefore not necessarily meaningful, as annealing strongly affects the grain-boundary conductivity of RuSr2GdCu2O8McCrone1999 . In general the resistivity of sintered high-$`T_c`$ materials is also affected by sample densityCooper1991 ; Fisher1996 .

For the 2.5% and 7.5% samples the resistivity (Figure 2) is metallic, and similar in magnitude to the undoped sample. The 20% sample has a higher resistivity and shows a small semiconducting upturn at low temperatures, while both the magnitude and upturn are far larger for the 40% sample. Estimating $`T_c`$ from the onset of the resistive transition reveals a gradual increase from 40.5K for the 2.5% sample to 43.5K for the 20% sample, while the 40% sample has a reduced $`T_c`$ of just 30K.

The TEP, $`S(T)`$, and Hall effect, $`R_H(T)`$, are much less affected by grain-boundaries than the resistivity. In conventional high-Tc materials they closely reflect bulk CuO2 layer properties in conventional high-$`T_c`$ materialsCarrington1994 ; Fisher1996 .

The Hall effect data show a slow and monotonic decrease in $`T_c`$ with increasing Sn concentration, but it should be remembered that these data were taken in a field of 10T and only partly reflect the zero-field $`T_c`$. The vertical arrows in Figure 2 (and in subsequent figures) show the location of the magnetic transition. The TEP data show that $`T_c`$ (defined by the maximum in the derivative) rises by $``$4K in going from the 2.5% sample to the 20% sample, in good agreement with the resistivity data. The 40% sample shows a much lower transition temperature, both in $`R_H`$ and $`S`$. The increase in $`T_c`$ with Sn concentration is attributed to a transfer of holes into the CuO2 layerMaclaughlin1999\_2 , though we observe a smaller increase than the $``$12K reported previouslyMaclaughlin1999\_2 ; Maclaughlin2000 . In the earlier studies $`T_c`$ was defined from the resistivity onset, and the $`T_c`$ values obtained for low doping levels were significantly lower, possibly due to granularity.

On examining the temperature and doping dependence of the normal-state properties, one immediately observes that the room temperature TEP $`S_{290}`$ is little changed by the addition of Sn. This result is strange given the rise of $`4`$K in $`T_c`$ as the doping level is increased to 20%. The change in the Hall effect is also counter-intuitive: the 30% increase in going from 0 to 20% Sn would normally indicate a decrease in hole concentration. This apparent paradox is resolved when it is noticed that the anomalous drop in $`R_H`$ below $`T_{Magnetic}`$ is diminished in the 2.5% and 7.5% samples, and is absent in the 20% sample: as the Sn concentration is increased the RuO2 layer becomes well localized below $`T_{Magnetic}`$, reflecting significantly reduced conductivity at all temperatures. The changes in $`R_H(290)`$ and $`S_{290}`$ may then be explained quite simply: the introduction of Sn dopes a few extra holes into the CuO2 layer, increasing $`p_{Cu}`$ and raising $`T_c`$ by $``$4K, but also drives the RuO2 layer more insulating. Thus while $`R_H^{Cu}`$ probably decreases slightly, the overall Hall effect increases as the RuO2 layer no longer provides a parallel conduction pathway. The slight increase in $`p_{Cu}`$ which would normally decrease the measured TEP, is balanced by the decreasing $`\sigma _{Ru}`$, which removes the suppression of the TEP by the RuO2 layer, leaving it relatively unchanged overall. Certainly, the increase in doping is far smaller than one would expect from substituting Sn4+ for Ru5+, suggesting that the mean valency of the Ru ion is less than 5+. This conclusion is supported by recent XANES measurementsLiu2000 .

The reduction in the room-temperature $`R_H`$ of pure RuSr2GdCu2O8, due to conductivity in the RuO2 layer, was estimated to be of the order of 30%. This is entirely consistent with the rise in $`R_H`$ observed as the Sn concentration is increased to 20%, assuming that $`\sigma _{Ru}0`$. The 40% Sn-doped sample does not fit well into this picture, having a much larger $`R_H`$ at all temperatures. Given the much larger resistivity of this sample and its drastically reduced $`T_c`$, it is possible that some Sn$``$Cu substitution has occurred, reducing the CuO2 layer doping state, or that there are significant impurities present.

#### Nb-Doped RuSr2GdCu2O8

Nb also substitutes for Ru in the RuSr2GdCu2O8 structure, but has a dramatically different effect on the transport properties. In contrast to the Sn ion, which has a charge of 4+, Nb is believed to substitute in its usual 5+ stateMaclaughlin2000 , and thus for an average Ru valency of less than five will remove holes from the system, further under-doping it. The room-temperature TEP bears this out, showing a large increase proportionate with Nb doping (see Figure 3) and confirming that the CuO2 layer is progressively under-doped by the substitution of Nb. This conclusion is supported by the commensurate increase in the Hall effect and the rapid reduction of $`T_c`$, which is 19K for the 10% sample and below 1.5K (if present at all) in the 20% sample.

The effect of Nb doping on the Ru layer is less clear. The Hall effect of the 10% sample displays a maximum near $`T_{Magnetic}`$, suggesting increased RuO2 conductivity below this temperature, but the drop is not as clear as in the undoped sample. For 20% Nb$``$Ru substitution $`R_H`$ rises to a more or less constant value of 1.2$`\times 10^8`$m3C-1 below 200K, and there is no sign of a significant change at $`T_{Magnetic}`$. The conclusion from the TEP and Hall effect data, then, is that the transition from localized to itinerant behavior of the RuO2 layer is suppressed by the addition of Nb, as it is by the addition of Sn.

The resistivity, on the other hand, shows surprisingly little difference between the 10% and 20% samples – in fact the residual resistivity (extrapolated from the linear high-temperature data) actually decreases. A possible scenario consistent with this result is that the RuO2 layer becomes more itinerant both above and below $`T_{Magnetic}`$ as the Nb level is increased. However, if this were the case, the increased $`\sigma _{Ru}`$ would be expected to suppress both $`R_H`$ and $`S`$ below the CuO2 plane values. In fact, for $`T_c`$=19K and $`T_c^{max}100`$K, the universal relationship between $`S_{290}`$ and $`T_c`$ predicts $`S_{290}^{Cu}100\mu `$VK-1, as observed. Thus, while the increased $`S`$ and $`R_H`$, and the reduced $`T_c`$ are consistent with a reduced hole concentration in the CuO2 layer and a localized RuO2 layer, the relatively good conductivity of the 20% Nb-doped sample is not. One possible explanation is that the behavior of the resistivity is extrinsic to the bulk in the 20% Nb sample, resulting from either increased grain-boundary conductivity, or increased sample density.

#### Ce-Doped RuSr2GdCu2O8

Unlike Nb and Sn, which substitute for Ru, Ce substitutes for Gd in the layer separating the two CuO2 planes, and so would be expected to affect these more than the RuO2 layers from which it is relatively remote. The Ce ion is expected to be in the 4+ state in RuSr2GdCu2O8, as it is in the structurally similar compound RuSr2(Gd1+xCe1-x)Cu2O10Felner1999 ; hence its substitution for Gd3+ should further underdope the material.

Two samples (10% and 20% Ce$``$Gd) were measured, and of all the doped samples studied these exhibit the most remarkable and revealing transport properties: a large drop in $`R_H`$ below $`T_{Magnetic}`$ (in fact becoming negative in the 10% sample below $``$30K), and a large TEP at room temperature which, like the Hall effect, drops very rapidly below $`T_{Magnetic}`$. These data are shown in Figure 4, along with the resistivities of the two samples.

We note first that, as with other electron doping substitutions (Ce for Gd, La for Sr and hydrogen doping), $`T_{Magnetic}`$ is driven upwards. This appears to reflect an increasing Ru4+ fraction. The 10% Ce sample will be dealt with first. As with the undoped sample, the departure from cuprate-like properties below $`T_{Magnetic}`$ indicates a transition from localized to itinerant behavior in the RuO2 layer. In this case, however, the room temperature TEP, $`S_{290}=110\mu `$VK-1, indicates a much lower CuO2 layer carrier concentration of $`p0.03`$ holes/Cu. This is consistent with the increased Hall coefficient, which is probably still depressed from the true CuO2 value by residual conductivity in the RuO2 layer, and the large resistivity with its insulating upturn at low temperature. Having concluded that $`p`$, and hence $`\sigma _{Cu}`$, is much lower than in the undoped sample, the reason for the dramatic effects seen in $`R_H`$ and $`S`$ below $`T_{Magnetic}`$ becomes clear: the ratio $`\sigma _{Ru}/\sigma _{Cu}`$ is much larger in the Ce-doped sample at low temperature, allowing the intrinsic RuO2 layer properties to dominate the behavior.

The effect of changes in $`\sigma _{Ru}/\sigma _{Cu}`$ is greater for $`R_H`$ than $`S`$, but for these samples the increased $`\sigma _{Ru}`$ depresses the $`S_{Cu}`$ contribution to the total TEP so much that $`S_{Ru}`$ dominates below 50K. The basis for this assertion is the double maximum in $`\frac{dS}{dT}`$: initially, at low temperatures, $`S(T)`$ follows a curve reasonably consistent with the TEP of SrRuO3. This contribution would appear to saturate at a value of $`40\mu `$VK-1; however, above 50K increasing $`\sigma _{Cu}`$ allows $`S_{Cu}`$ to contribute, and the overall TEP then rises more rapidly.

The same qualitative treatment may be applied successfully to the Hall effect data, though in order to explain the negative values below $``$30K it is necessary to assume a negative Hall coefficient for the RuO2 layer of around $`1\times 10^9`$m3C-1. Examining typical data from the Srn+1RunO3n+1 series one finds that $`R_H`$ of Sr3Ru2O7 remains positive at all temperatures, while that of Sr2RuO4 becomes negative below 20K, but reaches just $`1\times 10^{10}`$m3C-1 near 1K. However, SrRuO3, which has the most similar ferromagnetic RuO2 layer to RuSr2GdCu2O8, has a negative $`R_H`$ below 100K, reaching a field-dependent value of $`1\times 10^9`$m3C-1 below 60KYang1999 . Thus the value of $`R_H`$ observed in the Ce-doped sample at low temperature is the same order of magnitude as that in SrRuO3, confirming that the RuO2 layer dominates the transport properties. It is interesting to note that, though it may not be a large effect, Ce substitution for Gd should drive the mean Ru valence closer to 4+, as it is in SrRuO3.

Turning now to the resistivity, one encounters a problem: if the RuO2 layer is indeed metallic below $`T_{Magnetic}`$, why does the resistivity increase so dramatically as $`T0`$? There are two possible answers to this question: either both the RuO2 and CuO2 layers are at least semiconducting, but such that $`\sigma _{Ru}/\sigma _{Cu}>1`$, or it may be that insulating grain boundaries cause the upturn. The second of these scenarios seems more likely. In this case the TEP and Hall effect, being much less sensitive to inter-grain connectivity, are determined by a weakly metallic intrinsic $`\sigma _{Ru}`$. Support for this conclusion is provided by close examination of the resistivity (Figure 5) which shows an extended metallic region below $`T_{Magnetic}`$.

This type of behavior is not uncommon in RuSr2GdCu2O8 – in fact extrinsic upturns in resistivity are observed in poorly annealed undoped samples. Interestingly though, transport measurements on SrRuO3 also show a minimum in resistivity below TCurie in samples where there is some disorder in the RuO2 layerKlein1996 . The temperature at which the minimum occurs, and the magnitude of the upturn below it both increase with RuO2 layer disorder: in good quality films the highest temperature minimum observed is 40K, coincident with the maximum residual resistivityKlein1996 . As the RuO2 layer may be considerably disordered in these doped RuSr2GdCu2O8 samples, it is possible that the behavior shown in Figure 5 is intrinsic to the RuO2 layer, and not a product of grain boundary resistivity. An extended anneal of the Ce-doped sample could be used to test the origin of the effect, as it should improve the quality of the grain boundaries in the material.

The 20% Ce-doped sample appears initially to be inconsistent with the interpretation of the 10% Ce data presented above: $`R_H`$ is higher at room temperature, as one would expect for even greater under-doping caused by the increase in Ce content, but $`S_{290}`$ is actually lower than that of the 10% sample, apparently implying an increased hole concentration. The paradox may be resolved straightforwardly if one accepts that $`\sigma _{Ru}\sigma _{Cu}`$, even at room temperature, in this sample. Then the values of $`S_{290}`$ and $`R_H(290)`$ for the material are no longer dominated by the CuO2 layer; $`S_{290}^{Cu}`$ and $`R_H^{Cu}(290)`$ may both in fact be much larger than the overall measured values so long as $`S_{290}^{Ru}`$ and $`R_H^{Ru}(290)`$ are lower, as expected. From the peak in the Hall effect data, $`T_{Magnetic}`$ of this sample appears to be depressed by about 30K from its value in the undoped sample, while the resistivity at 50K is more than two orders of magnitude greater, and the sample is insulating, probably due to extreme granularity.

#### Calculation of $`\beta (T)`$ in Ce-doped RuSr2GdCu2O8

The ratio $`\beta (T)=\sigma _{Ru}(T)/\sigma _{Cu}(T)`$ may be extracted from the data for the 10% Ce doped sample using the two-band model, as described for the undoped material in the previous Chapter. As in the previous calculation, typical $`S^{Cu}(T)`$ and $`R_H^{Cu}(T)`$ data are matched to the high temperature RuSr2GdCu2O8 data, where $`\sigma _{Ru}`$ is assumed to be small compared with $`\sigma _{Cu}`$ and the overall properties reflect those of the CuO2 layer most strongly. The deviation from cuprate-like behavior at lower temperatures is then used to extract the ratio $`\beta (T)`$.

For this sample, $`S^{Cu}(T)`$ data were taken as 1.05 times $`S(T)`$ measured on a sample of under-doped sintered YBa2Cu3O7-δ, with $`\delta 0.6`$Cooper1996 . $`R_H^{Cu}(T)`$ data were taken as 1.16 times $`R_H`$ measured on a similar sample with $`\delta 0.62`$Carringtonthesis . $`R_H`$ and $`S`$ are particularly strong functions of doping in this region of the phase diagram: the good agreement in the values of $`\delta `$ required for the two sets of data to match those of RuSr2GdCu2O8 suggests that the assumption of negligible $`\sigma _{Ru}`$ is reasonable.

The TEP of the RuO2 layer is approximated by that of sintered SrRuO3, as before. As the Hall effect becomes negative at low temperatures in Ce-doped RuSr2GdCu2O8 taking $`R_H^{Ru}0`$, as was done for the undoped material, will not work. Instead a rough approximation to data for SrRuO3 is usedYang1999 , which shows a field-dependent value of $`1\times 10^9`$ at 20K.

The measured and estimated data, together with the results of the calculations are shown in Figure 6. Above 50K there is remarkable agreement between $`\beta (T)`$ calculated from the TEP data ($`\beta _{TEP}`$) and that calculated independently from the Hall effect data ($`\beta _{Hall}`$), lending confidence both to the model and to the estimated $`R_H(T)`$ and $`S(T)`$ data for the RuO2 and CuO2 layers. Below 50K the agreement is not so good: $`\beta _{Hall}`$ carries on increasing, a direct result of $`R_H^{1212}`$ becoming very close to the estimated $`R_H^{Ru}`$ at low temperatures. $`S^{1212}`$ does not approach the estimated $`S^{Ru}`$ as closely, and hence $`\beta _{TEP}`$ does not continue to increase. Emerging clearly from these data is a large increase in $`\beta `$ below $`T_{Magnetic}`$. At 50K $`\sigma _{Ru}/\sigma _{Cu}1.9`$, whereas for the undoped material the increase in $`\sigma _{Ru}/\sigma _{Cu}`$ is just 0.3. The properties of the RuO2 layer dominate the overall transport properties of RuSr2GdCu2O8 below $`90`$K in this 10% Ce-doped sample.

#### Other Doped Samples

The remainder of the doped samples studied contained Y, Dy and Eu on the Gd site, plus a 5% Li-doped sample, in which Cu is substituted. The transport data for all these samples are shown in Figure 7.

The Hall effect shows the ‘usual’ anomalous downturn below $`T_{Magnetic}`$ in all these samples. The magnitude of the downturn, due to the transition to a more itinerant Ru layer, is approximately constant, leading to the conclusion that doping the Cu and Gd sites does not greatly affect the localization of carriers in the RuO2 layer.

Substituting a small amount of Li+ for Cu2+ causes virtually no change in the TEP, but depresses $`T_c`$ by approximately 20K. The Hall effect of this sample is slightly larger than that of the undoped sample, possibly due to some cross-substitution of Li with Ru, depressing $`\sigma _{Ru}`$, or a slight decrease in the CuO2 layer carrier concentration. These results are consistent with Li+ acting as a pair-breaker in the CuO2 layer, but with little other effect on transport properties. The rate of suppression of $`T_c`$ with Li substitution in RuSr2GdCu2O8, $`4`$K/%, is about one quarter of that observed in under-doped YBa2Cu3O7-δ ($`\delta =0.4`$) when either Li or Zn is substituted for Cu Bobroff1999 . However the concentration of Li in the RuSr2GdCu2O8 sample studied is nominal, and the difference in the rate of suppression may simply reflect loss of Li by vaporization during the long synthesis and anneal.

The isovalent substitution of Y or Dy for Gd actually causes a slight decrease in the TEP of RuSr2GdCu2O8, these being the only substitutions studied to do so. The implied increase in the doping level of the CuO2 layers, presumed to arise from an ion-size effect, is confirmed by the increased $`T_c`$ in these samples - approximately 2K (10% Y) and 6K (40% Dy) higher than in the undoped sample as seen by both resistivity and TEP measurements. However, the magnitude of the Hall effect is larger for these samples than for the undoped ones. Having argued that the CuO2 layer is less under-doped in these two samples, this effect may only arise from a decrease in the conductivity of the RuO2 layer, partially removing the ‘shorting’ of the CuO2 layer Hall effect.

Full substitution of Eu for Gd causes an increase in $`S_{290}`$ to 90$`\mu `$VK-1, $`T_c`$ as measured by the TEP or resistivity drops significantly, and $`R_H`$ is greatly increased. All these results suggest a drop in the CuO2 layer hole concentration, again consistent with the above-noted ion-size effect, perhaps coupled with a decrease in $`\sigma _{Ru}`$. The resistivity of this sample is much less metallic than that of the others (which are all metallic, with magnitudes two to three times the well-annealed undoped sample), consistent with this interpretation.

## Discussion

The Hall effect and TEP data described in this paper provide strong evidence for a transition from very poorly metallic to more itinerant behavior in the RuO2 layer below $`T_{Magnetic}`$. Results from substituted RuSr2GdCu2O8 samples confirm this picture. The universal relationship between $`S_{290}`$ and $`p_{Cu}`$ appears to hold in RuSr2GdCu2O8 as a result of the low $`\sigma _{Ru}`$ at room temperature, though below $`T_{Magnetic}`$ both $`S`$ and $`R_H`$ are reduced by ‘shorting’ of $`S_{Cu}`$ by the RuO2 layer.

The two-band model proposed is successful in explaining most of the existing data qualitatively: the anomalies, which lie mainly in resistivity data, are most likely due to grain boundary and density effects. The quantitative agreement is also reasonably good. The results support a picture in which the RuO2 layer in the pure compound is localized above $`T_{Magnetic}`$, with $`\sigma _{Ru}0.1\sigma _{Cu}`$, but becomes more conducting below $`T_{Magnetic}`$, mirroring the behavior of other ruthenates.

The transition from localized to itinerant RuO2 layer behavior at $`T_{Magnetic}`$ in the undoped compound may be modified by substituting Ru with Sn or Nb. Sn increases the doping level of the CuO2 layers, raising $`T_c`$ and suppressing $`T_{Magnetic}`$, and simultaneously drives the RuO2 layer more insulating. Nb under-dopes the CuO2 layers, lowering $`T_c`$, and also appears to drive the RuO2 layer insulating, though the 20% sample does not show the expected semiconducting resistivity. These results imply an initial Ru valence lying between 4+ and 5+, in agreement with X-ray Absorption Near Edge Spectroscopy (XANES) data which may be modelled as an admixture of Ru4+ and Ru5+ spectra in a 40:60 ratioLiu2000 .

As might be expected, doping of the Cu site has little effect on $`T_{Magnetic}`$ or the transport properties of the RuO2 layer. Li+ acts as a pair-breaking impurity in the CuO2 layer and causes a depression of $`T_c`$ in line with its behavior in other cuprates. Isovalent doping of the Gd site with other lanthanide elements changes the CuO2 layer doping level, with a remarkably strong variation in $`T_c`$. This appears to be a ion-size doping effect. Altervalent substitution of Ce for Gd rapidly reduces the doping level of the CuO2 layers and drives the material non-superconducting. In all but the Ce-doped samples, the conductivity of the RuO2 layer only ever reaches a modest fraction of that of the CuO2 layer. In the 10% Ce-doped sample the more heavily under-doped CuO2 layer has an insulating upturn at low temperature, whilst the RuO2 layer remains more metallic, and so the ratio of their conductivities reaches at least 1.9.

## Conclusions

To a first approximation the electronic properties of the CuO2 layer in RuSr2GdCu2O8 are the same as those of similar CuO2 layers in other high-$`T_c`$ cuprate superconductors in all respects. This conclusion is supported by the resistivity, TEP and Hall effect data presented here, and by results on the specific heat jump at $`T_c`$. On a more detailed level, magneto-transport measurements reveal an interaction between the carriers in these layers and the magnetization of the RuO2 layer. This interaction, with an energy which would seem to be of the same order as the SC energy gapMcCrone1999 , is not sufficient to destroy superconductivity.

The electronic properties of the RuO2 layer appear to bear a remarkable similarity to those observed in the ruthenate SrRuO3. At room temperature the conductivity of the layer is perhaps 10% of that of the CuO2 layer, with $`k_Fl_{Ru}0.2`$, indicating very badly metallic or localized behavior. Below $`T_{Magnetic}`$ the conductivity of the layer rises significantly - by at least 0.3$`\sigma _{Cu}`$. This increase raises the weighting of the RuO2 layer properties relative to those of the CuO2 layer in the admixture that determines the overall transport properties of RuSr2GdCu2O8. As the Hall effect and TEP of the RuO2 layer are both considerably smaller than those in the CuO2 layer the result is a drop in both $`R_H`$ and $`S`$ below $`T_{Magnetic}`$. In pure RuSr2GdCu2O8, and also in most of the substituted variants studied, $`\sigma _{Ru}`$ remains lower than $`\sigma _{Cu}`$ over the whole temperature range. For the Ce-doped samples studied, however, the CuO2 layer becomes insulating at low temperatures, allowing the poorly metallic RuO2 layer to dominate the conductivity, and its intrinsic transport properties to show strongly in the overall $`R_H`$ and $`S`$ of the material.

The two-band model of parallel conduction in the RuO2 and CuO2 layers has been very successful in modelling the transport properties observed in all the RuSr2GdCu2O8 samples studied, and it has been possible to describe well the effects of doping the different atomic sites. The inferred mixed valency of Ru, together with the onset of itinerancy at the magnetic transition suggests a possible role of a double-exchange mechanism in the magnetic interactions but also raises the possibility of charge ordering in these compounds at appropriate doping levels.

This work was supported by the U.K. Engineering and Physical Sciences Research Council and the New Zealand Marsden Fund for research (JLT) and travel funds (CB).


format: text

text:
# Discrete gap solitons in a diffraction-managed waveguide array

## I Introduction

### A Objectives of the work

Solitary-wave excitations in discrete nonlinear dynamical models (lattices) is a subject of great current interest, which was strongly bolstered by experimental observation of solitons in arrays of linearly coupled optical waveguides  and development of the diffraction management (DM) technique, which makes it possible to effectively control the discrete diffraction in the array, including a possibility to reverse its sign (make the diffraction anomalous) . It has recently been shown that a lattice subject to periodically modulated DM can also support stable solitons, both single-component ones  and two-component solitons with nonlinear coupling between the components via the cross-phase-modulation (XPM) .

Two-component nonlinear-wave systems, both continuum and discrete, which feature a linear coupling between the components, constitute a class of media which can support gap solitons (GSs). A commonly known example of a continuum medium that gives rise to GSs is a nonlinear optical fiber carrying a Bragg grating , whose standard model is based on the equations

$$\begin{array}{c}i\mathrm{\Psi }_t+i\mathrm{\Psi }_x+\left(|\mathrm{\Psi }|^2+2|\mathrm{\Phi }|^2\right)\mathrm{\Psi }+\mathrm{\Phi }=0,\hfill \\ i\mathrm{\Phi }_ti\mathrm{\Phi }_x+\left(|\mathrm{\Phi }|^2+2|\mathrm{\Psi }|^2\right)\mathrm{\Phi }+\mathrm{\Psi }=0,\hfill \end{array}$$
(1)
where $`\mathrm{\Psi }(x,t)`$ and $`\mathrm{\Phi }(x,t)`$ are amplitudes of the right- and left-propagating waves, and the Bragg-reflection coefficient is normalized to be $`1`$. Another optical system that may give rise to GSs is a dual-core optical fiber with asymmetric cores, in which the dispersion coefficients have opposite signs .

In this work, we demonstrate that the use of the DM technique provides for an opportunity to build a double lattice in which two discrete subsystems with opposite signs of the effective diffraction are linearly coupled, thus opening a way to theoretical and experimental study of discrete GSs, as well as of solitons of different types (solitons in linearly coupled lattices with identical discrete diffraction in the two subsystems have recently been considered in Ref. ; a possibility of the existence of discrete GSs in a model of a nonlinear-waveguide array consisting of alternating cores with two different values of the propagation constant was also considered recently ). The objective of the work is to introduce this class of systems and find fundamental solitons in them, including the investigation of their stability. We will also consider, in a more concise form, another physically relevant possibility, viz., a discrete system with a second-harmonic-generating (SHG) nonlinearity, in which the diffraction has opposite signs at the fundamental and second harmonics. Solitons will be found and investigated in the latter system too.

It is relevant to start with equations on which our $`\chi ^{(3)}`$ model (the one with the cubic nonlinearity) is based,

$`i{\displaystyle \frac{d\psi _n}{dt}}`$ $`=`$ $`\left(C\mathrm{\Delta }_2+q\right)\psi _n\left(\left|\psi _n\right|^2+\beta \left|\varphi _n\right|^2\right)\psi _n\kappa \varphi _n=0,`$ (2)
$`i{\displaystyle \frac{d\varphi _n}{dt}}`$ $`=`$ $`\delta \left(C\mathrm{\Delta }_2+q\right)\varphi _n\left(\left|\varphi _n\right|^2+\beta \left|\psi _n\right|^2\right)\varphi _n\kappa \psi _n=0,`$ (3)
where $`\psi _n(t)`$ and $`\varphi _n(t)`$ are complex dynamical variables in the two arrays (sublattices), $`\kappa `$ and $`\beta `$ being coefficients of the linear and XPM coupling between them, and $`t`$ is actually not time, but the propagation distance along waveguides, in the case of the most physically relevant optical interpretation of the model. The operators $`C\mathrm{\Delta }_2\psi _nC\left(\psi _{n+1}+\psi _{n1}2\psi _n\right)`$ and $`\left(\delta \right)C\mathrm{\Delta }_2\varphi _n\left(\delta \right)C\left(\varphi _{n+1}+\varphi _{n1}2\varphi _n\right)`$ represent discrete diffraction induced by the linear coupling between waveguides inside each array, the diffraction being normal in the first sublattice and anomalous in the second, with a negative relative diffraction coefficient $`\delta `$ and intersite coupling constant $`C`$ (one may always set $`C>0`$, which we assume below). Physical reasons for having $`\delta <0`$ are explained below. Finally, the real coefficient $`q`$ accounts for a wavenumber mismatch between the sublattices.

We also choose a similar SHG model, following the well-known pattern of discrete SHG systems with normal diffraction at both harmonics , :

$`i{\displaystyle \frac{d\psi _n}{dt}}`$ $`=`$ $`C\mathrm{\Delta }_2\psi _n\psi _n^{}\varphi _n,`$ (4)
$`2i{\displaystyle \frac{d\varphi _n}{dt}}`$ $`=`$ $`\delta C\mathrm{\Delta }_2\varphi _n\psi _n^2\kappa \varphi _n,`$ (5)
where the asterisk stands for the complex conjugation and $`\kappa `$ is a real mismatch parameter. In this case too, we assume $`\delta <0`$.

There are at least two different physical realizations of the $`\chi ^{(3)}`$ model based on Eqs. (2) and (3). First, one may consider two parallel arrays of nonlinear waveguides with different effective values $`n^{(1)}`$ and $`n^{(2)}`$ of the refractive index in them corresponding to a given (oblique) direction of the light propagation. To this end, the waveguides belonging to the two arrays may be fabricated from different materials; alternatively, they may simply differ by the transverse size of waveguiding cores, or by the refractive index of the filling between the cores, see e.g., Fig. 1. The difference in the effective refractive index gives rise to the mismatch parameter $`q`$ in Eqs. (2) and (3). More importantly, it may also give rise to different coefficients of the discrete diffraction. Indeed, the DM technique assumes launching light into the array obliquely, the effective diffraction coefficient in each array being

$$D^{(1,2)}=2Cd^2\mathrm{cos}\left(k_{}^{(1,2)}d\right),$$
(6)
where $`d`$ is the spacing of both arrays, and $`k_{}^{(1,2)}`$ are transverse components of the two optical wave vectors. As it follows from Eq. (6), the diffraction coefficients are different if $`k_{}^{(1)}k_{}^{(2)}`$.

Despite the fact that $`k_{}^{(1)}`$ and $`k_{}^{(2)}`$ are assumed different, we assume that the propagation directions of the light beams are parallel in the two arrays, as a conspicuous walkoff (misalignment) between them will easily destroy any coherent pattern. On the other hand, the light coupled into both arrays has the same frequency, hence the absolute values of the two wave vectors are related as follows: $`k^{(1)}/k^{(2)}=n^{(1)}/n^{(2)}`$, where $`n^{(1,2)}`$ are the above-mentioned effective refractive indices. Combining the latter relation and the classical refraction law, and taking into regard the condition that the propagation directions are parallel inside the arrays, one readily arrives at the conclusion that

$$k_{}^{(1)}/k_{}^{(2)}=n^{(1)}/n^{(2)}.$$
(7)
Note that the two incidence angles $`\theta ^{(1,2)}`$ (at the interface between the arrays and air) are related in a similar way, $`\left(\mathrm{sin}\theta ^{(1)}\right)/\left(\mathrm{sin}\theta ^{(2)}\right)=n^{(1)}/n^{(2)}`$, hence the incident beams (in air) must be misaligned, in order to be aligned in the arrays.

Equation (6) shows that there is a critical direction of the beam in each array, corresponding to $`k_{}^{(1,2)}d=\pi /2`$, at which the effective diffraction coefficient changes its sign . Due to the difference between $`k_{}^{(1)}`$ and $`k_{}^{(2)}`$, the critical directions are different in the two arrays. Then, if the common propagation direction in the arrays is chosen to be between the two critical directions, Eq. (6) gives different signs of the two diffractive coefficients. Note that this interpretation of the model implies no XPM coupling between the arrays, i.e., $`\beta =0`$ in Eqs. (2) and (3).

An alternative realization is possible in a single array of bimodal optical fibers, into which two parallel beams with orthogonal polarizations, $`u`$ and $`v`$, are launched obliquely. If the two polarizations are circular ones, then $`\beta =2`$ in Eqs. (2) and (3), and the asymmetry between the beams, which makes it possible to have different signs of the coefficient (6) for them, may be induced by birefringence, which, in turn, can be easily generated by twist applied to the fibers . The birefringence also gives rise to the mismatch $`q`$. As for the linear mixing between the two polarizations, which is assumed in the model, it can be easily induced if the fibers are, additionally, slightly deformed, having an elliptic cross section . If the two polarizations are linear, then the birefringence is induced by the elliptic deformation, and the linear mixing is induced by the twist, the XPM coefficient being $`2/3`$ in this case (assuming that, as usual, the birefringence makes it possible to neglect four-wave mixing nonlinear terms ).

It is interesting to note that the discrete model based on Eqs. (2) and (3) with $`\kappa =0`$ is exactly tantamount to a formal discretization of the above-mentioned continuum model which was introduced in Ref.  to describe a dual-core optical fiber with opposite signs of dispersion in the cores. Another quite noteworthy feature of the present model is that, if $`\beta =2`$, it turns out to be formally equivalent to a discretization of the standard Bragg-grating model (1), which is produced by replacing $`\mathrm{\Psi }_x\left(\mathrm{\Psi }_{n+1}\mathrm{\Psi }_{n1}\right)/2`$ and $`\mathrm{\Phi }_x\left(\mathrm{\Phi }_{n+1}\mathrm{\Phi }_{n1}\right)/2`$. Indeed, making the substitution (“staggering transformation”)

$$\mathrm{\Psi }_ni^n\varphi _n,\mathrm{\Phi }_ni^n\psi _n,$$
(8)
one concludes that the discrete version of Eqs. (1) takes precisely the form of Eqs. (2) and (3) with $`\delta =1,q=2C,\kappa =1`$, and $`\beta =2`$.

### B The linear spectrum

Before proceeding to the presentation of numerical results for solitons found in the system of Eqs. (2) and (3), it is relevant to understand at which values of the propagation constant $`\mathrm{\Lambda }`$ (spatial frequency) solitons with exponentially decaying tails may exist in this model. There are two regions in which they may be found. Firstly, inside the gap of the system’s linear spectrum one may find discrete gap solitons, i.e., counterparts of the GSs found in the continuum version of the model in Ref. . Secondly, solitons specific to the discrete model may be found above the phonon band. To analyze these possibilities, an asymptotic expression for the tail,

$$\psi _n,\varphi _n\mathrm{exp}\left(i\mathrm{\Lambda }t\lambda \left|n\right|\right)$$
(9)
is to be substituted into the linearized version of Eqs. (2) and (3).

Investigating the possibility of the existence of solitons above the phonon gap, it is sufficient to focus on the particular case $`\delta =1`$ and $`q=0`$, when the system’s spectrum takes a simple form (we have also considered more general cases with positive $`\delta `$ different from $`1`$ and $`q0`$, concluding that they do not yield anything essentially different from this case). The final result, produced by a straightforward algebra, is that solitons are possible in the region

$$\mathrm{\Lambda }^2>\mathrm{\Lambda }_{\mathrm{edge}}^216C^2+\kappa ^2\text{,}$$
(10)
$`\pm \mathrm{\Lambda }_{\mathrm{edge}}`$ being edges of the phonon band. In what follows below, we will assume $`\mathrm{\Lambda }>0`$, as in this case positive and negative values of $`\mathrm{\Lambda }`$ are equivalent.

To understand the possibility of the existence of the discrete GSs, we, first, set $`\delta =1`$ as above, but keep the mismatch $`q`$ as an arbitrary parameter. Then, the gap is easily found to be

$$\mathrm{\Lambda }^2<\mathrm{\Lambda }_{\mathrm{gap}}^2\{\begin{array}{cc}q^2+\kappa ^2& \mathrm{if}\text{q¡0},\\ \kappa ^2& \mathrm{if}\mathrm{\hspace{0.17em}\hspace{0.17em}0}q4C,\\ \left(q4C\right)^2+\kappa ^2& \mathrm{if}\text{q¿4C}\end{array}$$
(11)
(recall that, by definition, $`C>0`$). An essential role of the mismatch parameter is that it makes the gap broader if it is negative.

In the more general case, $`\delta 1`$, two different layers can be identified in the gap, similar to what was found in the continuum limit . For instance, if $`q=0`$, the inner and outer layers are

$$0<\mathrm{\Lambda }^2<\frac{4\delta }{\left(\delta +1\right)^2}\kappa ^2,\mathrm{and}\frac{4\delta }{\left(\delta +1\right)^2}\kappa ^2<\mathrm{\Lambda }^2<\kappa ^2$$
(12)
(in the case $`\delta =1`$, the outer layer disappears). The difference between the layers is the same as in the continuum limit : in the outer layer, solitons, if any, have monotonically decaying tails, i.e., real $`\lambda `$ in Eq. (9), while in the inner layer $`\lambda `$ is complex, and, accordingly, soliton tails are expected to decay with oscillations.

### C The structure of the work

The rest of the paper is organized as follows. In section II, we display results for solitons found above the phonon band, i.e., in the region (10). The evolution of the solitons is monitored, starting from the anti-continuum (AC) limit $`C=0`$, and gradually increasing $`C`$. Any branch of soliton solutions in this region must disappear, approaching the continuum limit. Indeed, as the radiation band (frequently called “phonon band”, referring to linear phonon modes in the lattice dynamics) becomes infinitely broad in this limit, see Eq. (10), the solution branch with $`\mathrm{\Lambda }=\mathrm{const}`$ will crash hitting the swelling phonon band. However, in many cases the soliton of this type is found to remain stable as long as it exists, so it may be easily observed experimentally in the optical array.

In section III we present results for solitons existing inside the gap. In the outer layer \[which is defined as per Eq. (12), provided that $`\delta 1`$\], we were able to find only solitons of an “antidark” type, that sat on top of a nonvanishing background. However, in the inner layer \[recall it occupies the entire gap in the case $`\delta =1`$, according to Eq. (12)\], true solitons are easily found (in accord with the prediction, their tails decay with oscillations). In the case $`q=0`$, these solutions appear as stable ones in the AC limit, get destabilized at some finite critical value of $`C`$, and continue, as unstable solutions, all the way up to the continuum limit, never disappearing. It is quite interesting that sufficiently large negative mismatch strongly extends the stability range for these solitons.

As was mentioned above, the $`\chi ^{(3)}`$ model based on Eqs. (2) and (3) may be considered as a discretization of the standard gap-soliton system (1). In this connection, it is natural to search for discrete counterparts of the usual GSs in the latter system. However, the discrete GSs found in section III do not have any counterpart in the continuum system (1), as the staggering transformation (8) makes direct transition from the discrete equations (2) and (3) to the continuum system (1) impossible. At the end of section III, we specially consider discrete solitons which are directly related to GSs in the system (1). We find that such solitons exist indeed at all the values of $`C`$, their drastic difference from those found in sections II and III is that they are essentially complex solutions to the stationary version of Eqs. (2) and (3). At all finite values of $`C`$, they are unstable, but the instability asymptotically vanishes in the AC and continuum limits, $`C0`$ and $`C\mathrm{}`$.

In section IV, we briefly consider the SHG model (4), (5). Solitons are found in this model too, and their stability is investigated. When the solitons are linearly unstable, the development of their instability is examined (in all the sections II, III, and IV) by means of direct numerical simulations, which show that the instability may initiate a transition to a localized breather, or to lattice turbulence, or, sometimes, complete decay of the soliton into lattice phonon waves.

## II Solitons above the phonon band

### A General considerations

Stationary solutions to Eqs. (2)-(3) are sought for the form

$$\psi _n=e^{i\mathrm{\Lambda }t}u_n,\varphi _n=e^{i\mathrm{\Lambda }t}v_n,$$
(13)
where $`\mathrm{\Lambda }`$ is the propagation constant defined above. In figures displayed below, the stationary solutions will be characterized by the norms of their two components,

$$P_u^2\underset{n=\mathrm{}}{\overset{+\mathrm{}}{}}u_n^2,P_v^2\underset{n=\mathrm{}}{\overset{+\mathrm{}}{}}v_n^2.$$
(14)
Once such solutions are numerically identified by means of a Newton-type numerical scheme, we then proceed to investigate their stability, assuming that the solution is perturbed as follows:

$`\psi _n`$ $`=`$ $`\left[u_n+ϵa_n\mathrm{exp}(i\omega t)+ϵb_n\mathrm{exp}(i\omega ^{}t)\right]\mathrm{exp}(i\mathrm{\Lambda }t),`$ (15)
$`\varphi _n`$ $`=`$ $`\left[v_n+ϵc_n\mathrm{exp}(i\omega t)+d_n\mathrm{exp}(i\omega ^{}t)\right]\mathrm{exp}(i\mathrm{\Lambda }t),`$ (16)
where $`ϵ`$ is an infinitesimal amplitude of the perturbation, and $`\omega `$ is the eigenvalue corresponding to the linear (in)stability mode. The set of the resulting linearized equations for the perturbations $`\{a,b^{},c,d^{};\omega \}`$ is subsequently solved as an eigenvalue problem. This is done by using standard numerical linear algebra subroutines built into mathematical software packages . If all the eigenvalues $`\omega `$ are purely real, the solution is marginally stable; on the contrary, the presence of a nonzero imaginary part of $`\omega `$ indicates that the soliton is unstable. When the solutions were unstable, their dynamical evolution was followed by means of fourth-order Runge-Kutta numerical integrators, to identify the development and outcome of the corresponding instabilities.

In what follows below, we describe different classes of soliton solutions, which are generated, in the AC limit, by expressions with different symmetries. Still another class of solitons, which carries over into the usual GSs in the continuum system (1), will be considered in the next section.

### B Solution families which are symmetric in the anti-continuum limit

As it was said above, in this section we set $`\delta =1`$ and $`q=0`$, since comparison with more general numerically found results has demonstrated that this case adequately represents the general situation, as concerns the existence and stability of solitons. Figure 2 shows a family of soliton solutions found for $`\kappa =0.1`$, $`\mathrm{\Lambda }=2`$ and $`\beta =0`$, as a function of the coupling constant $`C`$. In this case, the family starts, in the AC limit ($`C=0`$), with a solution that consists of a symmetric excitation localized at a single lattice site $`n_0`$, with

$$u_{n_0}=v_{n_0}=\pm \sqrt{\frac{\mathrm{\Lambda }\kappa }{1+\beta }},$$
(17)
and terminates at finite $`C`$. Figure 2 demonstrates that this branch is always unstable. The termination of the branch happens when it comes close to the phonon band, that swells with the increase of $`C`$. The branch terminates at $`C=0.464`$, when the upper edge of the band is at $`\mathrm{\Lambda }_{\mathrm{edge}}=\sqrt{\kappa ^2+16C^2}1.859`$, according to Eq. (10). This value is still smaller than the fixed value of the soliton’s propagation constant, $`\mathrm{\Lambda }=2`$, for which the soliton branch is displayed in Fig. 2. The branch, if it could be continued, would crash into the upper edge of the phonon band at $`C=0.499`$. The slightly premature termination of this soliton family is a consequence of the nonlinear character of the solutions, as the above prediction for the termination point was based on the linear approximation.

An example of the development of the instability of this solution, as found from direct simulations of the full equations (2) and (3), is given in Fig. 3 for $`C=0.1`$. It is clearly seen that the unstable soliton turns into a stable breather.

On the contrary, in the presence of XPM with the physically relevant value of $`\beta =2`$, a similar solution branch, found for the same values $`\kappa =0.1`$ and $`\mathrm{\Lambda }=2`$, is stable for all $`C`$, until it terminates at $`C=0.499`$. Note that, at this point, the upper edge (10) of the phonon band is $`\mathrm{\Lambda }_{\mathrm{edge}}=1.999`$, which is extremely close to $`\mathrm{\Lambda }=2`$, i.e., the termination of the solution family is indeed accounted for by its crash into the swelling phonon band. Details of this stable branch are shown in Fig. 4.

Direct simulations of this solution have corroborated its stability (details are not shown here). In fact, in all the cases when solitons are found to be stable in terms of the linearization eigenvalues (see other cases below), direct simulations fully confirm their dynamical stability.

### C Solution families which are anti-symmetric in the anti-continuum limit

Another branch of solutions is initiated, in the AC limit, by an anti-symmetric excitation localized at a single lattice site, cf. Eq. (18):

$$u_{n_0}=v_{n_0}=\pm \sqrt{\frac{\mathrm{\Lambda }+\kappa }{1+\beta }}.$$
(18)
The solution belonging to this branch is shown in Fig. 5 for the same values of parameters as in Fig. 2, i.e., $`\delta =1`$, $`\kappa =0.1`$, $`\mathrm{\Lambda }=2`$, and $`\beta =0`$. With the increase of $`C`$, this branch picks up an oscillatory instability at $`C0.257`$, and terminates at $`C0.407`$. Unlike the solutions displayed above, the termination of this branch occurs not through its crash into the phonon band, but via a saddle-node bifurcation. The latter bifurcation implies a collision with another branch of solutions. That additional branch (which is strongly unstable) was found but is not shown here.

In fact, the numerical algorithm is able to capture other solutions (unstable ones) past the point $`C0.407`$ at which the present solution terminates. The newly found solutions are shown in the bottom part of Fig. 5. However, the new family cannot be continued beyond $`C=0.467`$ \[cf. the termination point $`C=0.464`$ for the solutions initiated in the AC limit by the expression (17)\].

The development of the oscillatory instability of the solution shown in Fig. 5 was also studied in direct simulations. It leads to onset of a state where one component of the soliton is fully destroyed \[it cannot completely disappear, due to the presence of the linear couplings in Eqs. (2) and (3), but it is reduced to a level of small random noise\]. An example of this is given in Fig. 6 for $`C=0.3`$, $`\beta =0`$ and $`\kappa =0.1`$. The instability (with the initial growth rate $`0.07`$ in this case) develops after $`t60`$, destroying one component of the soliton in favor of further growth of the other one. In this case, a uniformly distributed noise perturbation of an amplitude $`10^4`$ was added to accelerate the onset of the instability, as the initial instability is very weak (which implies that the unstable soliton may be observed in experiment).

A counterpart of the solution from Fig. 5, but with $`\beta =2`$, rather than $`\beta =0`$, is shown in Fig. 7. This branch is always unstable (i.e., in the case of the solutions starting from the anti-symmetric expression in the AC limit, the XPM nonlinearity destabilizes the solitons, while in the case of the branch that was initiated by the symmetric expression in the AC limit, the same XPM nonlinearity was stabilizing). It terminates at $`C0.219`$, again through a saddle-node bifurcation. As in the previous case, a new family of solutions can be captured by the numerical algorithm past the termination point. The new family is found for $`0.22<C<0.498`$, and it is also shown in Fig. 7. Comparing the value $`\mathrm{\Lambda }_{\mathrm{edge}}=1.995`$ given by Eq. (10) in this case with the actual value $`\mathrm{\Lambda }=2`$ of the soliton’s propagation constant, we conclude that the termination of the latter branch is caused by its collision with the phonon band. Notice also that the latter branch becomes unstable only very close to its termination point, at $`C>0.494`$.

In the case of $`\beta =2`$, direct simulations show that the instability of the anti-symmetric branch gives rise to rearrangement of the solution into a very regular breather shown in Fig. 8 for $`C=0.1`$ and $`\kappa =0.1`$.

### D Solution families which are asymmetric in the anti-continuum limit

Additional branches of the solutions may start in the AC limit from asymmetric configurations, provided that $`\mathrm{\Lambda }`$ is still larger, namely for $`\mathrm{\Lambda }>2\kappa `$. In particular, such an extra branch can be initiated by the following AC-limit solution excited at a single site $`n=n_0`$ (here, $`\beta =0`$), cf. Eqs. (17) and (18):

$`u_{n_0}^2`$ $`=`$ $`{\displaystyle \frac{1}{2}}\left[\mathrm{\Lambda }\pm \sqrt{\mathrm{\Lambda }^24\kappa ^2}\right],`$ (19)
$`v_{n_0}`$ $`=`$ $`\kappa ^1(\mathrm{\Lambda }u_{n_0}u_{n_0}^3).`$ (20)
An example of this solution for the upper sign in Eq. (19) is shown, for $`\mathrm{\Lambda }=2`$, $`\kappa =0.5`$ and $`\delta =1`$, in Fig. 9. Such asymmetric branches may be stable for sufficiently weak coupling (in this case, for $`C<0.204`$), but they eventually become unstable, and disappear soon thereafter (at $`C0.213`$, in this case).

The evolution of the instability (for $`C>0.204`$) for this asymmetric branch is strongly reminiscent of that shown in Fig. 3, resulting in a persistent breathing state.

The branch that commences from the AC expression (19) with the lower sign is shown for $`\mathrm{\Lambda }=2`$, $`\kappa =0.75`$ and $`\delta =1`$ in Fig. 10. The branch remains stable as long as it exists, i.e., for $`C<0.46`$. At this point, it disappears colliding with the phonon band, whose upper edge is located, according to Eq. (10), at $`\mathrm{\Lambda }_{\mathrm{edge}}1.987`$, which is very close to the family’s fixed propagation constant, $`\mathrm{\Lambda }=2`$.

## III Gap solitons

### A Solitons in the inner layer of the gap

All the solutions that were examined in the previous section had their propagation constant above the upper edge of the phonon spectrum. Another issue of obvious interest is to study possible gap solitons (GSs), whose propagation constant is located inside the gap (11), i.e., below the lower edge of the phonon band. Unlike the solitons found above the band, GSs may persist up to the continuum limit.

An example of such a solution for $`\mathrm{\Lambda }=0.75`$, $`\kappa =1`$, $`\delta =0.9`$, and $`\beta =0`$ is shown in Fig. 11. In the AC limit, this branch starts with the expression (17). The branch is stable for small $`C`$, but then it becomes unstable due to oscillatory instabilities. The first two instabilities occur at $`C=0.242`$ and $`C=0.349`$, as is shown in Fig. 11. Past the onset of the instabilities, this branch continues to exist (as an unstable one) indefinitely with the increase of $`C`$, and carries over into an (unstable) GS in the continuum limit. At large values of $`C`$, the distinct phonon bands, which are clearly seen in the example of the eigenvalue spectrum shown for $`C=0.4`$ in Fig. 11, eventually collide and, due to their opposite Krein signs (see the definition and discussion of these in Ref. ), which gives rise to a whole set of oscillatory instabilities. The result is clearly seen in the example of the eigenvalue spectrum shown in the bottom panel of Fig. 11 for a large value of the coupling constant, $`C=4`$. The characteristic size of the instability growth rate (largest imaginary part of the eigenvalue) is nearly the same for $`C=0.4`$ and $`C=4`$, in the latter case it being $`0.09`$. Notice, however, that, as the continuum limit is approached, the instabilities may be suppressed, in a part or completely, by finite-size effects (for an example of such finite-size restabilization, see Ref. ).

The development of the oscillatory instability of GS belonging to the inner layer is displayed, for $`C>0.242`$, in Fig. 12 for $`C=0.4`$, $`\delta =0.9`$, $`\kappa =1`$ and $`\beta =0`$. In this particular case, there are two oscillatory instabilities whose growth rates are in the interval $`0.1<\omega _i<0.2`$. As a result, symmetry breaking occurs, resulting in a shift of the central position of the soliton (from the site $`n=50`$ to $`n=49`$). Oscillatory features in the dynamics are also observed in the latter case, and a small amount of energy is emitted as radiation.

Similar results were obtained for smaller values of $`\mathrm{\Lambda }`$, for instance, $`\mathrm{\Lambda }=0.25`$. It was verified too that this scenario persists in the presence of the XPM nonlinearity (i.e., for $`\beta =2`$), as it is shown in Fig. 13. In the latter case, the evolution of the instability with the increase of $`C`$ is quite interesting, as it is nonmonotonic. The instability first arises at $`C0.11`$ (due to a collision between discrete eigenvalues with opposite Krein signs). Subsequent restabilization takes place at $`C0.16`$, but the solutions are unstable again for $`C>0.36,`$ and remain unstable thereafter, up to the continuum limit.

In the case of $`\beta =2`$, the dynamical development of the oscillatory instabilities is similar to the $`\beta =0`$ case, again demonstrating symmetry-breaking effects.

### B Solutions in the outer layer of the gap

In all the cases considered in the previous subsections, the soliton’s propagation constant $`\mathrm{\Lambda }`$ belonged to the inner layer of the gap, see Eq. (12). We have also examined the situation when $`\mathrm{\Lambda }`$ belongs to the outer layer defined in Eq. (12) (the outer layer exists unless $`\delta =1`$). An example is shown in Fig. 14, where $`\beta =0`$, $`\kappa =1`$, $`\delta =0.1`$, and $`\mathrm{\Lambda }0.787`$ is chosen to be in the middle of the outer layer. In this case, we typically obtained delocalized solitons, sitting on top of a finite background (they are sometimes called “antidark” solitons). As can be observed from Fig. 14, such solutions may be stable for sufficiently weak coupling, but become unstable as the continuum limit is approached, although they do not disappear in this limit (in Ref.  such delocalized solitons were found in the continuum counterpart of the present model).

The instability development in the case of the outer-layer GSs is demonstrated, for $`C=0.549`$, $`\beta =0`$, $`\kappa =1`$ and $`\delta =0.1`$, in Fig. 15. In this case, the non-vanishing background is also perturbed by the instability, resulting in, plausibly, chaotic oscillations throughout the lattice. Symmetry-breaking effects, which shift the central peak from its original position, are observed too in this case.

### C Stabilization of the gap solitons by mismatch

The above considerations show that, inside the inner layer of the gap, it is easy to identify families of soliton solutions that persist in the continuum limit as $`C\mathrm{}`$. However, all the examples considered above showed that the solutions get destabilized at finite $`C`$ and remain unstable with the subsequent increase of $`C`$. Therefore, a challenging problem is to find solution families that would remain stable for large values of $`C`$.

In fact, the introduction of a finite mismatch $`q`$ (recall it was set equal to zero in all the examples considered above) may easily stabilize the discrete GSs. To this end, we pick up a typical example, with $`C=0.5`$, $`\kappa =1`$, $`\delta =0.5`$, $`\mathrm{\Lambda }=0.75`$, and $`\beta =0`$, when the GS exists but is definitely unstable in the absence of the mismatch. Figures 16 and 18 show the effect of positive and negative values of the mismatch on the solitons. As is seen, large values of the positive mismatch can make the instability very weak, but cannot completely eliminate it. However, sufficiently large negative mismatch readily makes the solitons truly stable. Thus, adding the negative mismatch is the simplest way to stabilize the solitons at large $`C`$, which is not surprising, as Eq. (11) demonstrates that the negative mismatch makes the gap broader.

As an example of the dynamical evolution of unstable solitons in the case of positive mismatch, in Fig. 17 we display the case of $`C=0.5`$, $`\kappa =1`$, $`\beta =0`$, $`\delta =0.5`$ and $`q=1`$. In this case, the evolution leads to the establishment of a breather with a rather complex dynamical behavior.

The instability development in the case of negative mismatch, $`q=2.5`$, is demonstrated in Fig. 19. A localized breather with quasi-periodic intrinsic dynamics is observed in this case as an eventual state.

### D Discrete counterparts of gap solitons from the Bragg-grating model

As was shown in the introduction, the particular case of Eqs. (2) and (3) with $`\kappa =1`$, $`\delta =1`$, and $`\beta =2`$ may be interpreted, with regard to the transformation (8), as a discretization of the standard Bragg-grating system (1). This continuum model gives rise to a family of exact GS solutions ,

$$\mathrm{\Psi }=U(x)\mathrm{exp}\left(it\mathrm{cos}\theta \right),\mathrm{\Phi }=V(x)\mathrm{exp}\left(it\mathrm{cos}\theta \right),$$
(21)
$$U(x)=\frac{\mathrm{sin}\theta }{\sqrt{3}}\mathrm{sech}\left(x\mathrm{sin}\theta \frac{i}{2}\theta \right),V=U^{},$$
(22)
where the real parameter $`\theta `$ takes values $`0<\theta <\pi `$. A part of this interval, $`0<\theta <\theta _{\mathrm{cr}}1.01\left(\pi /2\right)`$, is filled with stable solitons , while the remaining part contains unstable ones.

All the discrete GSs considered above are not counterparts of the continuum solitons given by Eqs. (21) and (22). Establishing a direct correspondence between the latter ones and discrete solitons of Eqs. (2) and (3) is complicated by two problems: the transformation (8) does not have a continuum limit, and real symmetric or anti-symmetric GSs with $`\left|\mathrm{\Lambda }\right|<\kappa `$ do not exist in the AC limit, as is seen from Eqs. (17) and (18), i.e., the usual starting point of the analysis is not available in this case.

We have considered the discrete analogs of the Bragg-grating GSs in the following way. First, we took a formal discrete counterpart of the waveforms (21) and (22), and used them to obtain exact solutions of the (formal) discretization of the Bragg-grating model of Eq. (1). Then the transformation (8) was applied to these solutions, and the thus obtained expressions were used as an initial guess for finding a numerically exact stationary solution of Eqs. (2)-(3). This procedure naturally generates new solitons, a crucial difference of which from all the types considered above is that they are truly complex solutions, see examples in Figs. 20 and 21. In this case, we have used $`\theta =\pi /4`$.

Then, the solution was numerically continued, decreasing $`C`$, back to the AC limit, in order to identify its AC “stem”. The result is shown in Fig. 21. Obviously, this AC state is very different from all those considered above (in particular, it is complex).

Finally, linear stability eigenvalues were calculated for this new branch of the discrete GSs. The result (see Fig. 22) is that this branch is unstable for all finite values of $`C`$, getting asymptotically stable in both limits $`C0`$ and $`C\mathrm{}`$ (large values of $`C`$ are not shown in Fig. 22); the stability regained in the latter limit complies with the above-mentioned finding that a subfamily of the continuum Bragg-grating gap solitons are dynamically stable. Notice that the natural norm of the continuum soliton differs from that of the discrete one, given by Eq. (14), by an additional multiplier $`C^{1/2}`$ (which is proportional to the effective lattice spacing). We have checked that the thus renormalized norm of the soliton converges as $`C\mathrm{}`$, although data for large $`C`$ is not displayed in Fig. 22.

## IV Solitons in the model with the quadratic nonlinearity

Stationary solutions of the SHG system (4) and (5) are looked for in an obvious form, cf. Eqs. (13):

$$\psi _n=e^{i\mathrm{\Lambda }t}u_n,\varphi _n=e^{2i\mathrm{\Lambda }t}v_n,$$
(23)
and in this case we only consider the (most characteristic) case $`\delta =1`$. The linearization of Eqs. (4) and (5) demonstrates that one may expect termination of a soliton-solution branch, due to its collision with the phonon band, at (or close to) the point

$$\mathrm{\Lambda }=\kappa /4+C,$$
(24)
and the gap between two phonon bands is

$$0<\mathrm{\Lambda }<\kappa /4$$
(25)
(it exists only if $`\kappa >0`$).

Stationary solutions were constructed, again, by means of continuation starting from the AC limit, where the excitation localized on a single site of the lattice assumes the form

$`v_{n_0}`$ $`=`$ $`\mathrm{\Lambda },`$ (26)
$`u_{n_0}`$ $`=`$ $`\pm \sqrt{v_{n_0}(4\mathrm{\Lambda }\kappa )}.`$ (27)
Note that solutions with the propagation constant belonging to the gap (25) do not exist close to the AC limit. Indeed, the AC expression (27) shows that a necessary condition for its existence is $`4\mathrm{\Lambda }>\kappa `$. On the other hand, $`\mathrm{\Lambda }`$ stays in the gap (25) if $`4\mathrm{\Lambda }<\kappa `$, so the two conditions are incompatible.

A typical example of a numerically found soliton branch is shown, for the value of the mismatch parameter $`\kappa =0.9`$, in Fig. 23. This solution family is fixed by choosing $`\mathrm{\Lambda }=0.25`$, and starting from the expressions (26) and (27) with the upper sign. It is seen that the branch is stable for $`C<0.015`$, but then it becomes unstable, and eventually terminates at $`C0.025`$, in very good agreement with the prediction of Eq. (24).

Development of the instability of this SHG soliton branch for $`C>0.015`$ was numerically examined through direct simulations, results of which are presented in Fig. 24, for the case of $`C=0.02`$, the corresponding instability growth rate being $`0.03`$. A random uniformly distributed perturbation of an amplitude $`10^4`$ was added to the initial condition to accelerate the onset of the instability. The eventual result of the instability is the appearance of a breather with very regular periodic intrinsic vibrations.

A similar result for the case without a gap in the phonon spectrum \[see Eq. (25)\] is shown in Fig. 25 for $`\kappa =0.1`$, the other parameters being the same as in the previous case. This time, the branch terminates at $`C0.275`$, once again in complete agreement with the prediction of Eq. (24).

Lastly, another characteristic branch of solutions can be constructed starting from the pattern given by Eqs. (27) with the lower sign. This solution family is displayed in Fig. 26, for $`\kappa =0.75`$, $`\delta =1`$, and $`\mathrm{\Lambda }=0.25`$. The branch is stable for sufficiently weak coupling, but then it becomes unstable for $`C>0.046`$. The branch disappears colliding with the phonon band at $`C0.062`$, once again in full agreement with the prediction of Eq. (24).

An example of the development of instability of the present solution, that takes place at $`C>0.046`$, is shown in Fig. 27 for $`C=0.055`$ ($`\kappa =0.75`$; $`\delta =1`$). A random initial perturbation with an amplitude $`10^4`$ was added to the initial condition in this case. As is seen, the evolution results in complete destruction of the pulse into small-amplitude radiation waves.

## V Conclusion

In this work, we have introduced a model which includes two nonlinear dynamical chains with linear and nonlinear couplings between them, and opposite signs of the discrete diffraction inside the chains. In the case of the cubic nonlinearity, the model finds two distinct interpretations in terms of nonlinear optical waveguide arrays, based on the diffraction-management concept. A continuum limit of the model is tantamount to a dual-core nonlinear optical fiber with opposite signs of dispersion in the two cores. Simultaneously, the system is equivalent to a formal discretization of the standard model of Bragg-grating solitons. A straightforward discrete second-harmonic-generation \[$`\chi ^{(2)}`$\] model, with opposite signs of the diffractions at the fundamental and second harmonics, was introduced too. Starting from the anti-continuum (AC) limit and gradually increasing the coupling constant, soliton solutions in the $`\chi ^{(3)}`$ model were found, both above the phonon band and inside the gap. Above the gap, the solitons may be stable as long as they exist, but with transition to the continuum limit they inevitably disappear. On the contrary, solitons in the gap persist all the way up to the continuum limit. In the zero-mismatch case, they always become unstable before reaching the continuum limit, but finite mismatch may strongly stabilize them. A separate procedure had to be developed to search for discrete counterparts of the well-known Bragg-grating gap solitons. As a result, it was found that discrete solitons of this type exist at all values of the coupling constant $`C`$, but they appear to be stable solely in the limit cases $`C=0`$ and $`C=\mathrm{}`$. Solitons were also found in the $`\chi ^{(2)}`$ model. They too start as stable solutions, but then lose their stability.

In the cases when the solitons were found to be unstable, simulations of their dynamical evolution reveal a variety of different scenarios. These include establishment of localized breathers featuring periodic, quasi-periodic, or very complex intrinsic dynamics, or destruction of one component of the soliton, as well as symmetry-breaking effects, and even complete decay of both components into small-amplitude radiation. The outcome depends on the type of the nonlinearity (cubic or quadratic), and on the nature of the unstable solution.

## Acknowledgements

B.A.M. acknowledges hospitality of the Department of Applied Mathematics at the University of Colorado, Boulder, and of the Center for Nonlinear Studies at the Los Alamos National Laboratory. P.G.K. gratefully acknowledges the hospitality of the Center for Nonlinear Studies of the Los Alamos National Laboratory, as well as partial support from the University of Massachusetts through a Faculty Research Grant, from the Clay Mathematics Institute through a Special Project Prize Fellowship and from NSF-DMS-0204585. Work at Los Alamos is supported by the U.S. Department of Energy, under contract W-7405-ENG-36.

We appreciate valuable discussions with M.J. Ablowitz and A. Aceves.