text:
# Heavy Residues with A<90 from the Asymmetric Reaction of 20 AMeV 124Sn+27Al as a Sensitive Probe of the Onset of Multifragmentation

## 1 Introduction

The yields of the heavy residues, the large remnants of the heavy member of an asymmetric reacting pair of nuclei, are known to comprise a large fraction of the reaction cross section for intermediate energy nuclear collisions. The studies carried out in inverse kinematics at projectile energies of several tens of AMeV by Bazin et al., Faure-Ramstein et al. , Pfaff et al. , Hanold et al.  and Souliotis et al.  have shown the utility of the fragment separator approach in studying the heavy reaction products from Kr+X, Xe+X and Au+X mass asymmetric collisions at intermediate energies. Compared to asymmetric reactions in normal kinematics , where the low energies of residues ($``$0.015 AMeV) cause the loss of substantial portions ( 50%) of the product distributions, due to experimental thresholds , the studies in inverse kinematics using a high resolution spectrometer/detector system allow the observation of a wide range of nuclei including the projectile-like fragments (PLFs) and the heavy remnants of the hot nuclei originating from violent collisions. From such inclusive measurement, one can gain important information about the reaction mechanism, complementary to that obtained in exclusive studies, where only the light reaction partners are observed with high isotopic resolution.

For asymmetric reactions of a heavy beam with a light target nucleus ( typically C, Al ) the momentum transfer data obtained ( see e.g.  ) demonstrate the presence of reaction mechanisms ranging from quasi-elastic peripheral collisions to processes analogous to the complete fusion and incomplete fusion observed at low energies. However, since asymmetric reactions at several tens of AMeV are complex processes, the description of the data by the simple low energy concepts such as massive-transfer hypothesis or the high-energy geometric abrasion model have shown to be, despite qualitative agreement, far from quantitative.

With increasing target mass, the reactions become more symmetric and the observed yields of heavy residues are dominated by peripheral processes such as deep-inelastic transfer as, for instance, shown in detail in our recent work  in the reaction of 25 AMeV 86Kr + 64Ni. As it follows from the model analysis ( for details see the original work ), the projectile-like fragments observed in this reaction are heated by intense exchange of nucleons with the target. The excitation energy is sufficient for emission of complex fragments. Previous studies  of the collisions of heavy nuclei at intermediate energies have also shown that such reactions are predominantly binary. Evidence has been presented  for a sequential decay of one of the initial binary fragments leading to a three (or more) body final state. In some of the binary encounters the projectile-like fragments have been found  to have very high temperatures (T$``$7 MeV).

The focus of the present work is the formation of heavy residues in the reaction of 20 AMeV 124Sn on 27Al. The main motivation was the observation of heavy remnants of hot nuclei created in violent collisions at low impact parameters where a single hot source is created. In particular, it is of importance to establish in detail to what extent the observed properties of such heavy residues are influenced by the dynamics of the entrance channel ( pre-equilibrium emission, incomplete fusion ) and/or by the process of de-excitation ( emission of complex fragments, multifragmentation ). In the reaction 124Sn + 27Al one can assume that the residues with masses much lower than that of the beam originate from violent collisions. Such an assumption is supported by the work presented in , where it was shown that for the damped peripheral collisions 28Si + 112,124Sn at energies 30 – 50 AMeV the heavy target-like fragment remains rather cold, while the light projectile-like fragment is hot enough to undergo multifragmentation. This is caused by approximately equal sharing of the excitation energy imparted to the reaction partners due to nucleon exchange. In order to heat the heavy fragment to higher excitation energies, part of the lighter partner should fuse with the heavier one, thus converting the kinetic energy of the relative motion into heat. The recoil spectrometer MARS at the Cyclotron Institute of Texas A&M University offers the possibility to carry out this study with an appropriate angular and momentum acceptance, under high-resolution conditions and in the appropriate region of nuclei.

## 2 Experimental Method and Data Analysis

The present study was performed at the Cyclotron Institute of Texas A&M University. A 20 AMeV 124Sn27+ beam from the K500 superconducting cyclotron, with a typical current of $``$0.5 pnA, interacted with a 27Al target of thickness 2.0 mg/cm2. The reaction products were analyzed with the MARS recoil separator . The primary beam struck the target at 0o relative to the optical axis of the spectrometer. The direct beam was collected in a small square Faraday cup approx. 30 cm after the target, blocking the angular range 0.0–1.0o. The fragments were accepted in the remaining angular opening of MARS: 1.0–2.7o (the angular acceptance of MARS is 9 msr ). This angular range encompasses the grazing angle of 2.1o  for the present reaction. MARS optics  provides one intermediate dispersive image and a final achromatic image (focal plane) and offers a momentum acceptance of 4%.

At the focal plane, the fragments were collected in a large area (5$`\times `$5 cm) three-element ($`\mathrm{\Delta }`$E1, $`\mathrm{\Delta }`$E2, E) Si detector telescope. The $`\mathrm{\Delta }`$E1 detector was a position-sensitive Si strip detector of 63 $`\mu `$m thickness, whereas the $`\mathrm{\Delta }`$E2 and the E detector were single-element Si detectors of 150 and 950 $`\mu `$m, respectively. The position information from the $`\mathrm{\Delta }`$E1 strips provided a continuous monitoring of the focusing and collection of the fragments at the various settings of the separator. Time of flight was measured between two parallel plate avalanche counters (PPACs)  positioned at the dispersive image and at the focal plane, respectively, and separated by a distance of 13.2 m. The PPAC at the dispersive image was also X–Y position sensitive and used to record the position of the reaction products. The horizontal position, along with NMR measurements of the field of the MARS first dipole, was used to determine the magnetic rigidity $`B\rho `$ of the particles. Thus, the reaction products were characterized by an event-by-event measurement of the energy loss, residual energy, time of flight, and magnetic rigidity. The response of the spectrometer and detector system to ions of known atomic number Z, mass number A, ionic charge q and velocity was calibrated using low intensity primary beams of 124Sn at 20 AMeV and 40Ar, 44Ca and 86Kr at 25 AMeV. To cover the N/Z and velocity range of the fragments, a series of measurements was performed at overlapping magnetic rigidity settings in the range 1.3–1.6 Tesla-meters.

The determination of the atomic number Z was based on the energy loss of the particles in the first $`\mathrm{\Delta }E`$ detector  and their velocity, with a resulting resolution (FWHM) of 1.0 Z units for near-projectile fragments and 0.6 Z units for A$`<`$90. The ionic charge $`q`$ of the particles entering the spectrometer after the Al stripper, was obtained from the total energy Etot=$`\mathrm{\Delta }`$E1+$`\mathrm{\Delta }`$E2+E, the velocity and the magnetic rigidity according to the expression:

$$q=\frac{3.107}{931.5}\frac{E_{tot}}{B\rho (\gamma 1)}\beta \gamma $$
(1)
where Etot is in MeV, B$`\rho `$ in Tm, $`\beta =\upsilon /c`$ and $`\gamma =1/(1\beta ^2)^{\frac{1}{2}}`$. The measurement of the ionic charge q had a resolution of 0.8 Q units (FWHM) for near-projectile fragments and 0.5 Q units for A$`<`$90. Since the ionic charge must be an integer, we assigned integer values of q for each event by putting appropriate windows on each peak of the q spectrum at each magnetic rigidity setting of the spectrometer. Using the magnetic rigidity and velocity measurement, the mass-to-charge A/q ratio of each ion was obtained from the expression:

$$A/q=\frac{B\rho }{3.107\beta \gamma }$$
(2)
Combining the q determination with the A/q measurement, the mass A was obtained as:

$$A=q_{int}\times A/q$$
(3)
(qint is the integer ionic charge determined as above) with an overall resolution (FWHM) of 1.0 A units for near-projectile fragments and about 0.6 A units for A $`<`$ 90. We refer to our previous work in ref. , carried out using the same experimental setup within the same run, for more details. The reconstruction of the Z, q and A and the gating procedure were applied to the calibration beam data to ensure the reproduction of the expected Z, q and A values and the elimination of spurious yield contributions from neighboring Z and q values. For the heavier products with masses A $`>`$ 90 the experimental resolutions did not allow the mass to be resolved unambiguously. However, the gross features for such products can still be obtained. Close to the beam, part of the yield could not be detected due to background from the quasi-elastically scattered beam.

Combination and appropriate normalization of the data at various magnetic rigidity settings of the spectrometer provided fragment distributions with respect to Z, A, q and velocity. Correction of missing yields caused by charge changing at the PPAC (positioned at the dispersive image) was performed based on the equilibrium charge state prescriptions of Leon et. al. . The overall data reduction procedure was similar to that followed in earlier work on 197Au fragmentation  and 238U projectile fission  at 20 AMeV. The isotope distributions were subsequently summed over all values of q. It should be pointed out that the resulting distributions in Z, A and velocity are the fragment yield distributions measured in the angle interval 1.0–2.7o and in the magnetic rigidity range 1.3–1.6 T m.

## 3 Results and Discussion

The gross features of the measured data from the reaction of 20 AMeV 124Sn + 27Al are displayed in Figs. 1,2. In Fig. 1 isotopic yields ( contour plot ) are presented as a function of mass and atomic number. The masses of observed products cover the range from the projectile-like fragments to the border of the intermediate mass fragment domain at Z=20. The yields of projectile-like nuclei peak around the nuclide 116Te. The position of this peak is caused by missing yield from the region close to the beam where the reaction products cannot be resolved using the existing detector set-up due to significant background from quasi-elastic processes. At masses below A=90 the distribution becomes rather flat and follows approximately the corridor of stable isotopes ( thick dashed line ). These are the reaction products which are expected to originate from the hot excited systems. Fig. 2 presents the velocity distributions of residues with selected atomic numbers. In the region close to beam the velocity distributions appear to be convolutions of multiple contributions from different reaction channels ranging from quasielastic products close to the beam velocity to various incomplete fusion channels where the velocity decreases with the amount of mass transferred from the target. With decreasing mass of the reaction product, the velocity distributions develop toward a single contribution of a Gaussian shape with a low velocity tail. The centroids of these distributions shift with decreasing mass toward higher velocity. Such a trend is caused by the kinematic selection of the spectrometer.

The experimental heavy residue data presented in Figs. 1,2 contain useful information on the production mechanism. The properties of the final products are the result of a complex process including an early stage dominated by the dynamics of the entrance channel which is governed by the impact parameter and a de-excitation stage where the final partitions of the reaction products are generated. In such a case, a viable method of reaction mechanism analysis appears to be the use of various model frameworks for both the initial stage of the collision and the de-excitation. The model which eventually proves superior to the others can be considered as reflecting the physical process in most detail.

As a first choice for the description of the initial stage we use the model framework described in . The basic features of the reaction mechanism model for violent collisions are the pre-equilibrium emission and the incomplete fusion (ICF). The pre-equilibrium emission is treated using a simplified variant of the exciton model employing a phenomenological parametrization of emission probability as a function of exciton number and angular momentum. The incomplete fusion model is based on the concept of geometrical fragmentation refined for the Fermi energy domain where incomplete fusion occurs by fusion of the participant zone with one of the spectators. It is applied to the reconstructed projectile-like and target-like prefragments formed in the pre-equilibrium stage. For the dissipative peripheral collisions the model of deep inelastic transfer (DIT) is used as implemented by Tassan-Got and Stefan . The DIT scenario is employed for events where the overlap of projectile and target nuclei does not exceed 3 fm. Again the DIT stage is preceded by the pre-equilibrium stage. Such a hybrid framework proved rather successful  in the description of a wide range of data obtained in experiments ranging from inclusive measurements to highly exclusive ones in 4$`\pi `$ geometry and can be considered appropriate for the reaction investigated in the present work, with a capability of quantitative description. Throughout the paper, this calculation will be named for simplicity as DIT+ICF. Nevertheless, one has to keep in mind that pre-equilibrium emission is an essential part of this model.

In order to describe the de-excitation stage of hot nuclei, there exist several concepts implemented in various codes. The statistical de-excitation code GEMINI  uses Monte Carlo techniques and the Hauser-Feshbach formalism to calculate the probabilities for fragment emission with Z$``$2. Heavier fragment emission probabilities are calculated using the transition state formalism of Moretto . Within such a model, the final partition of products is generated by a succession of fragment emissions ( binary decays ). Alternative to this scenario is the model of statistical multifragmentation, where the fragment partition is generated at once in the so-called freeze-out configuration. In this work, we use the code SMM  as a representative implementation of the concept of prompt multifragmentation. Thermally equilibrated partitions of hot fragments are generated in the hot stage, which is followed by propagation of fragments in their mutual Coulomb field and secondary de-excitation of hot fragments flying in their asymptotic directions.

The results of the DIT+ICF/GEMINI calculation, compared to experimental observables are given in Fig. 3. In the GEMINI calculations, we used essentially the parameter set recommended by the author, featuring Lestone’s temperature dependent level density parameter , a fading of shell corrections with excitation energy and enabled IMF emission. This parameter set proved successful in our recent work  on the reaction 86Kr+64Ni at 25 AMeV.

In Fig. 3a the mass yield curve is presented. The measured data, normalized for beam current and target thickness are given in mb and presented as solid circles. The result of the DIT+ICF/GEMINI calculation, filtered by the spectrometer angular and momentum acceptance is given by the full line, whereas the dashed line gives the total (unfiltered) yield. A comparison of the measured yields to the calculated filtered yields shows reasonable agreement for the heavier projectile-like fragments (A$`>`$90). The missing experimental cross section at the masses close to the beam is caused by the limitations of the experimental set-up as explained in the discussion of Fig. 1. The yields of residues with A$`<`$90 are increasingly underestimated by the calculation despite the fact that the unfiltered calculated yields are rather flat in this region, which appears to reflect the trend of the experimental data. Since, according to Fig. 3a, the yield of unfiltered residues with A$`<`$90 is rather flat, the missing filtered yield below A=90 in the calculation appears to be caused by the fact that the kinematic properties of the simulated residues with A$`<`$90 increasingly miss the spectrometer acceptance with decreasing mass, leading to increasing losses in the filtering procedure.

In Fig. 3b, the calculated and measured yield distributions as a function of Z (relative to the line of $`\beta `$-stability, Zβ) and A are presented. The line of stability is calculated as: Zβ = A/(1.98 + 0.0155A2/3) . The calculated average values from DIT+ICF/GEMINI are shown as a thick dashed line (without acceptance cut) and as a thick full line (with spectrometer acceptance cut). The measured yield distribution is represented by a contour plot. The DIT+ICF/GEMINI calculation describes reasonably well the centroids of the experimental data.

Finally, in Fig. 3c, the velocity vs. mass distributions are given. The data are again shown as contours. The thick dashed line is from the DIT+ICF/GEMINI calculation without acceptance cut and the full line is with acceptance cut. In this case the filtered calculated data appear to follow the experimental trend at masses above A=90, while they exhibit slight increase of the mean velocity with decreasing masses below A=90.

An alternative calculation was carried out, where the SMM code was used for the de-excitation of the hot nuclei produced by the DIT+ICF simulation. In the SMM calculation, a freeze-out configuration with hot primary fragments was assumed. Hot fragments are propagated in the Coulomb field and de-excited by secondary emission. Only thermal excitation energy is used as input while the rotational energy ( typically not exceeding 10 MeV ) is subtracted. The results are presented in Fig. 4 in a fashion analogous to Fig. 3. As one can see in all panels of Fig.4, the DIT+ICF/SMM calculation provides very consistent description of experimental observables for A$`<`$105 ( the discrepancies in the region close to the beam are analogous to previous cases ). Using the ratio of filtered to unfiltered calculated yield for each mass, correction factors (whose magnitude are inferred from Fig. 4a) for the acceptance of the spectrometer can be obtained as a function of mass. When applied to the measured yield data, an estimate of the total yield, shown by the open circles in Fig. 4a, could be obtained.

Thus, implementation of the prompt multifragmentation scenario appears to lead to production of heavy residues with proper kinematics. Compared to light particles or intermediate mass fragments (IMFs) used for imaging of the emitting source via particle-particle correlations, the experimentally detected heavy residues possess direct information on the properties of hot multifragment partitions. The process of secondary emission, as can be concluded from the simulations, does not influence significantly the kinematic properties of the heavy residues with masses A=40–90, since emission of nucleons is a dominant channel of secondary de-excitation.

Detailed insight into the different kinematic properties of residues simulated using GEMINI and SMM codes can be obtained from Figs. 5,6 and 7. In Fig. 5 the calculated angular distributions as a function of residue mass for GEMINI (a) and SMM (b) are presented. Two horizontal lines mark the angular acceptance of the MARS separator. It is remarkable to notice that, in fact, the gross features of both distributions are very similar and the experimental effect appears to be caused by the distant tail of the distribution, which extends much further toward zero angle in the scenario where the hot nucleus disintegrates at once into more pieces. Nevertheless, the essential feature of the zero angle region is that it is highly selective toward products from incomplete or complete fusion channel since the hot source is flying essentially along the beam direction and the final angular distribution is determined by the de-excitation process. In the case when the hot nucleus disintegrates by sequential binary decay, the recoil from emitted fragments causes a shift of the residue angle away from the zero angle. In Figs. 6,7 the calculated velocity distributions as a function of mass obtained using the GEMINI and SMM codes are given, respectively. In both cases, unfiltered yields are presented in panel (a), while the filtered yields are given in panel (b). As one can observe in panel (a), the SMM calculation leads to a much wider velocity distribution for A$`<`$90 which corresponds to a flatter angular distribution extending further toward zero angle. In the GEMINI calculation, the velocity distribution is concentrated close to the center of the distribution for a given mass, while the tails toward higher/lower velocities are suppressed since such residues, due to kinematics, can be produced only close to zero angle. On the other hand, when comparing the fastest residues for a given mass at A$`<`$90, which are emitted close to zero angle, the GEMINI calculation leads to higher maximum velocities, due to a larger recoil caused by the scenario of sequential binary decay. Nevertheless, due to the momentum coverage of the whole measurement which is close to 20 %, such a difference does not have a strong influence on the filtered yield. The effect of angular acceptance of MARS is dominant and determines the filtered yields which are much higher in the SMM simulation than in the GEMINI simulation.

In Fig. 8 the normalized unfiltered residue mass distributions at various impact parameters are given for both calculations. For DIT+ICF/SMM, one can notice a rather strong dependence of the production rate for residues with A$`<`$90 on impact parameter. For the DIT+ICF/GEMINI calculation, the dependence is rather weak. Significant contribution of A$`<`$90 residues is produced already at b=8 fm and toward lower impact parameters the shape of the distribution changes only moderately. As one can see in Fig. 9a, the average excitation energy of the hot source obtained using DIT+ICF model changes with impact parameter rather strongly. At b=8 fm, the average excitation energy amounts only to 80 MeV, which translates into relative excitation below 1 AMeV. At such excitation energy, the emission of IMFs is rather improbable and the simulated residues with A$`<`$90 can be expected to be fission fragments in the traditional sense. However, the mass distribution of such residues is virtually flat due to a weak asymmetry dependence of fission barriers of nuclei with fissility $`x`$ = 0.4 - 0.5, as suggested by theory . In the GEMINI calculation, such fission fragments are produced with a probability of about two orders of magnitude larger than in the SMM calculation. This rather large difference in yield reflects a different treatment of fission in each code. The code GEMINI uses the asymmetry-dependent fission barriers  normalized to angular momentum dependent finite range fission barriers of Sierk . On the other hand, the SMM code, where fission and multifragmentation are treated separately, uses the fission barriers of Barashenkov  without explicit angular momentum dependence ( the mass distributions in the SMM are generated according to the parametrization of Adeev  ). For example, the nucleus 124Sn at angular momentum $`J`$=0 is supposed to have the fission barrier of 41.5 MeV in GEMINI and 44.8 MeV in SMM. At $`J`$=25, which is a typical value of angular momentum at b=8 fm ( see Fig. 9b ), the fission barrier of Sierk drops further to 36.9 MeV. Assuming a level density parameter $`a`$=A/10, the difference of level density indeed amounts to slightly less than two orders of magnitude, consistent with Fig. 8. Thus the production rates of residues with A$`<`$90 at b=8 fm are consistent with the production rates of fission fragments. An analogous conclusion applies to impact parameter b=6 fm. At central impact parameters b=4 fm and b=2 fm, where the average relative excitation energy approaches and exceeds 2 AMeV, the production rate of residues with A$`<`$90 increases dramatically in the SMM calculation and exceeds the production rate of GEMINI calculation, since the multifragmentation threshold is exceeded by an increasing part of the excitation energy distributions at central impact parameters b=4 fm and b=2 fm, as one can deduce from Fig. 10.

For the heavier asymmetric system 136Xe+48Ti at 18.5 AMeV, where, due to higher fissility, the fission mass distribution is dominated by symmetric mass splits, the experimental investigation of 3-body events carried out by Gui et al.  revealed that the sequential fission mode with symmetric mass distribution dominates at lower excitation energies, while at higher excitation energies the de-excitation mode with wider asymmetry range and short emission time was observed. The asymmetric de-excitation mode with short timescale can, in principle, be related to the onset of multifragmentation. The decrease of the emission time with increasing mass asymmetry was observed also for 3-body events in the symmetric reactions 100Mo+100Mo and 120Sn+120Sn around 20 AMeV , where, on the other hand, the typical spatial configuration with the small fragment located between the two large fragments suggests the occurence of dynamical effects specific to emission from the symmetric binary configuration.

At excitation energies above the multifragmentation threshold, the SMM scenario leads to competition of the hot multifragment partitions with the single hot thermally equilibrated fragment ( practically identical to the compound nucleus in traditional sense ). Within the sequential binary decay scenario with finite range fragment emission barriers, as implemented in GEMINI, the onset of multifragment channels is much smoother. Apart from the reduction of the fission barrier, the angular momentum of the hot residue leads to a reduction of the thermal excitation energy due to the energy of rotational motion. However, the values of angular momentum in Fig. 9b imply that in the present work the rotational energy typically does not exceed 5 MeV and thus does not influence the results, especially at mid-central to central collisions which lead to production of heavy residues with A$`<`$90. An interesting aspect of the excitation energy versus mass distributions, shown in Fig. 10, is the gradual change of the governing pattern from the essentially Wilczynski plots at peripheral impact parameters to the pattern typical for incomplete fusion scenario at central impact parameters. Thus, the success of SMM implies that the residues with A$`<`$90 observed at forward angles originate from incomplete fusion.

For further comparison, the calculations using only the model of deep inelastic transfer combined with both GEMINI and SMM, analogous to DIT/GEMINI simulation used in our recent work , have been carried out. In the case of the nearly symmetric reaction 86Kr+64Ni at 25 AMeV , the DIT/GEMINI simulations proved to be rather successful. Similar success was achieved in the reaction of 20 AMeV 124Sn + 124Sn , especially for the residues with A$`<`$90 on which this work focuses. It was verified using both DIT/GEMINI and DIT/SMM that, in the reaction 124Sn + 27Al at 20 AMeV, the yields of residues with A$`<`$90 can not be reproduced by taking into account peripheral processes only, independent of which de-excitation code was used. For instance, in the DIT/GEMINI simulations, the total unfiltered yields are roughly analogous to the simulation presented in Fig. 3 and the filtered yields again do not reproduce the observed yields of residues with A$`<`$90. In this case, the filtered distribution exhibits similar decrease in the A$`<`$90 region but extends further than in DIT+ICF/GEMINI case, due to higher excitation energy and different angular and momentum distribution of the generated hot nuclei. The difference can be attributed mostly to the absence of cooling via pre-equilibrium emission, since the products of violent collisions are not included in the DIT/GEMINI simulation and do not pass the filter procedure in the latter case. On the other hand, the DIT/SMM simulation, when compared to DIT+ICF/SMM simulation shown in Fig. 4, significantly underpredicts the yields of the residues with A$`<`$90 ( both unfiltered and filtered ), due to the fact that for peripheral collisions, the excitation energy of projectile-like fragments exceeds the value of 2 AMeV, where the multifragmentation threshold can be anticipated, only very rarely, and due to orbiting, such sources move typically at angles away from zero.

Furthermore, a backtracing procedure was used for the DIT+ICF/SMM simulation. In Fig. 11a,b,c the mass, charge, excitation energy distributions of hot heavy sources contributing to filtered yield of residues with A$`=`$65-75 are shown. Fig. 11d shows the distribution of contributing impact parameters. On average, a typical contributing hot source can be characterized as a nucleus 144Nd with excitation energy about 310 MeV ( 2.2 AMeV ). The average impact parameter is about 2 fm ( $`l`$ = 44 ) and the average spin of the hot source is $`J`$=30. Thus the results of the backtracing procedure are in good agreement with the results of the analysis presented above. The observed ( symbols ) and filtered ( solid histograms ) velocity and N/Z-distributions of residues with A$`=`$65-75 along with the velocity and N/Z-distributions of contributing hot heavy sources ( dashed histograms ) are compared in Fig. 12. The simulated velocity and N/Z-distributions reproduce in high detail the experimental ones, thus further demonstrating the appropriate description of reaction mechanism within the DIT+ICF/SMM simulation and the appropriate treatment of the spectrometer in the filtering procedure. The comparison of the observed and filtered velocity distribution with the backtraced velocity distributions shows that indeed the observed residues correspond to the forward kinematic solution, selected by the spectrometer angular acceptance. The momentum acceptance of the measurement does not appear to play a role in the filtering procedure. The width of the backtraced velocity distribution is rather small and corresponds to the incomplete fusion scenario. The comparison of the observed and filtered N/Z-distribution with the backtraced N/Z-distributions shows a strong shift of the N/Z-ratio towards the $`\beta `$-stability line which is typical for hot heavy sources where N/Z-equilibration is achieved during de-excitation. The widths of the distributions differ only slightly and the narrow N/Z-distributions correspond to the incomplete fusion scenario.

Since the experimental effect is determined by residues originating from violent collisions it is interesting to assess to what extent the description of incomplete/complete fusion can be considered adequate. As already mentioned above, the ICF calculation  proved rather successful in quantitative description of the available data on properties of hot sources created in the mid-central to central collisions. Thus, we believe that the description of the entrance channel does not bring uncertainties influencing the conclusions of the present study and that the experimental data provide information about the details of the de-excitation stage, allowing to distinguish between different scenarios based on quantitative observables, such as the yields of heavy residues.

In order to illustrate in detail the difference in the de-excitation scenarios of sequential binary decay and statistical multifragmentation, we carried out simplified calculations. In agreement with the results of backtracing procedure, the hot nucleus 144Nd with excitation energy 310 MeV flying along the beam direction with velocity $`\beta `$=0.17 was de-excited by both GEMINI and SMM. The spin of the hot nucleus was chosen $`J`$=30, in accord with backtracing procedure. The resulting angle versus mass contour plots are presented in Fig. 13 in a fashion analogous to Fig. 5. The behavior observed in Fig. 5 is essentially reproduced. Again, the SMM leads to larger yields of residues with A$`<`$90 at angles corresponding to the MARS acceptance. The effect of MARS angular acceptance is demonstrated in Fig. 14 where the velocity versus mass contour plots are shown for the simulated residues with $`\theta `$=1-3 deg acceptance. Again the results for the heavy residues corresponding to the forward kinematic solution ( upper arm in Fig. 14 ) are analogous to the results of full simulations as shown in Figs. 6,7. The heavy residues corresponding to the backward kinematic solution ( lower arm ) were not detected in the experiment and thus are rejected by the full MARS filtering procedure for the current experiment. As suggested above, the selection caused by the angular acceptance of MARS indeed seems to be the most crucial criterion of the filtering procedure for the heavy residues investigated in the present work. In Fig. 15 both unfiltered yields and yields filtered by the MARS angular acceptance ( upper and lower panels, respectively ) for GEMINI and SMM simulations are shown. As one can see, the unfiltered yield of heavy residues with A$`<`$90 is higher in the SMM simulation and this is essentially preserved after filtering by MARS angular acceptance. The essential difference between the results of SMM and GEMINI in the simplified calculation appears to be in the production rate of heavy residues with A$`<`$90, in good agreement with the conclusions set in the discussion of Fig. 8. Using a backtracing procedure of this simplified SMM simulation it was determined that the most probable hot fragment partition contributing to the filtered yield is essentially a ternary configuration resembling the ternary fission with two heavier fragments accompanied by a lighter third particle. The same conclusion applies to backtracing of the full simulation. In general, the simplified and full simulations lead essentially to the same results, which is even more convincing when taking into account that the source characteristics have been chosen using backtracing of SMM simulation but nevertheless the result is equivalent also for the simulations using GEMINI. Based on the simplified simulation, we assure that the reason for different production rates derives from the difference of physical scenarios used by the de-excitation codes.

In the particular case of the GEMINI code it is of importance to establish to what extent an adjustment of the parameter set could be justified in order to improve the agreement of the simulation with experiment. As already illustrated in the discussion of Fig. 8, the parameter essentially determining the decay rate of channels with heavy residues with A$`<`$90 is the fission barrier height. The comparison with the recent high-quality experimental data on asymmetric fission barriers  suggests that the finite range fission barriers used in the GEMINI code are in reasonable agreement with experiment. The asymmetry dependence of the experimental barriers is reproduced well, while the experimental values are typically higher by 1-2 MeV. Thus, there is no physical justification for an enhancement of the decay rate by lowering the fission barriers. In a similar fashion, the choice of the level density parameter of Lestone  is in agreement with the systematics of effective level density parameters obtained using the experimental data on particle emission spectra . The value of the level density parameter in the saddle point configuration is typically larger ( by 3-4 % ) than in the particle evaporation channels, due to the temperature dependence of the level density parameter of Lestone. Such a calculation is essentially equivalent to the calculation with temperature indepedent level density parameter and af/an ratio according to models taking into account the surface area increase in the saddle configuration . However, several experimental works  suggest that a value of af/an=1 is good approximation at excitation energies in the saddle configuration above 40 MeV, mostly based on statistical model analysis of evaporation residue cross sections and pre-fission neutron multiplicity. In any case, the further increase of af/an ratio does not appear to be justified. Thus, the parameter set used can be considered appropriate as demonstrated in other works  where heavy residues originating from peripheral collisions have been investigated. In these deep-inelastic collisions, the details of de-excitation are overshadowed by the broad excitation energy, angle and velocity distributions.

The insight into the difference between the sequential binary and simultaneous decay can be obtained using the conclusions of the paper of Lopez and Randrup , where a statistical model of sequential binary decay is compared to a more general statistical model where a transitional state formalism is extended to the multifragment partition with any number of fragments. The authors observe that the transitional state model of multifragmentation leads also to a more general expression for the decay rate of binary channels leading to additional enhancement at high temperatures due to additional available degrees of freedom when compared to the binary decay formula used in the GEMINI code. In a similar fashion we can conclude that the decay rate of channels with higher fragment multiplicities will be enhanced in the more general model of Lopez and Randrup when compared to the corresponding de-excitation chain with binary splits, due to still higher number of available degrees of freedom. When considering, as suggested by the authors, that the transitional state model of multifragmentation presented in  is an alternative formulation to the models of nuclear disassembly such as statistical multifragmentation model , on which the SMM code is based, we can attribute the success of the SMM simulation in describing the experimental yields to the proper treatment of the available degrees of freedom in de-excitation of hot sources with the highest excitation energies created in the most violent collisions.

From the comparison of the measured yields with the results of simulations, one can conclude that the heavy residues with A$`<`$90 observed in the experiment are produced in violent collisions at low impact parameters, where a single hot source is created by incomplete fusion which further de-excites via simultaneous multifragmentation. The production of such residues can not be explained either by peripheral collisions or by fission and sequential decay of the hot source. The DIT+ICF/SMM simulation is shown to be a realistic model framework for the description of asymmetric nucleus-nucleus collisions at projectile energies around the Fermi energy, offering a quantitative description of the experimental data.

## 4 Summary and conclusions

In the present work, the cross sections and velocity distributions of heavy residues from the reaction of 20 AMeV 124Sn on 27Al have been measured at forward angles using the MARS recoil separator at Texas A&M. A consistent overall description of the measured cross sections and velocity distributions was achieved using a model calculation employing the concept of deep-inelastic transfer for the primary stage of the peripheral collisions, pre-equilibrium emission and incomplete fusion for the primary stage of the more violent central collisions and the statistical model of multifragmentation for the deexcitation stage. An alternative calculation employing the model of the sequential binary decay could not reproduce the observed yields of the residues from violent collisions (A$`<`$90) since the angular distribution of the generated residues does not extend far enough toward zero angle. The success of the statistical multifragmentation model can be attributed to a better treatment of decay rates and kinematic properties of de-excitation channels, leading to heavy residues with A$`<`$90, at excitation energies around the multifragmentation threshold, within the scenario where various multifragment partitions compete according to their statistical decay widths. The decay widths for sequential binary decay exhibit restrictions inherent to the physical picture of the asymmetric binary fission. As demonstrated in the present work, a high resolution measurement using a kinematic separator is a sensitive method allowing to investigate the details of de-excitation of hot nuclei created in violent collisions.

## 5 Acknowledgment

We gratefully acknowledge the support of the operations staff of the Cyclotron Institute during the measurements. We express our gratitude to L. Tassan-Got for the opportunity to use his DIT code, to R. Charity for the use of the GEMINI code and to A.S. Botvina for the use of the SMM code. Financial support for this work was given, in part, by the U.S. Department of Energy under Grant No. DE-FG03-93ER40773, by the Robert A. Welch Foundation under Grant No. A-1266 and by Slovak Scientific Grant Agency under Grant No. VEGA-2/1132/21. M.V. would like to express his gratitude to the Cyclotron Institute of Texas A&M University for hospitality during his stay.


format: text

text:
# Frequency range and explicit expressions for negative permittivity and permeability for an isotropic medium formed by a lattice of perfectly conducting Ω particles

## Abstract

An analytical model is presented for a rectangular lattice of isotropic scatterers with electric and magnetic resonances. Each isotropic scatterer is formed by putting appropriately 6 $`\mathrm{\Omega }`$-shaped perfectly conducting particles on the faces of a cubic unit cell. A self-consistent dispersion equation is derived and then used to calculate correctly the effective permittivity and permeability in the frequency band where the lattice can be homogenized. The frequency range in which both the effective permittivity and permeability are negative corresponds to the mini-band of backward waves within the resonant band of the individual isotropic scatterer.

PACS: 41.20.Jb, 78.20.Ci, 42.25.Bs, 42.70.Qs.

## 1 Introduction

In 1967 Veselago considered an isotropic continuous medium with both negative permittivity and negative permeability and proved that a negative refraction would occur when a monochromatic wave impinges on the surface of such a medium . As a result he obtained theoretically a quasi-lens (the lens does not focus parallel rays) from a slab of such a hypothetic material. In the Veselago theory this quasi-lens must form an image of any source without distortion (since each interface is a parallel plate). The phase of the radiated field at the focal point should be the same as the phase of the source dipole. These results are the consequences of backward waves in such a medium. The phase velocity and the Pointing vector are in opposite direction for a backward wave. In 2000 Pendry claimed that the Veselago quasi-lens must possess an extraordinary property to amplify and focus the evanescent waves in the spatial spectrum of a point source and therefore this lens must be perfect . Since then much attention has been attracted to this research area both theoretically and experimentally (see e.g. $``$). Such negative media do not exist in nature but have already got different names in the literature such as backward-wave media, the Veselago media, left-handed media and negative meta-materials. Since it is not clear yet how to make such media in the optical range of frequencies, the experimental efforts have been concentrated on creating negative meta-materials at microwave frequencies. A structure allowing the propagation of backward waves within a certain frequency band has been suggested in  and then experimentally studied in , where the negative refraction phenomenon was primarily observed. However, the structure of  (a combination of two lattices: a lattice of infinitely long parallel wires and a lattice of the so-called split-ring resonators) is not isotropic. This structure can be described with negative permittivity and permeability only if the propagation direction is orthogonal to the axes of wires .

To the best of our knowledge, the only isotropic medium with possible negative permittivity and permeability (within a certain frequency band) that has been reported is formed by many isotropic cubic cells of $`\mathrm{\Omega }`$ particles . Each isotropic cubic cell is made by putting 6 $`\mathrm{\Omega }`$-shaped perfectly conducting particles on its faces (as shown in Fig. 1). Fig. 1 shows how to make an isotropic unit cell by putting successively 3 pairs of $`\mathrm{\Omega }`$-particles on the opposite faces of a cubic unit cell in an appropriate way. Then each cubic cell can be described approximately as an isotropic resonant scatterer. The frequencies at which both the effective permittivity and permeability of such a composite medium are negative are within the resonant band of each individual isotropic scatterer. Note that conducting $`\mathrm{\Omega }`$-particles were first suggested in  for creating the bianisotropic composites.

In the present paper we consider a negative meta-material formed by a rectangular lattice of isotropic cubic unit cells of $`\mathrm{\Omega }`$ particles (as shown in Fig. 1). We introduce an analytical model to describe the dispersion properties of the lattice and calculate correctly the effective permittivity $`ϵ_{eff}(\omega )`$ and permeability $`\mu _{eff}(\omega )`$ in the frequency band where the lattice can be homogenized. We compare the results obtained by the present model with the previous results obtained by the Maxwell Garnett model, and show that the latter can not give reliable values for either the real parts or the imaginary parts of the effective material parameters $`\mu _{eff}`$ and $`ϵ_{eff}`$.

From the negative dispersion of the lattice at low frequencies we obtain the frequency range of backward waves within the resonant band of the isotropic scatterers. In this frequency range the phase velocity is in an opposite direction of the group velocity. Our calculations show that this backward wave range is exactly the same as the frequency range within which both $`ϵ_{eff}`$ and $`\mu _{eff}`$ are negative.

The case when the lattice is formed by many parallel bianisotropic $`\mathrm{\Omega }`$-particles is also considered in an appendix. It is shown that the bianisotropy of the lattice particles reduces the possibility for the existence of backward waves.

## 2 Polarizability of an isotropic cubic unit cell

For bianisotropic particles (a $`\mathrm{\Omega }`$-particle is a special case), one has the following four dyadic polarizability matrices relating the electric dipole moment $`𝐩`$ and the magnetic dipole moment $`𝐦`$ with the local electric and magnetic fields:

$`𝐩`$ $`=`$ $`\overline{\overline{a}}_{ee}𝐄^{\mathrm{loc}}+\overline{\overline{a}}_{em}𝐇^{\mathrm{loc}},`$ (1)
$`𝐦`$ $`=`$ $`\overline{\overline{a}}_{me}𝐄^{\mathrm{loc}}+\overline{\overline{a}}_{mm}𝐇^{\mathrm{loc}}.`$ (2)
For reciprocal particles, one has $`\overline{\overline{a}}_{me}=\overline{\overline{a}}_{me}^T`$.

For a cubic unit cell of $`\mathrm{\Omega }`$ particles as shown in Fig. 1, the situation is quite different. The bianisotropy cancels out since e.g. each pair of opposite $`\mathrm{\Omega }`$ particles shown in Fig. 1 have equal absolute value but opposite sign for $`\overline{\overline{a}}_{em}`$. The polarizabilities $`\overline{\overline{a}}_{ee}`$ and $`\overline{\overline{a}}_{mm}`$ of a cubic unit cell of 6 $`\mathrm{\Omega }`$-particles (see Fig. 1) are the sums of the polarizabilities $`\overline{\overline{a}}_{ee}^\mathrm{\Omega }`$ and $`\overline{\overline{a}}_{mm}^\mathrm{\Omega }`$ for individual particles . The influence of the mutual coupling of the particles is small and the resonant frequency for such a cubic unit cell of 6 $`\mathrm{\Omega }`$-particles has only a small shift from that for an isolated $`\mathrm{\Omega }`$-particle. The resonant excitation by an electric field directed along the arms of $`\mathrm{\Omega }`$ particles is mainly due to the presence of the arms. The resonant behavior retains for the total electric polarization of the two opposite omega particles. Within the resonant band, the quasi-static polarizability component $`a_{ee}^{zz\mathrm{\Omega }}`$ (which describes the response of the particle to an electric field normal to the direction of the arms; see Fig. 1(a)) of the $`\mathrm{\Omega }`$ particle is small as compared to $`a_{ee}^{yy\mathrm{\Omega }}`$. Such a cubic unit cell can be therefore considered as an isotropic scatterer and thus one has $`a_{ee}^{xx}=a_{ee}^{yy}=a_{ee}^{zz}=a_{ee}`$ and $`a_{mm}^{xx}=a_{mm}^{yy}=a_{mm}^{zz}=a_{mm}`$ . One can easily show that

$$a_{ee}2a_{ee}^{xx\mathrm{\Omega }},a_{mm}=2a_{mm}^{yy\mathrm{\Omega }}.$$
(3)
The expressions for the parameters $`a_{ee}^{xx\mathrm{\Omega }}`$ and $`a_{mm}^{yy\mathrm{\Omega }}`$ have been given in . In the present paper we consider a material formed by a rectangular lattice of such scatterers (see Fig. 1(b)).

The electric and magnetic polarizabilities for a cubic cell of $`\mathrm{\Omega }`$-particles are calculated (cf. ) and shown in Fig. 2 as a numerical example. The geometric parameters for the $`\mathrm{\Omega }`$-particles are chosen as $`r=1.5`$ mm, $`w=0.4`$ mm, $`h=0.2`$ mm and $`l=2`$ mm (see Fig. 1(a)). The size of the cubic unit cell is $`4`$ mm. The relative permittivity for the background medium is chosen as $`ϵ_b=1.5`$ (the permeability of the background medium is assumed be the same as the one for vacuum in the present paper). Resonances of $`a_{ee}`$ and $`a_{mm}`$ (frequencies at which the real parts of these parameters become zero) occur at the frequency of 8.05 GHz and 8.14 GHz, respectively. We have chosen the time dependence $`e^{j\omega t}`$, and thus the imaginary parts of these parameters are negative. Below we consider a material formed by a rectangular lattice of such isotropic scatterers (see Fig. 1(b)).

## 3 Wrong results predicted by the Maxwell Garnett model for the regular lattice

In , the well-known Maxwell Garnett model was used and the following expressions for the material parameters were obtained,

$$ϵ_{eff}=ϵ_b+\frac{1}{F}\left(\frac{Na_{ee}}{ϵ_0}\frac{N^2a_{ee}a_{mm}}{3ϵ_0ϵ_b\mu _0}\right),$$
(4)
$$\mu _{eff}=1+\frac{1}{F}\left(\frac{Na_{mm}}{\mu _0}\frac{N^2a_{ee}a_{mm}}{3ϵ_0ϵ_b\mu _0}\right),$$
(5)
where

$$F=1\frac{Na_{ee}}{3ϵ_0}\frac{Na_{mm}}{3\mu _0}\frac{N^2a_{ee}a_{mm}}{9ϵ_0ϵ_b\mu _0}.$$

In the Maxwell Garnett model, one only needs to know the scatterer polarizabilities and the density $`N`$ of the scatterers. Thus the obtained result for a regular lattice of scatterers will be the same as that for a random distribution of scatterers if the density of the scatterers is the same. For the case considered in the present paper (i.e., a rectangular lattice of isotropic cubic unit cells as shown in Fig. 1(b) with periods $`d_x,d_y,d_z`$ along the Cartesian axes), one has $`N=1/d_xd_yd_z`$. The frequency dependencies of $`ϵ_{eff}`$ and $`\mu _{eff}`$ predicted by the Maxwell Garnett model are shown in Fig. 3 when $`d_x=d_y=d_z=d=8`$ mm (i.e. $`N=(1/8^3)10^9`$ m<sup>-3</sup>). From this figure one sees that both $`\mathrm{Re}(\mu _{\mathrm{eff}})`$ and $`\mathrm{Re}(ϵ_{\mathrm{eff}})`$ are negative within a certain frequency band (8.34$``$8.62 GHz). However, within this frequency band the imaginary parts of $`\mu _{\mathrm{eff}}`$ and $`ϵ_{\mathrm{eff}}`$ have high values (about $`46`$). The radiation losses are so large that the obtained negative effective permittivity $`ϵ_{eff}`$ and permeability $`\mu _{eff}`$ have no practical importance for applications.

The above results predicted by the Maxwell Garnett model (in  for a random arrangement of cubic unit cells) are not correct for the present case of the regular lattice. Within the resonant band (7.6-8.6 GHz) of each $`\mathrm{\Omega }`$-particle, the wavelength ($`\lambda =2833`$ mm) in the background medium ($`ϵ_b=1.5`$) is much larger than the distance $`d`$ (8 mm) between the isotropic cubic scatterers, and this distance is larger than the size (4 mm) of each scatterer. From this point of view one may expect that the homogenization is possible. However, the results predicted by the Maxwell Garnett model are inconsistent near the resonance of the scatterer. The Maxwell Garnett model is based on the Clausius-Mossotti relations for the electric and magnetic fields. These relations for cubic lattices of static dipoles have been found as good approximations in the works of Sivukhin  (for an infinite lattice and for a half-space) and McPhedran  (for arrays of finite sizes). The Clausius-Mossotti relations were also confirmed to be accurate enough for cubic lattices of dipoles in a time-harmonic case  under the condition that the wavelength is large compared to the lattice period $`d`$. However, the resonant case is not considered in these works. In the resonant case, the Maxwell Garnett model may give very high absolute values for the real parts of $`ϵ_{eff}`$ and $`\mu _{eff}`$ within the resonant band, and this would lead to a dramatic shortening of the wavelength in the homogenized medium. The wavelength in the homogenized medium becomes comparable with the particle distance $`d`$ and thus the Maxwell Garnett model becomes contradictory at some frequencies.

In the present case (see Fig. 3) the absolute values of real parts of $`ϵ_{eff}`$ and $`\mu _{eff}`$ are so high within the narrow band 8.37-8.40 GHz where the Maxwell Garnett model becomes contradictory. Outside this narrow band $`|\mathrm{Re}(ϵ_{eff})|`$ and $`|\mathrm{Re}(\mu _{eff})|`$ are not very high. However, the inconsistency of the Maxwell Garnett model for the present case of regular lattice is evident over the whole resonant band of $`\mathrm{\Omega }`$-particles (8.1-8.5 GHz). In Fig. 3 one can see that the imaginary parts of $`ϵ_{eff}`$ and $`\mu _{eff}`$ are rather significant at those frequencies. On the other hand, we know that $`\mathrm{Im}(\mu _{\mathrm{eff}})`$ and $`\mathrm{Im}(ϵ_{\mathrm{eff}})`$ must be identically zeros physically if the lossless scatterers are arranged in a regular lattice. The cancelation of the radiation resistances of the scatterers (due to their electromagnetic interaction) is obvious for an infinite regular array and thus there is no radiation loss at all (see e.g. ).

In the present paper, we introduce another homogenization model which gives correctly zero value for the imaginary parts of $`ϵ_{eff}`$ and $`\mu _{eff}`$. A comparison in Section 5 will show that the real parts of $`ϵ_{eff}`$ and $`\mu _{eff}`$ predicted by the Maxwell Garnett model are also wrong within the resonant band of the scatterers.

## 4 Dispersion equation for the lattice of the isotropic scatterers

Consider the rectangular lattice of the isotropic scatterers (shown in Fig. 1) with periods $`d_x,d_y`$ along the $`x`$ and $`y`$ axes and $`d_z=d`$ along the $`z`$ axis. In the present case, the electric dipole moment is orthogonal to the magnetic dipole moment for each scatterer. Let $`𝐩`$ and $`𝐦`$ denote the dipole moments of a reference scatterer located at the origin. Thus one has

$`𝐩`$ $`=`$ $`a_{ee}𝐄^{\mathrm{loc}},`$ (6)
$`𝐦`$ $`=`$ $`a_{mm}𝐇^{\mathrm{loc}}.`$ (7)
We need to find the eigenmodes for such a lattice of isotropic scatterers. Assume the wave propagates along the $`z`$-axis. Let $`\beta `$ denote the propagation constant of an eigenmode. The first Brillouin zone is the interval $`0\beta d\pi `$. Then $`𝐩(n_x,n_y,n_z)`$ and $`𝐦(n_x,n_y,n_z)`$ (the electric and magnetic dipole moments of a scatterer located at a lattice node with coordinates $`x=n_xd_x,y=n_yd_y,z=n_zd`$) can be expressed through $`𝐩`$ and $`𝐦`$ (for the reference scatterer) as

$`𝐩(n_x,n_y,n_z)`$ $`=`$ $`𝐩e^{jn_z\beta d},`$ (8)
$`𝐦(n_x,n_y,n_z)`$ $`=`$ $`𝐦e^{jn_z\beta d}.`$ (9)
Since eigenwaves are linearly polarized in the present case, we can assume without loss of generality that $`𝐩=p𝐱_0`$ and $`𝐦=m𝐲_0`$ (see also Fig. 1). Then the local electric and magnetic fields are directed along the $`x`$ and $`y`$ directions, respectively (the local field has the same polarization as the eigenmode). This allows the following scalar expressions for the local field amplitudes in terms of $`p`$ and $`m`$:

$`E^{\mathrm{loc}}`$ $`=`$ $`Ap+Bm,`$ (10)
$`H^{\mathrm{loc}}`$ $`=`$ $`Cm+Dp,`$ (11)
where $`A,B,C,D`$ are the so-called interaction factors of the lattice, which depend only on the lattice geometry, the frequency $`\omega `$ and the propagation factor $`\beta `$ (they do not depend on the polarizabilities of the scatterers). From the reciprocity it follows that $`B=D`$. From the duality it follows that $`C=A/\eta ^2`$, where $`\eta =\sqrt{\mu _0/ϵ_bϵ_0}`$ is the wave impedance of the background medium.

From Eqs. (6),(7), (10) and (11), one obtains the following two important relations,

$`p(1a_{ee}A)`$ $`=`$ $`a_{ee}Dm,`$ (12)
$`m(1{\displaystyle \frac{a_{mm}A}{\eta ^2}})`$ $`=`$ $`a_{mm}Dp.`$ (13)
Denote the dimensionless ratio $`\eta p/m`$ as $`\alpha `$, then we have

$$\alpha \eta \frac{p}{m}=\frac{\left(\frac{\eta }{a_{mm}}\frac{A}{\eta }\right)}{D}.$$
(14)
This formula will be used in the next section to find the effective material parameters of the lattice at low frequencies.

From Eqs. (12) and (13) one obtains

$$\left(\frac{1}{a_{ee}}A\right)\left(\frac{\eta ^2}{a_{mm}}A\right)=\eta ^2D.$$
(15)
The relation

$$\mathrm{Im}\left(\frac{1}{a_{ee}}\right)=\frac{\omega ^3(ϵ_0\mu _0)^{\frac{3}{2}}ϵ_b^{\frac{1}{2}}}{6\pi ϵ_0}$$
(16)
was first obtained in  for a particle re-radiating the light and later was reproduced as the consequence of the energy conservation for a lossless dipole scatterer (see e.g. ). For lossless magnetic scatterers, one has the following corresponding form,

$$\mathrm{Im}\left(\frac{1}{a_{mm}}\right)=\frac{\omega ^3(ϵ_0ϵ_b\mu _0)^{\frac{3}{2}}}{6\pi \mu _0}.$$

For the interaction factor $`A`$, an explicit approximation has been found in :

$$A=\frac{\eta \omega }{2d_xd_y}\left(q_0+\frac{\mathrm{sin}kd}{\mathrm{cos}kdcos\beta d}\right)+j\frac{k^3}{6\pi ϵ_0ϵ_b},$$
(17)
where $`k=\omega \sqrt{ϵ_0\mu _0ϵ_b}`$ is the wavenumber in the background medium, and $`q_0`$ is the real part of the dimensionless interaction factor of a 2D grid of dipoles with periods $`d_x,d_y`$. $`q_0`$ has the following explicit expression when $`d_x=d_y=a`$ ,

$$q_0=\frac{1}{2}\left(\frac{\mathrm{cos}kas}{kas}\mathrm{sin}kas\right),$$
(18)
where $`s1/1.4380=0.6954`$. Relation (17) is almost exact when $`dd_x,d_y`$, and has only a small error when $`d=d_x=d_y`$ .

Since

$$\mathrm{Im}A=\mathrm{Im}\left(\frac{1}{a_{ee}}\right)=\mathrm{Im}\left(\frac{\eta ^2}{a_{mm}}\right),$$

the left-hand side of (15) is real. Thus we have

$$\mathrm{Re}\left(\frac{1}{a_{ee}}A\right)\mathrm{Re}\left(\frac{\eta ^2}{a_{mm}}A\right)=\eta ^2D.$$
(19)
Below we show that $`D`$ is also real and we derive its explicit expression. From the definition of $`D`$, one has

$$H_y^{\mathrm{loc}}(x=0,y=0,z=0)=\underset{n_x}{}\underset{n_y}{}\underset{n_z}{}H_y(n_xd_x,n_yd_y,n_zd)Dp,$$
(20)
where $`H_y(n_xd_x,n_yd_y,n_zd)`$ is the $`y`$-component of the magnetic field at the origin produced by the $`x`$-polarized electric dipole with dipole moment $`𝐩(n_x,n_y,n_z)=p(n_x,n_y,n_z)𝐱_0`$. From the reciprocity we know that $`H_y(n_xd_x,n_yd_y,n_zd)`$ is equal to the field at the point $`(x=n_xd_x,y=n_yd_y,z=n_zd)`$ produced by the reference dipole $`𝐩`$. In Eq. (20) the summation is over all integers $`(n_x,n_y,n_z)`$ from $`\mathrm{}`$ to $`\mathrm{}`$ except $`n_x=n_y=n_z=0`$. Here we use the following plane-wave representation for the magnetic field produced by an electric dipole $`𝐩=p𝐱_0`$ located at the origin :

$$H_y(x,y,z)=\frac{j\omega }{8\pi ^2}\mathrm{\Psi }(z)p\underset{\mathrm{}}{\overset{\mathrm{}}{}}\underset{\mathrm{}}{\overset{\mathrm{}}{}}e^{j(q_xx+q_yy)\mathrm{\Psi }(z)\sqrt{k^2q_x^2q_y^2}z}𝑑q_x𝑑q_y,$$
(21)
where $`\mathrm{\Psi }(z)=+1`$ for $`z>0`$, $`\mathrm{\Psi }(z)=1`$ for $`z<0`$ and $`\mathrm{\Psi }(z)=0`$ for $`z=0`$. Substituting Eqs. (21) and (8) into definition (20) and changing the order of the summation over $`n_z`$ and the integration over $`q_x`$ and $`q_y`$, one obtains

$$D=\frac{2j\omega }{8\pi ^2}\underset{n_x}{}\underset{n_y}{}\underset{\mathrm{}}{\overset{\mathrm{}}{}}\underset{\mathrm{}}{\overset{\mathrm{}}{}}e^{j(q_xx+q_yy)}\underset{n_z=1}{\overset{\mathrm{}}{}}e^{jn_z\beta d}e^{n_zd\sqrt{k^2q_x^2q_y^2}}dq_xdq_y.$$

Summarizing the geometrical series and using known explicit formulas for the integration over $`q_x`$ and $`q_y`$, one can obtain the following final result for $`D`$,

$$D=\frac{\omega }{2d_xd_y}\underset{n_x}{}\underset{n_y}{}\frac{\mathrm{sin}\beta d}{\mathrm{cos}\sqrt{k^2\left(\frac{2\pi n_x}{d_x}\right)^2\left(\frac{2\pi n_y}{d_y}\right)^2}d\mathrm{cos}\beta d}.$$
(22)
Numerical calculations have shown that the contribution from the terms with $`n_x,n_y=\pm 1,\pm 2,\mathrm{}`$ becomes significant only if $`d_x,d_yd`$. These terms correspond to the effect of the evanescent Floquet spatial harmonics generated by the 2D grids of dipoles with $`n_z=\pm 1,\pm 2\mathrm{}`$ (the reference dipole is located in the 2D grid with $`n_z=0`$) . If we consider the case $`d_x=d_y=ad`$, then we can neglect these terms and take into account only the plane-wave interaction between the 2D grids of dipoles. For such a case, one has the following real-valued expression for $`D`$,

$$D=\frac{\omega }{2d_xd_y}\frac{\mathrm{sin}\beta d}{\mathrm{cos}kd\mathrm{cos}\beta d}.$$
(23)
Note that relation (17) was derived in  under the same condition. Formulas (17) and (23) are both approximate but self-consistent. Therefore, these relations lead to the real-valued dispersion equation (19).

Substituting (17) and (23) into (19), one obtains

$$\mathrm{cos}^2\beta d(1+\gamma _1\gamma _2)\mathrm{cos}\beta d\left(2\gamma _1\gamma _2\mathrm{cos}kd(\gamma _1+\gamma _2)\mathrm{sin}kd\right)$$

$$(1\gamma _1\gamma _2)\mathrm{cos}^2kd(\gamma _1+\gamma _2)\mathrm{sin}kd\mathrm{cos}kd=0,$$
(24)
where

$$\gamma _1=\frac{\mathrm{Re}\left(\frac{1}{a_{ee}}\right)2d_xd_y}{\eta \omega },\gamma _2=\frac{\mathrm{Re}\left(\frac{\eta ^2}{a_{mm}}\right)2d_xd_y}{\eta \omega }.$$

The dispersion relation (24) has two roots for a fixed frequency. One of the roots satisfies the following condition (obtained from Eqs. (12)and (13))

$$\frac{\left(\frac{\eta }{a_{mm}}\frac{A}{\eta }\right)}{D}=\frac{\eta D}{\left(\frac{1}{a_{ee}}A\right)}.$$

This root is $`\beta `$ that we show in our dispersion plots below. The other root does not satisfy the above condition and thus is spurious.

Fig. 4 shows the dispersion for the same cubic lattice of isotropic scatterers as used for Fig. 3. From this figure one can see that there is a band (between 8.14 and 8.37 GHz) in which the group velocity is negative. This is the backward-wave band. Since the Poynting vector in a lossless medium must be directed along the group velocity, this band corresponds to the case when the energy transports oppositely to the phase. Of course, this situation would be trivial for some high-frequency dispersion branches. For a lattice of dipoles, negative dispersion is inherent for every even dispersion branch ($`\pi (2n+1)<kd<\pi 2(n+1)\pi `$, where $`n=0,1,2\mathrm{}`$). However, only in the first frequency zone (where $`0<kd<\pi `$) the lattice period is smaller than $`\lambda /2`$ (half of the wavelength in the background medium) and thus the homogenization of the lattice is possible. Only if the lattice of scatterers can be homogenized, the material parameters $`ϵ_{eff}`$ and $`\mu _{eff}`$ can be introduced and the lattice of scatterers can be interpreted as a continuous medium. It the Veselago theory  the backward wave in a continuous medium should correspond to a case when both $`ϵ_{eff}`$ and $`\mu _{eff}`$ are negative. The negative dispersion in the homogenized structure implies the backward wave and therefore our $`ϵ_{eff}`$ and $`\mu _{eff}`$ should be negative if the Veselago theory is applicable for dispersive media.

Fig. 5 gives an enlarged view of Fig. 4 for the frequency dependence of $`\mathrm{Re}\beta `$ and $`\mathrm{Im}\beta `$ in the frequency range 7-10 GHz. From this figure one sees that there are two stopbands within the resonant band of each scatterer (cf. Fig. 2). The lower one (stopband 1) is very narrow (8.06-8.14 GHz), and the higher (stopband 2) is wider (8.37-8.52 GHz). Between these two stopbands there is the band of backward waves. Stopband 2 is a conventional lattice stopband which has $`\mathrm{Re}\beta =0`$. Stopband 1 corresponds to the so-called complex mode which is known in the theory for lattices of infinite wires (see e.g. ,). This complex mode is decaying and the real part of the propagation factor is identically equal to $`\pi /d`$. In our case, the relation $`\mathrm{Re}\beta =\pi /d`$ reflects the fact that the directions of the dipole moments of the scatterers are alternating along the propagation axis (i.e., two adjacent isotropic scatterers have opposite polarizations).

Substituting relations (17) and (23) into Eq. (14), one obtains the following explicit formula for the parameter $`\alpha `$,

$$\alpha (\omega )=\frac{\left[q_0\frac{2d_xd_y\eta }{\omega }\mathrm{Re}\left(\frac{1}{a_{mm}}\right)\right](\mathrm{cos}kd\mathrm{cos}\beta d)+\mathrm{sin}kd}{\mathrm{sin}\beta d}.$$
(25)
Here $`a_{mm}`$ is a function of $`\omega `$ for given geometrical parameters of scatterers and given $`ϵ_b`$. Since the frequency dependence of $`\beta (\omega )`$ is known from the solution of the dispersion equation (24), the frequency dependence of $`\alpha (\omega )`$ can be calculated from the above equation. In the next section we will use this frequency dependence of $`\alpha (\omega )`$ to find the explicit expressions for the effective material parameters $`ϵ_{eff}(\omega )`$ and $`\mu _{eff}(\omega )`$ of the homogenized lattice.

Our calculations show that $`\alpha `$ is also a crucial parameter for the frequency bounds of the backward-wave range. The central frequencies of both stopbands are frequencies at which $`\alpha (\omega )`$ changes the sign. It is positive below 8.1 GHz, negative in 8.1-8.45 GHz (backward-wave range) and positive again above 8.45 GHz. Both stopbands are centered by the 2 zero-points of $`\alpha (\omega )`$. Of course, both stopbands are within the resonant band of each isotropic scatterer.

## 5 Homogenization of the lattice of the isotropic scatterers

Consider the frequency band $`0<kd\pi `$ within which homogenization is possible and effective parameters $`ϵ_{eff}`$ and $`\mu _{eff}`$ can be, perhaps, introduced. To determine the frequency dependence of the two parameters $`ϵ_{eff}(\omega )`$ and $`\mu _{eff}(\omega )`$, we needs to find two relations (involving the propagation constant $`\beta (\omega )`$) between these two parameters. The first relation can be obtained by fitting the value $`\beta /k`$ with the effective refraction index $`n_{eff}=\sqrt{ϵ_{eff}\mu _{eff}}`$, i.e., one has

$$\mu _{eff}=\frac{\beta ^2}{k^2ϵ_{eff}}.$$
(26)
Note that $`k`$ is proportional to $`\omega `$ (i.e., $`k=\omega \sqrt{ϵ_0\mu _0ϵ_b}`$) and $`\beta =\beta (\omega )`$ is given by the dispersion curve calculated in the previous section.

To find another relation, we define the following averaged (over a cubic unit cell) electric and magnetic polarizations per unit volume:

$$P=\frac{p}{V},M=\frac{m}{V},$$

where $`V=d^3`$ (in this section we assume that the lattice periods are the same in all the three directions, i.e., $`d_x=d_y=d_z=d`$; otherwise the lattice will be anisotropic). In the same way we can introduce the following averaged electric and magnetic fields,

$$<E>=\frac{1}{V}\underset{V}{}E𝑑V,<H>=\frac{1}{V}\underset{V}{}H𝑑V.$$

Since the material considered here is isotropic, its wave impedance $`Z`$ can be expressed as the ratio $`<E>/<H>`$, i.e.,

$$\frac{<E>}{<H>}Z.$$
(27)
If we use the conventional expression $`Z=\sqrt{\frac{\mu _0\mu _{eff}}{ϵ_0ϵ_{eff}}}`$, the square root can be chosen either positive (in the range of forward waves) or negative (in the range of backward waves; the ratio $`\mu _{eff}/ϵ_{eff}`$ is always positive). For a backward wave in a continuous isotropic medium, the vectors $`<𝐄>=E𝐱_0`$ and $`<𝐇>=H𝐲_0`$ form a left-hand triad with the wave vector $`𝐊=\beta 𝐳_0`$ , and thus the wave impedance should be negative. Substituting (26) into $`Z=\sqrt{\frac{\mu _0\mu _{eff}}{ϵ_0ϵ_{eff}}}`$, one obtains

$$Z=\sqrt{\frac{\mu _0}{ϵ_0}}\frac{\beta }{kϵ_{eff}}.$$
(28)
The above expression is more convenient for use since it gives automatically a negative impedance for a backward wave (when $`ϵ_{eff}<0`$).

From the definitions for the effective material parameters, one has

$$ϵ_0ϵ_{eff}<E>=ϵ_0ϵ_b<E>+P,\mu _0\mu _{eff}<H>=\mu _0<H>+M.$$

Using Eqs. (27) and (28) and taking into account the above relations, one can express the dimensionless parameter $`\alpha =\eta P/M`$ (determined by Eq. (25)) as

$$\alpha =\eta \frac{P}{M}=\eta \frac{ϵ_0(ϵ_{eff}ϵ_b)}{\mu _0(\mu _{eff}1)}\frac{<E>}{<H>}=$$

$$\eta \frac{ϵ_0(ϵ_{eff}ϵ_b)}{\mu _0(\mu _{eff}1)}Z=\frac{(ϵ_{eff}ϵ_b)}{\sqrt{ϵ_b}(\mu _{eff}1)}\frac{\beta }{kϵ_{eff}}.$$
(29)
The above equation gives the second required relation between the two effective parameters.

Finally, from the two relations (26) and (29) we obtain the following expression for determining $`ϵ_{eff}(\omega )`$ from the dispersion relation $`\beta (\omega )`$,

$$ϵ_{eff}(\omega )=\frac{\frac{\beta ^2}{k^2}+\frac{\beta }{k\alpha }\sqrt{ϵ_b}}{1+\frac{\beta }{k\alpha \sqrt{ϵ_b}}}.$$
(30)
Note that $`k=\omega \sqrt{ϵ_0\mu _0ϵ_b}`$ and $`\alpha =\alpha (\omega )`$ is given by formula (25). After $`ϵ_{eff}(\omega )`$ is calculated, $`\mu _{eff}(\omega )`$ is then calculated from Eq. (26).

Fig. 6 shows the effective permittivity and permeability for the same lattice of isotropic scatterers as used before in Figs. 3$``$5. From this figure one sees that both $`ϵ_{eff}`$ and $`\mu _{eff}`$ are negative in the band of backward waves. This is consistent to the Veselago theory.

In any of the two stopbands, the value $`\beta `$ contains non-zero imaginary part and thus Eqs. (30) and (26) give complex values for $`ϵ_{eff}`$ and $`\mu _{eff}`$ at these frequencies. Since the structure is lossless, $`ϵ_{eff}`$ and $`\mu _{eff}`$ must be real at the frequencies for which the homogenization is possible. Therefore, we can conclude that the homogenization can not be performed within the stopbands. That’s why in Fig. 6 no value is given for $`ϵ_{eff}`$ or $`\mu _{eff}`$ in the stopbands. From Fig. 6 one sees that outside the two stopbands the absolute values of $`ϵ_{eff}`$ and $`\mu _{eff}`$ are not very high. The wavelength in the homogenized medium (which is equal to $`2\pi /\beta `$) is not dramatically shortened within the passbands (in which $`\beta d<\pi `$). This indicates that the present homogenization model is self-consistent (unlike expressions (4) and (5) derived from the Maxwell Garnett model). By comparison Fig. 6 with Fig. 3, we can conclude that for a lattice of resonant particles the Maxwell Garnett model can not give reliable values for either the real parts or the imaginary parts for the effective material parameters $`(\mu _{eff})`$ and $`(ϵ_{eff})`$. The Maxwell Garnett model can only predicts qualitatively the location of the frequency band where the effective permittivity and permeability have negative real parts.

## 6 Conclusion

In the present paper, we have presented an analytical model for an isotropic negative meta-material, which is formed by a rectangular lattice of isotropic scatterers with electric and magnetic resonances. Each isotropic scatterer is formed by putting 3 pairs of $`\mathrm{\Omega }`$-particles on the opposite faces of a cubic unit cell in an appropriate way. We have derived a self-consistent dispersion equation and studied the dispersion properties of the lattice. The obtained dispersion curves are then used to calculate correctly the effective permittivity and permeability in the frequency band where the lattice can be homogenized. It is interesting to see that the dispersion curves agree well with the Veselago theory which predicts backward waves when both permittivity and permeability are negative. The Veselago theory is non-complete (it was only for monochromatic waves and did not take into account the dispersion of the medium) and some questions are left open. It has been shown in the present paper that the frequency range in which both the effective permittivity and permeability are negative corresponds to the mini-band of negative dispersion of the lattice at low frequencies within the resonant band of the individual isotropic scatterer. It has been shown that the frequency range of backward waves is determined by the edges of the two stopbands at low frequencies. We have compared the results obtained by the present model with the previous results obtained by the Maxwell Garnett model, and shown that the latter can not give reliable values for either the real part or the imaginary part for any effective material parameter ($`\mu _{eff}`$ or $`ϵ_{eff}`$) of the lattice. For lattices of resonant scatterers, the Maxwell Garnett model can only predict approximately the location of the frequency band where the effective permittivity and permeability have negative real parts.

The case when the lattice is formed by parallel bianisotropic $`\mathrm{\Omega }`$-particles has also been studied (in the appendix). It has been shown that the bianisotropy of the lattice particles suppresses those terms in the dispersion equation that can allow a backward wave propagation and thus reduces the possibility for the existence of backward waves. Therefore, the present formulas for calculating correctly $`ϵ_{eff}`$ and $`\mu _{eff}`$ for the negative meta-material from resonant isotropic scatterers (the only isotropic negative meta-material that we know) is of importance.

## Acknowledgment

The partial support of the Royal Swedish Academy of Sciences is gratefully acknowledged. The authors are also grateful to Pavel A. Belov for some helpful discussion.

## Appendix: The case of a lattice of bianisotropic particles

In this appendix, we present briefly our studies for the case when the lattice is formed by an infinite number of parallel $`\mathrm{\Omega }`$-particles. We also explain why the magnetoelectric coupling of these bianisotropic particles reduce the possibility for the existence of backward waves.

We use a simplified wire-and-loop model for each $`\mathrm{\Omega }`$-particle, in which we assume that the current flowing around the loop is uniform (i.e., the electric polarization of the loop is neglected).

If the particle is excited by a $`x`$-directed electric field $`E^{loc}`$ (see Fig. 1), then the current in the loop and the straight wire arms will be equal to the induced voltage in the straight wire portion divided by the impedance for the series connection of the loop impedance $`Z_l`$ and the straight wire impedance $`Z_w`$, i.e.,

$$I=\frac{}{Z_l+Z_w}=\frac{E^{loc}l}{Z_l+Z_w}.$$

Thus one obtains the following electric and magnetoelectric polarizabilities,

$$a_{ee}\frac{Il}{j\omega E^{loc}}=j\frac{l^2}{Z_l+Z_w},$$

$$a_{me}\frac{I\mu _0S}{E^{loc}}=\frac{\mu _0Sl}{Z_l+Z_w},$$

where $`S=\pi r^2`$ is the loop area.

Now let the particle be excited by a $`y`$-directed local magnetic field $`H^{loc}`$. An electromotive force (i.e., voltage $``$) is induced in the loop, and one has

$$I=\frac{}{Z_l+Z_w}=\frac{j\omega \mu _0SH^{loc}}{Z_l+Z_w}.$$

Thus, one has electromagnetic and magnetic polarizabilities

$$a_{em}\frac{Il}{j\omega H^{loc}}=\frac{\mu _0lS}{Z_l+Z_w},$$

$$a_{mm}\frac{I\mu _0S}{H^{loc}}=\frac{j\omega \mu _0^2S^2}{Z_l+Z_w}.$$

From the above expressions for the polarizabilities, the following identities can be obtained:

$$a_{ee}a_{mm}=a_{em}a_{me},$$
(A1)
$$a_{me}=a_{em}.$$
(A2)
The second identity is due to the reciprocity property of the particle. For such a bianisotropic particle (local electric field is $`x`$-polarized and magnetic field is $`y`$polarized), we have (cf. Eqs. (6) and (7))

$`p`$ $`=`$ $`a_{ee}E^{\mathrm{loc}}+a_{em}H^{\mathrm{loc}},`$ (A3)
$`m`$ $`=`$ $`a_{me}E^{\mathrm{loc}}+a_{mm}H^{\mathrm{loc}}.`$ (A4)
In the rectangular lattice of parallel $`\mathrm{\Omega }`$-particles, the electric dipole moment of each particle is directed along the $`x`$-axis and the magnetic dipole moment is directed along the $`y`$-axis when the eigenmode with $`𝐄=E𝐱_0`$ and $`𝐇=H𝐲_0`$ propagates along the $`z`$axis (the local field has the same polarization as the eigenmode). Substituting relations (10) and (11) into Eqs. (A3) and (A4), one obtains

$`p(1a_{ee}Aa_{em}D)`$ $`=`$ $`m(a_{ee}D+a_{em}{\displaystyle \frac{A}{\eta ^2}}),`$ (A5)
$`m(1a_{mm}{\displaystyle \frac{A}{\eta ^2}}a_{me}D)`$ $`=`$ $`p(a_{mm}D+a_{em}A).`$ (A6)
The above two equations give

$$(1a_{ee}Aa_{em}D)(1a_{mm}\frac{A}{\eta ^2}a_{me}D)=(a_{ee}D+a_{em}\frac{A}{\eta ^2})(a_{mm}D+a_{em}A).$$

Since all terms containing $`D`$ cancel out due to identities (A1) and (A2), one obtains

$$1a_{ee}Aa_{mm}\frac{A}{\eta ^2}=0,$$

which can be rewritten in the following form

$$\frac{1}{a_{ee}+\frac{a_{mm}}{\eta ^2}}=A.$$
(A7)
For a bianisotropic particle, Eq. (16) is generalized to the following relation ,

$$\mathrm{Im}\left(\frac{1}{a_{ee}+\frac{a_{mm}}{\eta ^2}}\right)=\frac{\omega ^3(ϵ_0\mu _0)^{\frac{3}{2}}}{6\pi ϵ_0ϵ_b}.$$
(A8)
Substituting Eq. (17) into Eq. (A7) and using Eq. (A8), one obtains the following dispersion relation,

$$\mathrm{cos}\beta d=\mathrm{cos}kd+\frac{\mathrm{sin}\beta d}{X},$$
(A9)
where

$$X=\frac{2d_xd_y\mathrm{Re}\left(\frac{1}{a_{ee}+\frac{a_{mm}}{\eta ^2}}q_0\right)}{\omega \eta }.$$

Eq. (A9) has the same form as the classical dispersion equation for a periodically loaded line ), which is also known as the Kronig-Pennie equation in the solid-state physics. In the known cases described by this equation (see e.g. ) the negative dispersion exists only within the even zones of frequencies ($`2\pi n<kd<\pi (2n+1),n=1,2,3,\mathrm{}`$). We did not find any physically realizable $`\mathrm{\Omega }`$-particle for which the above dispersion relation would give a negative dispersion within the first frequency zone ($`0<kd<\pi `$).

In the above analysis, the model for the $`\mathrm{\Omega }`$-particle is rough (since the electric polarization of the loop is neglected) and identity (A1) is approximate. We have studied the lattice of parallel $`\mathrm{\Omega }`$-particles using a more accurate model (as used in ). This accurate model leads to another dispersion equation (different from Eq. (A9)) which is quadratic with respect to $`\mathrm{cos}\beta d`$ and contains the electromagnetic interaction factor $`D`$. Again, we did not find any physically realizable $`\mathrm{\Omega }`$-particle for which this dispersion relation would give a negative dispersion within the first frequency zone ($`0<kd<\pi `$). The error in Eq. (A1) for real $`\mathrm{\Omega }`$particles is quite small and the terms containing $`D`$ in this dispersion equation are relatively small. Neglecting these small terms and removing the spurious root, we obtained the same result as what Eq. (A9) leads to.

Therefore, we can conclude that the bianisotropy of the lattice particles suppresses those terms in the dispersion equation that can allow a backward wave propagation at low frequencies. The possibility for the existence of backward waves is thus reduced (it was reported with another reason in  that the bianisotropy is not good for the existence of backward waves).

## Figure captions

(a) Geometry for a cubic unit cell of $`\mathrm{\Omega }`$-particles. For graphic clearness, only 1 pair of $`\mathrm{\Omega }`$-shaped perfectly conducting particles are shown on 2 opposite faces of the cubic unit cell. (b) A rectangular lattice of isotropic cubic unit cells of $`\mathrm{\Omega }`$-particles.

Left: The real (solid lines) and imaginary (dashed lines) parts of the electric polarizability (left figure) and the magnetic polarizability (right figure) for an isotropic scatterer (a cubic unit cell of $`\mathrm{\Omega }`$-particles).

The frequency dependencies of $`ϵ_{eff}`$ and $`\mu _{eff}`$ predicted (incorrectly for the regular lattice) by the Maxwell Garnett model for the negative meta-material formed by a lattice of isotropic cubic unit cells of $`\mathrm{\Omega }`$-particles (as shown in Fig. 1(b)). $`\mathrm{Re}(ϵ_{\mathrm{eff}})`$ (solid line), $`\mathrm{Re}(\mu _{\mathrm{eff}})`$ (dashed line), $`\mathrm{Im}(ϵ_{\mathrm{eff}})`$ (circles), and $`\mathrm{Im}(\mu _{\mathrm{eff}})`$ (crosses) are shown.

The dispersion (indicated by the circles) for the same lattice of isotropic resonant scatterers as used for Fig. 3. The thin lines indicates the dispersion of the homogenous background medium.

An enlarged view of Fig. 4 for the frequency dependence of $`\mathrm{Re}\beta `$ (crosses) and $`\mathrm{Im}\beta `$ (circles) near the resonance. The thin line indicates the dispersion of the homogenous background medium.

The frequency dependencies of the effective permittivity $`ϵ_{eff}`$ (circles) and the effective permeability $`\mu _{eff}`$ (crosses) for the same lattice of isotropic resonant scatterers as considered in Figs. 3$``$5.