text:
# First-principles calculation of the band offset at BaO/BaTiO3 and SrO/SrTiO3 interfaces.

## I Introduction

The search for alternative gate dielectric materials to replace silica (SiO2) in microelectronic devices is one of the grand challenges that the materials science community and the Si-based semiconductor industry are facing at the current time MRS . The rapid scaling of the physical gate lengths of Metal-Oxide-Semiconductor Field-Effect-Transistors (MOSFET) requires a concomitant rapid reduction of the gate dielectric thickness in order to preserve a high gate oxide capacitance. This can no more be accomplished by lowering the size of the SiO2 layer because, together with problems in the thickness control, the leakage current would become inacceptably high. Indeed, the leakage current from the channel to the gate is due to the direct tunneling of carriers and increases exponentially with the decrease of both the thickness of the gate dielectric and the height of the electrostatic barrier for the electrons through the gate stack. The current roadmap projection (assessed by the the International Technology Roadmap for Semiconductors, ITRS ITR ) imposes the choice of an alternative gate dielectric with a good capacitance for a thick-enough layer and its full implementation into the production line by 2005.

The properties that the new dielectric should meet are well established and have been reviewed recently by Wilk and Wallace Wilk and Wallace (2001). They can be divided into fundamental material properties, and device processing and performance properties. Amongst the material properties, we can enumerate $`(i)`$ a higher dielectric constant than amorphous silica ($`\kappa _{\mathrm{SiO}_2}=3.9`$) in order to increase the capacitance without decreasing the thickness, $`(ii)`$ large band gaps and band offsets with Si to prevent tunneling currents, $`(iii)`$ a good thermodynamic stability in contact with the Si substrate, $`(iv)`$ a good quality of the interface with the Si channel, which means a small number of electrical defects and a low midgap interface state density, and $`(v)`$ film morphology avoiding the formation of polycristalline films and grain boundaries. Amongst the device properties, we can cite $`(vi)`$ a good compatibility with metallic gate electrodes, $`(vii)`$ a compatiblity with the deposition mechanism during the fabrication process, $`(viii)`$ reliability.

Many materials satisfy some subset of the previous criteria, but the identification of a dielectric that addresses simultaneously all of the requirements is a real challenge. Investigations on oxides like Al2O3, ZrO2, HfO2, Ta2O5, Y2O3, Gd2O3, and TiO2, have thrown encouraging results in the last few years Schlom et al. (2002). Amongst the most promising candidates, ABO3 perovskite oxides (where A stands for Ba or Sr and B stands for Ti) appear in good position.

The ABO3 compounds have a dielectric constant above 300, one order of magnitude higher than the other candidates. Although they are thermodynamically unstable in direct contact with Si (they react to form titanium silicide and alkaline-earth silicate Schlom et al. (2002); Demkov (2001)), they can be grown in perfect registry with the Si substrate by means of Molecular Beam Epitaxy (MBE) when including a silicon-compatible buffer layer. On one hand, this layer must be sufficiently thick to ensure the physical separation between the substrate and the perovskite. On the other hand, it must remain thin enough to keep the benefit of the high dielectric constant of the ABO3 compound (the capacitance of the lower-$`\kappa `$ buffer layer being in series with that of the perovskite).

In the McKee-Walker process McKee et al. (1998, 2001); McK , the buffer consists in few atomic layers of AO alkaline-earth oxide that can eventually be alloyed during the growth in order to accomodate the lattice mismatch with Si. The growth of AO on Si includes the presence, at the interface, of a (sub)monolayer of ASi2 silicide so that the final structure corresponds to the sequence Si/ASi2/AO/ABO3. The epitaxy is such that $`\mathrm{ABO}_3\left(001\right)\mathrm{AO}\left(001\right)\mathrm{Si}\left(001\right)`$, and $`\mathrm{ABO}_3110\mathrm{AO}100\mathrm{Si}100`$, i. e. the ABO3 atomic planes are rotated 45o around the $`\left(001\right)`$ AO direction McKee et al. (1991). The epitaxial crystalline growth at the oxide/semiconductor interface avoids the formation of defects and ensures the continuity of the dielectric displacement McKee et al. (2001). MBE techniques allow the control of the growing sequence at the submonolayer level preventing grain-boundaries and providing a good quality interface and extremly smooth surface morphology.

First attempts to make MOSFETs including perovskite oxides have been reported recently. Using a 110 Å-thick SrTiO3 layer as the gate dielectric, Eisenbeiser et al Eisenbeiser et al. (2000) have fabricated a transistor that behaves comparably to a 8 Å-thick SiO2/Si MOSFET. The improvement in transistor performance was very satisfactory, and the leakage currents was two order of magnitudes smaller than in a similar SiO2-based device.

As it was pointed out before, the barrier height of the dielectric with respect to the Si substrate should be large enough to minimize carrier injection into the conduction band states. A large value of the Conduction Band Offset, CBO, between Si and the gate dielectric is required, and typically materials with CBO smaller than 1.0 eV are rejected for further applications. Robertson and Chen Robertson and Chen (1999), aligning the Charge Neutrality Levels (CNL) Tersoff (1984) of both semiconductors, have estimated the CBO for a Si/SrTiO3 interface to -0.14 eV (SrTiO3 below, that is no barrier at all for the electrons) in very good agreement with experimental results Chambers et al. (2001). This prevents, in principle, the use of the titanate as the gate dielectric in electronic devices. However, the presence of the buffer alkaline-earth oxide in the heterostructure was missing in their approach. In this paper, we will show that in addition to providing a physical separation between Si and the perovskite, the presence of the alkaline-earth oxide also allows to monitor efficiently the band offset.

We report a study of the properties of BaO/BaTiO3 (from now on, we will refer to this heterostructure as the Ba-interface), and SrO/SrTiO3 (Sr-interface) structures from first-principles. The method on which the simulations are based is described in Section II. In Section III, we discuss the details of the atomic structure at the interfaces. The electronic structure is presented in Sections IV, where we analyze the density of states at the junctions. In Section V we study the band-offset at the interface. Finally, in Section VI, an estimate of the band alignment of the whole whole Si/SrO/SrTiO3/Pt and Si/BaO/BaTiO3/SrRuO3 structures will be given.

## II Technicalities

Our calculations have been performed within Density Functional Theory (DFT) Hohenberg and Kohn (1964) and the Local Density Approximation (LDA) Kohn and Sham (1965). We used a Numerical Atomic Orbital (NAO) method, as it is implemented in the Siesta code Soler et al. (2002); Sánchez-Portal et al. (1997); Ordejón et al. (1996). The exchange-correlation functional was approximated using the Perdew and Zunger Perdew and Zunger (1981) parametrization of Ceperley-Alder data Ceperley and Alder (1981).

Core electrons were replaced by ab-initio norm-conserving fully-separable Kleinman and Bylander (1982) Troullier-Martin Troullier and Martins (1991) pseudopotentials. Due to the large overlap between the semi-core and valence states, the $`3s`$ and $`3p`$ electrons of Ti, $`4s`$ and $`4p`$ electrons of Sr, and $`5s`$ and $`5p`$ electrons of Ba were explicitly included in the calculation. Ti, Sr and Ba pseudopotentials were generated scalar-relativistically. The reference configuration and cutoff radii for all the atoms we used are shown in Table-1.

The one-electron Kohn-Sham eigenstates were expanded in a basis of strictly-localized Sankey and Niklewsky (1989) Numerical Atomic Orbitals Artacho et al. (1999). Basis functions were obtained by finding the eigenfunctions of the isolated atoms confined within the new soft-confinement spherical potential proposed in Ref. Junquera et al., 2001. We used single-$`\zeta `$ basis set for the semicore states of Ti, Sr and Ba, and double-$`\zeta `$ plus polarization for the valence states of all the atoms. For Sr (respectively Ba) an extra shell of $`4d`$ (respectively $`5d`$) orbitals was added. All the parameters that define the shape and the range of the basis functions for Ba, Ti and O were obtained by a variational optimization in cubic bulk BaTiO3, following the procedure described in Ref. Junquera et al., 2001. For Sr, another optimization was performed in bulk SrTiO3, frozen in the atomic orbitals of Ti and O to these previously optimized in BaTiO3 Not .

The electronic density, Hartree and exchange-correlation potentials, as well as the corresponding matrix elements between the basis orbitals, were calculated in an uniform real-space grid Soler et al. (2002). An equivalent plane wave cutoff of 200 Ry was used to represent the charge density. Once self-consistency was reached, the grid was refined (reducing the distance between grid points by half) to compute the total energy, atomic forces and stress tensor.

The integrals in reciprocal space were well converged, using in all the cases a sampling in $`\stackrel{}{k}`$ of the same quality as the $`\left(6\times 6\times 6\right)`$ Monkhorst-Pack Monkhorst and Pack (1976) mesh in bulk $`\mathrm{BaTiO}_3`$. The equivalent cutoff-length Moreno and Soler (1992) for this sampling, 13 Å, was the one employed in all simulations. This represents a large number of $`\stackrel{}{k}`$-points thought that all the materials involved in the heterojunctions are insulators. However it has been proved that this fineness is mandatory while dealing with perovskites King-Smith and Vanderbilt (1994).

Test of the performance of the Siesta method on perovskites were done in bulk BaTiO3 Ordejón et al. (2001). Lattice constants, ferroelectric distorsions, Born effective charges, and phonon dispersion curves are in very good agreement with plane waves King-Smith and Vanderbilt (1994); Ghosez et al. (1999, 1995); Zhong et al. (1994) and full potential LAPW calculations Cohen and Krakauer (1990).

## III Atomic structure at the interface

In Table 2 we report the experimental and calculated lattice parameters of the different materials involved in our heterostructures, together with the lattice mismatch with respect to the Si substrate. The misfit is defined as $`f=100\times \left(aa_{Si}\right)/a_{Si}`$, where $`a`$ and $`a_{Si}`$ are, respectively, the lattice constant of the epilayer and Si. The value of $`f`$ is positive when the epilayer is compressed and negative when it is expanded. In Table 2, we observe that the LDA produces a systematic underestimate of the lattice constant (about 1$`\%`$). Nevertheless, the correct sequence of lattice mismatch is obtained so that the calculations will reproduce the experimental strain conditions when working at the theoretical lattice constants of the substrate.

Interfaces were simulated using a supercell approximation. The basic unit cell, periodically repeated in space corresponds to the generic (AO)n/(AO-BO2)m formula, where $`n`$ and $`m`$ are respectively the number of AO oxide atomic planes and the number of ABO3 unit cells McKee et al. (2001). For even $`n`$ and odd $`m`$ (the only cases studied in this work), this structure possesses two mirror symmetry planes located on the central AO and BO2 layers.

We considered pseudomorphic heterojunctions, so that the lattice constant parallel to the plane of the interface, $`a_{}`$, is assumed to remain the same on both sides of the structure. The choice of $`a_{}`$ allows to treat implicitly the mechanical effect of the substrate, which is not included explicitely in the calculations.

To establish the notation, we will call the plane parallel to the interface the $`(x,y)`$ plane, whereas the perpendicular direction will be referred to as the $`z`$ axis.

Under the strain conditions imposed by the Si substrate, the epitaxial layers will minimize the elastic energy by elongation or compression of the lattice constant along $`z`$, $`a_{}`$. To determine its value, strain relaxations of the bulk unit cells of AO and ABO3 were performed under the constraint of fixed $`a_{}`$. Since the lattice misfit between the substrate and the epilayers is small enough to remain in the linear regime, the different values of $`a_{}`$ with respect to the in-plane lattice constant can be predicted from the Macroscopic Theory of Elasticity (MTE), and therefore an estimation of the atomic structure of the interface can be done. Following the description of Ref. de Walle and Martin, 1986, and for an interface orientation along (001):

$`a_{i,}`$ $`=`$ $`\left[1D_i\epsilon _{i,}\right]a_i`$
$`\epsilon _{i,}`$ $`=`$ $`{\displaystyle \frac{a_{i,}}{a_i}}1`$
$`D_i`$ $`=`$ $`2{\displaystyle \frac{c_{12}^i}{c_{11}^i}}`$ (1)
where $`a_i`$, $`c_{11}^i`$ and $`c_{12}^i`$ stand for, respectively, the equilibrium lattice parameter and the elastic constants of material $`i`$. Theoretical values of the elastic constants are reported in Table 3. Bulk structures from the macroscopic theory are in excellent agreement with the first-principles results, as can be drawn from the results in Table 4 (relative errors within 1 % for all the cases).

The resulting bulk tetragonal unit cells were used as the building blocks of our supercell. However, as interplanar distances in the region close to the interface are not predicted properly from MTE Peressi et al. (1993), a full relaxation of the geometry using first principles methods was needed.

For each interface, a reference ionic configuration was defined by piling up truncated bulk strained materials. Atomic coordinates were then relaxed until the maximum component of the force on any atom was smaller than 10 meV/Å. The maximum component of the stress tensor along $`z`$ was smaller than $`5\times 10^3`$ eV/Å3 for the Ba-interface and than $`7\times 10^3`$ eV/Å3 for the Sr-interface. It has been confirmed that additional relaxation of $`a_{}`$ (neglected in this work) does not produce any significant change.

In order to characterize the atomic displacements induced by the relaxation, we define $`\delta _z\left(M_i\right)`$ (respectively $`\delta _z\left(O_i\right)`$) as the displacement of the cation (respectively oxygen) along $`z`$ at layer $`i`$, with respect to the initial reference configuration. We introduce the displacement of the mean position of each atomic plane as $`\beta _i=\left[\delta _z\left(M_i\right)+\delta _z\left(O_i\right)\right]/2`$, and the change in the interplanar distance between consecutive planes $`i`$ and $`j`$ as $`\mathrm{\Delta }d_{ij}=\beta _i\beta _j`$. The rumpling parameter of layer $`i`$ describes the movement of the ions with respect to the mean position of each atomic plane and corresponds to $`\eta _i=\left[\delta _z\left(M_i\right)\delta _z\left(O_i\right)\right]/2`$. It is positive when the cation $`M_i`$ is above the oxygen, and negative otherwise.

Fig. 1 shows a schematic view of the atomic relaxations for both Ba and Sr-interfaces, when the in-plane lattice constant was constrained to the theoretical one of Si. The most important features are: $`(i)`$ a compression of the interplanar distance at the interface layer ; $`(ii)`$ the appearance of an ionic interface dipole, due to the opposite motion of the anion and the cation at the AO-layer closest to the interface – the anion moves towards the AO region whereas the cations displaces inside the ABO3 part –; $`(iii)`$ a monotonic decay of the absolute value of the rumpling parameter as a function of the distance to the AO-interface layer, where the major relaxations are localized, and $`(iv)`$ the oscillatory behaviour of the sign of $`\eta _i`$ and $`\beta _i`$ inside the perovskite from layer to layer, as it happens also in ABO3 free-standing slabs Padilla and Vanderbilt (1997); Cheng et al. (2000).

The main difference between the Ba and the Sr-heterostuctures is the magnitude of the relaxations at the interface, larger in the last case. All these conclusions are independent of the in-plane lattice constant imposed in the calculation, and show very good agreement with the results obtained using the Abinit Gonze et al. (2002) plane-wave pseudopotential code. Zimmer et al. (2002).

## IV Electronic structure at the interface

In Figure 2, we report the energy band structure along a selected high symmetry line in the first Brillouin zone for the bulk alkaline-earth oxides ($`\mathrm{\Gamma }`$X-line) and cubic bulk perovskite structures ($`\mathrm{\Gamma }`$R-line). Only bands close to the Fermi level are represented. The valence bands are mainly composed of O $`2p`$ states that, in the case of the perovskites, show significant hybridizations with Ti $`3d`$ orbitals.

All the alkaline-earth and perovskite oxides we consider are insulators (experimental gaps above 3 eV). Within the cubic symmetry (in the absence of strains) and neglecting spin-orbit couplings, the computed band gap for both BaTiO3 and SrTiO3 is indirect ($`R\mathrm{\Gamma }`$). The top-most valence and the bottom-most conduction bands are three times degenerated at these high symmetry points. Under the same conditions, the gap is also indirect in SrO, with three degenerated upper most valence bands at $`\mathrm{\Gamma }`$ and a single lowest conduction band at $`X`$, whereas BaO exhibits a direct gap, between single bands at $`X`$.

In Table 5 we report experimental and theoretical bands gaps (within LDA) for all the materials involved in our study. We see that, due to the well known DFT “band gap problem”, the theoretical values are understimated by about 50 % in each case. Nevertheless, it is usually accepted that this error can be roughly compensated by an appropriate shift of the conduction bands which should not affect the conclusions of the character of the gap reported in this Section.

A uniaxial strain along (001) lowers the symmetry of the perovskites from Pm3m to P4mm. This translates into a splitting of the top of the valence bands into a singlet and a doublet. The singlet is above (below) for a compressive (tensile) strain. For the alkaline-earth oxides, the symmetry reduces from Fm$`\overline{3}`$m to I4/mmm. The top of the valence band of SrO is therefore split but, in this case, the doublet is above (below) the singlet for a compressive (tensile) strain. Spin-orbit couplings (not considered in this work) might introduce further splittings.

Fig. 3 shows the Projected Density Of States (PDOS) on the different atoms (sum of the projections of the DOS on all the atomic orbitals of the given atom) as a function of the depth of the layer inside the material for the BaO/BaTiO3 and SrO/SrTiO3 interfaces. The main conclusions that can be extracted are as follows: $`(i)`$ the absence of any interface induced gap states clearly demonstrates the semiconductor character of the heterostructures; $`(ii)`$ the features of the PDOS on the alkaline-earth and the O atom at the interface layer (labelled as 4) are much closer to the ones displayed in bulk-ABO3 than in bulk-AO, showing a dominant ABO3 character of the interface; $`(iii)`$ the PDOS converges very quickly to the bulk properties and many of the bulk features can be recovered even at the atomic layers closest to the interface; $`(iv)`$ atomic relaxations have small effects on the shape of the PDOS, as can be seen comparing the solid and dotted curves in the figure. Only a shift in the SrO layers towards the zero energy (chosen as the top of the valence band in each case) is noticeable. This effect is a direct consequence of the relaxation-induced interface dipole discussed in Section III. The different magnitude of the dipole between the Ba and Sr-interfaces explains why the shift is almost negligeable in the Ba-heterostructure.

## V Band offset

One of the most important physical quantities that characterize the interface between semiconductors or insulators is the band offset, i. e., the relative position of the energy levels on both sides of the interface. The valence-band offset, VBO (respectively conduction-band offset, CBO) is defined as the difference between the positions of the top of the valence bands (respectively the bottom of the conduction bands) of the two materials. These band discontinuities play a fundamental role in calculating the transport properties through heterojunction devices.

The determination of these offsets from first-principles cannot be achieved from a direct comparison of the corresponding band edges in the two compounds as obtained from two independent bulk band-structure calculations. The reason is the lack of an intrinsic energy scale to refer all the energies: in a first-principles simulation, the hamiltonian eigenvalues are referred to an average of the electrostatic potential that is ill-defined for infinite systems Kleinman (1981) (it is only defined to within an arbitrary constant). Consequently, together with the eigenvalue difference, we must consider the lineup of this average between the two materials. This potential shift depends on the dipole induced by the electronic charge transferred from one part of the interface to the other after the interfacial hybridization (the electronic charge density of each system will decay into the other in an, in principle, unknown way). The transfer of charge depends not only on the materials that constitute the interface, but also on the particular orientation, so the lineup can only be obtained from a self-consistent calculation on a supercell including both materials.

Therefore, from the theoretical point of view, the band offsets (BO) are usually split into two terms Colombo et al. (1991); Peressi et al. (1998):

$`\mathrm{BO}`$ $`=`$ $`\mathrm{\Delta }E_{v,c}+\mathrm{\Delta }V`$ (2)
The first contribution, $`\mathrm{\Delta }E_v`$ (resp. $`\mathrm{\Delta }E_c`$), is referred to as the band-structure term. It is defined as the difference between the top (resp. bottom) of the valence (resp. conduction) bands as obtained from two independent standard bulk band-structure calculations at the same strained geometries as in the supercell construction. Within LDA, only a first estimate of the band-structure term can be obtained, $`\mathrm{\Delta }E_{v,c}^{LDA}`$. To get more accurate results, a correction dealing with many-body effects in the quasiparticle spectra should be added:

$$\mathrm{\Delta }E_{v,c}=\mathrm{\Delta }E_{v,c}^{LDA}+\mathrm{\Delta }E_{v,c}^{corr}$$
(3)
Self-energies corrections are often obtained within the GW approximation Hybertsen and Louie (1986). They strongly modify the description of the conduction bands, and tend to solve the “band gap problem” mentioned in Section IV. Even the valence band energies might be subject of certain errors, specially in oxides Kralik et al. (1998). Unfortunately, no accurate GW data are currently available for AO and ABO3 compounds. Only model GW calculations have been performed recently for SrO and SrTiO3 and with limited success Capellini et al. (2000). To overcome the problem, we make the approximation that the errors in the valence bands are smaller than those for the conduction bands and of the same order of magnitude for the two compounds taking part in the heterostructures so that they tend to cancel each other ($`\mathrm{\Delta }E_v^{corr}=0`$). Knowing the relative position of the valence bands, we simply add the experimental band gaps (see Table 5) to obtain the discontinuities for the conduction bands ($`\mathrm{\Delta }E_c=\mathrm{\Delta }E_v^{LDA}+\mathrm{\Delta }E_{gap}^{expt}`$).

The second term, $`\mathrm{\Delta }V`$, is the lineup of the average of the electrostatic potential through the heterojunction. This macroscopic quantity summarizes all the intrinsic interface effects, such as the chemical composition, structural details and orientation. To obtain it, we start from the total (ionic plus electronic) microscopic electrostatic Hartree potential, output of the self-consistent supercell calculation (in this Section we will define the zero-energy level as the average of this potential in the unit cell del ). Then, we apply the double-macroscopic average Baldereschi et al. (1988); Colombo et al. (1991); Peressi et al. (1998) technique. It consists of performing first the average of the electrostatic potential over planes parallel to the interface, and then averaging the obtained quasiperiodic one-dimensional function with two step-like filter functions whose lenghts, $`l_1`$ and $`l_2`$, are determined by the periodicity of the constituents. Here $`l_1`$ and $`l_2`$ have been set up to the distance between equivalent AO and TiO2 planes in the alkaline-earth oxide and in the perovskite respectively. A full description of the method to the AO/ABO3 heterostructures can be found in Ref. Junquera and Ordejón, 2002. The resulting profile of the macroscopic potential is flat on both sides far enough from the interface (bulk-like regions). $`\mathrm{\Delta }V`$ is defined as the difference between these two plateau values (see Fig. 4). The lineup should be independent of the length used in the filter functions. However we have checked how doubling the size of the step-like functions introduce a numerical uncertainty in $`\mathrm{\Delta }V`$ of the order of 30 meV. This is the main source of inaccuracy in our calculations of the band offsets.

It is worth noticing that neither $`\mathrm{\Delta }E_v`$ nor $`\mathrm{\Delta }V`$ have any physical meaning by their own, being pseudopotential dependent numbers. Only the sum of both is physically significant and quite independent of the choice of the pseudopotential Peressi et al. (1998).

### V.1 First principles results

Fig. 4 shows a schematic representation of the band structure discontinuities for the BaO/BaTiO3 and SrO/SrTiO3 heterostructures, both of them calculated fixing the in-plane lattice constant at the theoretical one of Si. Band splittings steming from strain are taken into account in the figure. We use the same sign convention as Van de Walle and Martin in Ref. de Walle and Martin, 1987: a positive value of the band offset for the discontinuity at a junction A/B corresponds to an upward step in going from A to B.

From the figure, we conclude that Ba-interface is type-II, with both the valence and conduction bands of BaTiO3 falling in energy below the corresponding ones of BaO. Within the accuracy of our calculations, BaO and BaTiO3 topmost valence bands are almost aligned (an offset of only -0.06 eV is predicted), so the barrier in the conduction bands is mainly due to the difference in the band gaps of both materials and, inferred from the experimental values, amounts to -1.62 eV.

Sr-interface is type-I, meaning that the band gap of SrTiO3 lies completely inside the gap of SrO. An upward step of +0.18 eV for the valence bands is theoretically expected, which implies a CBO of -2.22 eV.

A rough estimate of the valence band offset was already accessible Rao et al. (1997) by identifying in Fig. 3 the position of the top of the valence band in the PDOS for the O atom at both symmetry planes, in the bulk-like regions of the materials that constitute the interface (numbered as 1 for AO and 9 for ABO3). The values deduced from the figure are -0.09 eV for the Ba-interface (BaO above) and +0.28 eV for the Sr-interface (SrTiO3 above), close to those obtained using the macroscopic average technique. However these numbers must be taken with care Peressi et al. (1998): this method to compute band offsets requires calculations with a higher number of special $`\stackrel{}{k}`$-points than those needed to converge the charge density or the potential lineup.

It is important to note here the crucial role played by the atomic relaxations at these polar interfaces. As was pointed out in section III, after the relaxation process an extra dipole appears at the junction that modifies the electrostatic lineup accross the interface Martin and Kunc (1981) and, consequently, the band offsets :

$`\delta \left(\mathrm{\Delta }V\right)`$ $`=`$ $`{\displaystyle \frac{4\pi }{a_{}^2}}{\displaystyle \underset{\kappa \alpha }{}}{\displaystyle \frac{Z_{\kappa ,\alpha z}^{(T)}}{ϵ_{\mathrm{}}}}\mathrm{\Delta }u_{\kappa \alpha }`$ (4)
where $`\delta \left(\mathrm{\Delta }V\right)`$ is the change in the electrostatic lineup along $`z`$-direction due to the atomic displacements, $`Z_{\kappa ,\alpha \beta }^{(T)}`$ is the Born effective tensor of atom $`\kappa `$, $`\mathrm{\Delta }u_{\kappa \alpha }`$ its displacement along cartesian direction $`\alpha `$ during the relaxation and $`ϵ_{\mathrm{}}`$ the optical dielectric constant. Looking at the magnitude of the atomic displacements, it is reasonable that the change should be more remarkable for the Sr-interface than for the Ba one. From our ab-initio calculations, and for the same supercell used to get results in Fig. 4, we observe a change in $`\mathrm{\Delta }V`$ of -0.67 eV for the Sr-interface (from 1.16 eV for the unrelaxed geometry to 0.49 eV after the relaxation), whereas in the Ba-interface the deviation amounts to -0.11 eV (from 0.44 eV to 0.33 eV). This emphasizes the importance of performing accurate first-principles atomic relaxations for correct predictions of the barriers.

To what extent do these discontinuities change with the in-plane lattice constant? This is an important question because a dependence with strain would allow us to tune the band offsets (for example, replacing the Si substrate by Ge McKee et al. (2001) in order to impose a different lattice parameter throughout the interface) depending on the required values for a given device. In order to check this point, we have carried out calculations at different in-plane lattice constants. In Table 6 and Fig. 5 we summarize the results for both, Ba and Sr-heterostructures. In both cases a variation by about 0.1 eV in VBO with the in-plane lattice constant is observed, mainly due to the band-structure term (consequence of the strain-induced splittings of the top valence-band manifold), as it happens for other lattice-mismatched, isovalent, common anion interfaces Ventra et al. (1996a). The change is almost linear, and tends to lower the energy of the valence bands of the ABO3 perovskite with respect to the AO alkaline-earth oxide.

The band-structure term displays a linear behaviour with strain for the Sr-interface (see Fig. 5-panel a). The anomalous behaviour of $`\mathrm{\Delta }E_v`$ for the Ba-interface is due to a modification in the character of the top of the valence band of BaO under strain. It changes from $`X`$ when BaO is compressed to $`Z`$ when it is expanded. This transformation occurs for a lattice constant around 5.43 $`\mathrm{\AA }`$ (theoretical lattice parameter of BaO). In Fig. 5(a) we plot the difference between the top of the valence band of BaTiO3 and the highest occupied state at $`X`$ and $`Z`$ of BaO. The crossing point is clearly identified in the figure. No extra changes in the linear beaviour of $`\mathrm{\Delta }E_v`$ are expected for longer lattice constants.

The almost-linear change in the lineup term can be explained according to an analytic scaling law proposed in Ref. Ventra et al., 1996b. Once $`\mathrm{\Delta }V`$ is known for a reference configuration with an in-plane lattice constant $`a_{}`$, then, supposing an uniform strain throughout the structure, $`\mathrm{\Delta }V^{}`$ for any other strained configuration $`a_{}^{^{}}`$ can be extrapolated from:

$$\mathrm{\Delta }V^{^{}}\frac{1}{\left(1+\epsilon _{}^2\right)}\mathrm{\Delta }V$$
(5)
where $`\epsilon _{}=(a_{}^{^{}}/a_{}1)`$. Fig. 5(b) shows a comparison of the first-principles and extrapolated values, where the heterostructure at the in-plane lattice constant of Si has been chosen as the reference configuration. Results are in good agreement within the numerical accuracy of the ab-initio results.

In summary, the VBO varies almost linearly for a large range of in-plane strains. The only deviation is observed for BaO/BaTiO3 and is explained by a change of character of the BaO gap under compression.

## VI Interface with Si

As it was pointed out in the Introduction, AO/ABO3 interface is only a part of the gate stack of technological interest for the semiconductor industry. AO acts as a buffer layer between the Si substrate and the high-$`\kappa `$ perovskite. The whole heterostructure epitaxially grown following the McKee-Walker process McK  is made of Si/ASi2/AO/ABO3 McKee et al. (1998). As it will be emphasized in this Section, the role of the buffer layer is not only the passivation of the Si substrate, but also the efficient tuning of the offsets between the perovskite and the channel.

Combining our results with various data available in the literature, we can estimate the band discontinuities along the whole heterostructures of technological interest as summarized in Fig. 6. Previous theoretical calculations of the band offsets between the alkaline-earth oxide AO and Si have been reported recently (Si/BaO Gul , Si/SrO Buongiorno-Nardelli et al. (2002). In addition, we can find in the literature theoretical estimations for the Schottky barriers between perovskites and prototypical metallic electrodes (SrTiO3/Pt Robertson and Chen (1999), or BaTiO3/SrRuO3 Jun ).

We must notice that, altough most of the previous works (except Ref. Robertson and Chen, 1999) have been done using the same basic approaches (DFT, pseudopotentials, $`\stackrel{}{k}`$-point samplings, supercells, etc), they differ in the details of the calculations (for example, they have been carried out at different in-plane lattice constant, and the size of the supercells, or the parameters used to generate the pseudopotentials might change from one to the other). So, only a rough estimate of the barriers can be deduced from the comparison and any quantitative conclusion is beyond the scope of this Section.

Within LDA, plus GW corrections, Boungiorno-Nardelli et al. Buongiorno-Nardelli et al. (2002) investigated the structural and electronic properties of the Si/SrSi2/SrO interface. They predicted a VBO between Si and SrO of -1.91 eV for the most stable interface configuration. Using the experimental gaps to locate the conduction bands, it translates in a CBO of 2.69 eV. The Schottky barrier $`\varphi _n`$ (difference between the Fermi level and the bottom of the conduction band) between SrTiO3 and Pt has been evaluated Robertson and Chen (1999) to -0.89 eV, which implies a barrier height $`\varphi _p`$ (difference between the Fermi level and the top of the valence band) of 2.41 eV. These results are summarized in the first panel of Figure 6.

Through first-principles gradient-density-functional calculations Gulleri et al. Gul  focused on the structure, energetic and band offsets of the Si/BaO interface. For the favoured geometry, they obtained a VBO of -2.0 eV. Some of us evaluated the Schottky barriers between BaTiO3 and SrRuO3 (a typical metallic electrode in ferroelectric devices Eom et al. (1992)) to be equal to $`\varphi _p`$ = +1.44 eV and $`\varphi _n`$ = -1.76 eV. These results are summarized in the second panel of Figure 6.

For both stacks, we can clearly see how the problem of the large carrier injection (expected for the perovskite in direct contact with Si Robertson and Chen (1999)) is overcome by the use of the AO layer. The electrostatic barriers for both the electrons and holes, between the gate electrode and the channel are large enough to prevent carrier injections and to push the use of ABO3 perovskites to a prominent position to replace silica as the gate dielectric oxide in MOSFETs.

## VII Conclusions

We have studied structural and electronic properties of BaO/BaTiO3 and SrO/SrTiO3 interfaces from first-principles. Atomic relaxations have been performed. Interface dipoles, due to the opposite motion of the anion and cation atoms at the interface, appear for both heterostructures. No interface electronic states are induced in the band gap. The character of the AO layer at the interface is mainly perovskite-like. Under the experimental strain conditions, the valence bands of BaO and BaTiO3 are almost aligned (within the accuracy of our calculations), whereas a slightly larger barrier is predicted for SrO/SrTiO3. Interface dipoles, induced by atomic relaxations, have a strong effect on the band alignments at the interface. A nearly linear variation of the VBO with in-plane strain is observed.

Gathering together our results and various data available in the literature, we make a guess for the band alignment of whole Si/SrO/SrTiO3/Pt and Si/BaO/BaTiO3/SrRuO3 structures. In both cases large enough electrostatic barriers for electrons and holes between the gate electrode and the channel are estimated, preventing the injection of carriers and suggesting that both perovskites compounds are promising candidates to replace silica in MOSFETs. Our results should be confirmed by more accurate calculations for the whole heterostructure.

## VIII Acknowledgments

The authors are indebted with Maria Peressi, Luciano Colombo, Rodney McKee, Marco Buongiorno-Nardelli, Hermann Kohlstedt, Darrel Schlom, Alex Demkov, and Jun Wang for useful discussions. Motorola PSRL supported the initial stages of this work. This work was supported by the Volkswagen-Stiftung (www.volkswagenstiftung.de) within the program “Complex Materials: Cooperative Projects of the Natural, Engineering, and Biosciences” with the title: “Nano-sized ferroelectric Hybrids” under project number I/77 737. Three of us (J. J., M. Z., and Ph. G.) acknowledge support from FNRS-Belgium (grants 9.4539.00 and 2.4562.03) and the Université de Liège (Impulsion grant). J. J. and P. O. acknowledge finantial support from the Fundación Ramón Areces and Spanish MCyT grant BFM2000-1312. Calculations have been performed on the NIC computer at the Université de Liège, and on computational resources from CESCA and CEPBA.


format: text

text:
# Untitled Document

Odd Order Pandiagonal Latin and Magic Cubes in Three and Four Dimensions
Solomon Gartenhaus, Purdue University, W. Lafayette, IN 47907

garten@physics.purdue.edu

I Introduction

By a magic square of order $`n`$ is here meant an arrangement, without repeats, of the integers $`\{0,1,2,\mathrm{}n^2`$ \- 1} into the $`n^2`$ cells of an $`n\times n`$ square in a way that the sum of the elements of each row, of each column and of each of the two diagonals is the same. Since the sum of $`1+2+\mathrm{}+(n^21)`$ is $`n^2(n^2`$-1)/2 and since for a magic square this must be equal to the sum of the sum of the $`n`$-integers in each of the $`n`$-rows (or columns) it follows that the common sum, $`\sigma _2`$ of these integers must be

$$\sigma _2=\frac{n(n^21)}{2}$$
$`(1.1)`$
Figure I shows two magic squares for $`n`$ = 4 and $`n`$ = 5

$$\begin{array}{ccccccccc}\hfill 0& 10\hfill & 15\hfill & \hfill 5& 13\hfill & 19\hfill & 20\hfill & 1\hfill & \hfill 7\\ \hfill 7& 13\hfill & 8\hfill & \hfill 2& 21\hfill & 2\hfill & 8\hfill & 14\hfill & \hfill 15\\ \hfill 9& 3\hfill & 6\hfill & \hfill 12& 9\hfill & 10\hfill & 16\hfill & 22\hfill & \hfill 3\\ \hfill 14& 4\hfill & 1\hfill & \hfill 11& 17\hfill & 23\hfill & 4\hfill & 5\hfill & \hfill 11\\ & & & & 0\hfill & 6\hfill & 12\hfill & 18\hfill & \hfill 24\\ & (a)\hfill & & & & & (b)\hfill & & \end{array}$$

Figure I

with sum values 30 and 60, respectively. Of particular interest is the magic square in Ib. Not only do the sums along the rows, the columns and the two diagonals have the value 60, but so do the eight broken diagonals that are obtained if all partial diagonals are completed by imagining the square bent around into a cylinder. Examples of broken diagonals in Figure Ib are: 21, 10, 4, 18, 7; 9, 23, 12, 15, 1; 17, 6, 3, 14, 20; 11, 22, 8, 19, 0. Magic $`n\times n`$ squares, such as this one, where the sums along the $`n`$-rows and columns, the two diagonals and the $`2n2`$ broken diagonals all have the same value are called pandiagonal or “diabolic” magic squares \[9, ch.$`10`$\],. To simplify, we shall use the term ROW in all of the following to mean any one of a row, a column, a diagonal or a broken diagonal. Thus in a pandiagonal magic square the sum along any ROW is the same and given by (1.1).

By a magic pandiagonal cube of order $`n`$ we shall mean an arrangement of the integers $`0,1,2,\mathrm{}n^3`$ \- 1 on the lattice points of an $`n\times n\times n`$ cube so that the 3n squares parallel to the faces of the cube as well as the six “diagonal” squares which bisect the cube and contain its 4 body-diagonals are all pandiagonal magic squares. Thus in each of the $`3n`$ \+ 6 (overlapping) squares contained in the cube, the sums along all of the the ROWS—which here includes also the files—are the same. It is easily seen that in this case the common sum $`\sigma _3`$ must be

$$\sigma _3=\frac{n(n^31)}{2}$$
$`(1.2)`$
By extension, we define a four dimensional (4-D) pandiagonal magic $`n\times n\times n\times n`$ cube as one in which the integers 0, 1, 2,…$`n^41`$ are placed, without repeats, at the $`n^4`$ lattice sites of a 4-D hypercube, so that the three dimensional (overlapping) $`n\times n\times n`$ cubes that can be formed within it, are pandiagonal magic cubes. The sum along each ROW of a 4-D cube is easily shown to be

$$\sigma _4=n(n^41)/2$$
$`(1.3)`$
In constructing pandiagonal magic squares, cubes, and hypercubes, we follow the idea of Euler —who was concerned exclusively with magic squares—and used latin squares (LS) in their construction. As defined by Denes and Keedwell  and by Laywine and Mullen  a latin square of order $`n`$ is an $`n\times n`$ array each of whose lattice points is occupied by one of $`n`$ given symbols, in a way so that no row or column contains any one of these symbols more than once. The number of LS’s grows rapidly with $`n`$ and for example for $`n=10`$ and 15 this number is $`10^{36}`$ and $`1.5\times 10^{86}`$, respectively \[9, p.$`5`$\]. Here we are interested in the relatively small subset of the LS’s, the pandiagonal latin squares, for which in addition each diagonal as well as each of the $`2n2`$ broken diagonals—in short it’s ROWS—also contains each of the $`n`$ elements precisely once.

For a general algebraic theory of pandiagonal (diabolic) magic squares see the analysis by Rosser and Walker . The books by Andrews , Kraitchik  and Rause Ball and Coxeter  contain a more empirical approach to magic squares and include various examples and practical rules for constructing magic squares. Pasles  describes some unusual magic squares constructed by Benjamin Franklin.

Some 5 $`\times `$ 5 pandiagonal magic squares and some of their important features are described by Gardner . In a companion column  he details a magic cube of order seven and notes its non-pandiagonal nature. This is reenforced by the results of Wynne  who studied magic cubes of order seven and showed that even if every square of the cube that is parallel to a cube-face is pandiagonal, not all of the six, diagonal squares of the cube can be pandiagonal. This is consistent with the present analysis, according to which only for $`n11`$ do magic pandiagonal cubes exist. (Similarly we find that only for $`n17`$ can 4-D pandiagonal magic cubes exist.) That non-pandiagonal magic cubes of order 7 do exist, however, is established by Andrews  and by Alspach and Heinrich , the latter, incidental to their discussion of cubes of order 4m.

II. Odd-$`n`$ Pandiagonal Latin Squares and Cubes in 3 and 4 Dimensions

A. Pandiagonal Latin Squares

To set the stage for our discussion of magic pandiagonal cubes in 3 and 4 dimensions, in this section we collect—and in some instances amplify on—properties of pandiagonal latin squares ,.

A latin square of odd order $`n`$, LS, is an $`n\times n`$ array involving $`n`$ distinct symbols with each symbol appearing once and only once in each row and column. We shall invariably take the $`n`$ symbols to be the integers comprising the set $`S`$ defined by

$$S=\{0,1,2,\mathrm{},n1\}$$
$`(2.1)`$
Also useful will be the set $`\overline{S}`$ defined by

$$\overline{S}=\{1,2,\mathrm{},n1\}$$
$`(2.2)`$
A diagonal LS has the additional property that each of its two diagonals also contains the $`n`$ chosen symbols exactly once. Finally a pandiagonal LS is a diagonal LS in which, in addition, each of the $`n`$ symbols appears once and only once in each of the $`2n2`$ “broken diagonals”, as defined in the preceding section.

Figure II shows 4 LS’s of order 3,4,5, and 5, respectively. The $`3\times 3`$ LS in part (a) is not a diagonal LS by virtue of the three zeros

$$\begin{array}{c}\begin{array}{ccc}0& 1& 2\\ 2& 0& 1\\ 1& 2& 0\\ & (a)& \end{array}\begin{array}{cccc}0& 1& 2& 3\\ 2& 3& 0& 1\\ 3& 2& 1& 0\\ 1& 0& 3& 2\\ & (b)& & \end{array}\\ \begin{array}{ccccc}3& 4& 0& 1& 2\\ 1& 2& 3& 4& 0\\ 4& 0& 1& 2& 3\\ 2& 3& 4& 0& 1\\ 0& 1& 2& 3& 4\\ & & (c)& & \end{array}\begin{array}{ccccc}4& 0& 1& 2& 3\\ 2& 3& 4& 0& 1\\ 0& 1& 2& 3& 4\\ 3& 4& 0& 1& 2\\ 1& 2& 3& 4& 0\\ & & (d)& & \end{array}\end{array}$$

Figure II

appearing in the upper left to lower right diagonal. By contrast the $`4\times 4`$ LS in part (b) is diagonal but not pandiagonal while the $`5\times 5`$, LS in part (c) is pandiagonal. Note that, given a pandiagonal LS, if we permute its symbols $`\{0,1,2,\mathrm{}n1\}`$ it remains pandiagonal. For example, if in IIc, we carry out the permutation 0 $``$ 1; 1 $``$ 2; 2 $``$ 3, 3 $``$ 4; 4 $``$ 0, the resultant LS remains pandiagonal and is given by Figure IId. Thus a given $`n\times n`$ pandiagonal LS is the basis for $`n`$! different ones that result from the $`n!`$ possible permutations of the $`n`$-symbols.

Consider for odd $`n`$, the $`n\times n`$ array $`L`$ with elements $`L_{ij}`$

$$L_{ij}\alpha _1i+\alpha _2j(modn)$$
$`(2.3)`$
where $`i,j`$ run over the elements of $`S`$ as do the elements $`L_{ij}`$ themselves and $`\alpha _1`$ and $`\alpha _2`$ are non-zero positive integer parameters and thus are elements of $`\overline{S}`$ in (2.2). When we wish to stress the dependence of $`L`$ on $`\alpha _1`$ and $`\alpha _2`$, we shall also use the notation $`L=<\alpha _1,\alpha _2>`$.

We now establish:

Theorem (2.1): If $`L=<\alpha _1,\alpha _2>`$ is the array in (2.3) and the greatest common divisor with $`n`$ of each of $`\alpha _1,\alpha _2,\alpha _1\pm \alpha _2`$ is 1, that is

$$(\alpha _i,n)=1;i=1,2$$
$`(2.4a)`$
$$(\alpha _1+\alpha _2,n)=1$$
$`(2.4b)`$
$$(\alpha _1\alpha _2,n)=1$$
$`(2.4c)`$
then $`L`$ is a pandiagonal LS for $`n5`$. (Note that (2.4) cannot be satisfied for even $`n`$.)

Proof: Firstly, since for $`n=3`$, the possible $`\alpha `$-values are 1 and 2 and since these violate (2.4b), it follows that (2.3) is not a pandiagonal LS for $`n=3`$. On the other hand, for $`n=5`$, the pairs $`[\alpha _1,\alpha _2]=[1,2]`$ and $`[1,3]`$, for example, do satisfy each of (2.4).

Secondly to establish first that $`L`$ in 2.3 for $`n5`$ is a diagonal LS, consider it for fixed $`j`$, say, as $`i`$ ranges over the $`j`$th row. As $`i`$ thus varies from 0 to $`n1`$, the $`j`$th row of $`L`$ varies over the same set by virtue of the hypothesis $`(\alpha _1,n)=1`$. Similarly for fixed $`i`$, as $`j`$ varies over the $`i`$th column, since $`(\alpha _2,n)=1`$, $`L_{ij}`$ varies over the same set S. Thus the rows and columns of $`L`$ satisfy the condition that $`L`$ be a latin square. Along the diagonal from $`(i,j)=(0,0)`$ to $`(n1,n1)`$, $`i=j`$, so that along here $`L_{ij}(\alpha _1+\alpha _2)i`$ ($`mod`$ $`n`$). Thus again as $`i`$ varies over the set S the $`n`$ diagonal elements of $`L`$ must be some permutation of S since $`(\alpha _1+\alpha _2,n)=1`$ according to (2.4b). Finally along the other the diagonal from $`(i,j)=(0,n1)`$ to $`(n1,0)`$, $`i+j=n1`$, so that $`L_{ij}(\alpha _1\alpha _2)i+\alpha _2(n1)`$ mod $`n`$ which has the same property by virtue of (2.4c), $`(\alpha _1\alpha _2,n)=1`$. Thus we have established that $`L`$ in (2.3) under the constraints in (2.4) is a diagonal LS; it remains only to establish that it is also pandiagonal.

To this end, consider, for example, the “split” diagonal just above—and parallel to—the lower left to upper right diagonal of $`L`$ and its appendage in the lower right hand corner, $`(i,j)=(n1,0)`$. (We assume $`i`$ increases from 0 to the right and $`j`$ increases upward from 0.) Its entries are defined by $`j=i+1`$, $`i=0,1,2,\mathrm{}n2`$ and $`(i,j)=(n1,0)`$. Substituting $`j=i+1`$ into (2.3) we find

$$L_{i,i+1}(\alpha _1+\alpha _2)i+\alpha _2(modn)$$

If we now let $`i`$ run over the complete set S, it is easy to see that since $`L_{in}L_{io}`$ $`(modn)`$ we obtain both parts of the broken diagonal and that no two of these $`n`$-elements are the same since they are simply a permutation of the elements of S by virtue of (2.4b). This argument is easily repeated for all broken diagonals and we conclude that each broken diagonal consists of a permutation of the elements of S. We have thus established the theorem.

Remark 2.1: If $`\alpha _1,\alpha _2`$ satisfy each of (2.4), then so do $`k\alpha _1,k\alpha _2`$ for any $`kϵ\overline{S}`$ for which $`(k,n)=1`$.

Remark 2.2: Reference to (2.3) shows that if $`L=<\alpha _1,\alpha _2>`$ is a pandiagonal LS, then so is $`<k\alpha _1,k\alpha _2>`$ for any positive integer $`kϵ\overline{S}`$ provided $`(k,n)=1`$. For according to (2.3), $`<k\alpha _1,k\alpha _2>`$, is simply $`<\alpha _1,\alpha _2>`$ with its elements permuted in some way.

Remark 2.3: If we add an integer $`xϵ\overline{S}`$ to the right hand side of (2.3), we obtain $`<\alpha _1,\alpha _2>`$ with its elements permuted in some way.

An important notion relating to LS’s is that of the orthogonality \[9, ch. 2\]. Two $`n\times n`$, LS’s are said to be orthogonal, if when they are superposed, none of the $`n^2`$ ordered pairs of elements that result occurs more than once. Thus if $`L_{ij}^{(1)}`$ and $`L_{ij}^{(2)}`$ $`ϵS`$ are any corresponding elements of the two LS’s, $`L^{(1)}`$ and $`L^{(2)}`$ then they are orthogonal if and only if the ordered pairs $`[L_{ij}^{(1)},L_{ij}^{(2)}]`$ and $`[L_{\alpha \beta }^{(1)},L_{\alpha \beta }^{(2)}]`$ differ $`(modn)`$ for any fixed values of $`i,jϵS`$ but for all choices of $`\alpha ,\beta ,ϵS`$.

Remark 2.4 If two LS’s are orthogonal, they remain so when the elements of either or both undergo arbitrary permutations.

We now establish:

Theorem (2.2) Let $`L^{(1)}=<\alpha _1,\alpha _2>`$ and $`L^{(2)}=<\beta _1,\beta _2>`$ be two pandiagonal LS’s with both pairs $`\alpha _1,\alpha _2`$ and $`\beta _1,\beta _2`$ each satisfying (2.4). If the determinant $`d_2`$ defined by

$$d_2=|\begin{array}{cc}\alpha _1& \alpha _2\\ \beta _1& \beta _2\end{array}|$$
$`(2.5)`$
is relatively prime to $`n`$, i.e.

$$(d_2,n)=1$$
$`(2.6)`$
then $`L^{(1)}`$ and $`L^{(2)}`$ are orthogonal pandiagonal LS’s.

Proof: Suppose for $`i,jϵS`$ there is a second pair $`k,\mathrm{}ϵS`$ for which the ordered pairs $`[L_ij^{(1)},L_ij^{(2)}][L_k\mathrm{}^{(1)},L_k\mathrm{}^{(2)}]`$ mod $`(n)`$. That is, suppose

$$\begin{array}{cccc}\hfill \alpha _1i+\alpha _2j& & \alpha _1k+\alpha _2\mathrm{}\hfill & \\ & & & (modn)\hfill \\ \hfill \beta _1i+\beta _2j& & \beta _1k+\beta _2\mathrm{}\hfill & \end{array}$$

which for convenience we express in matrix notation

$$\left(\begin{array}{cc}\alpha _1& \alpha _2\\ \beta _1& \beta _2\end{array}\right)\left(\begin{array}{c}(ik)\hfill \\ (j\mathrm{})\hfill \end{array}\right)0(modn)$$

Now since by hypothesis $`(d_2,n)=1`$, so that in particular $`d_20`$ $`(modn)`$, we can multiply this relation by

$$d_2(\begin{array}{cc}\alpha _1& \alpha _2\\ \beta _1& \beta _2\end{array})^1=(\begin{array}{cc}\hfill \beta _2& \hfill \alpha _2\\ \hfill \beta _1& \hfill \alpha _1\end{array})$$

so that it becomes

$$d_2(\begin{array}{c}(ik)\\ (j\mathrm{})\end{array})0modn$$

Finally, since $`(d_2,n)=1`$ it follows that $`i=k`$ and $`j=\mathrm{}`$. Thus $`L^{(1)}`$ and $`L^{(2)}`$ are orthogonal pandiagonal LS’s.

The orthogonality criterion of the two LS’s in (2.6) is very convenient and as we shall see is extendable to higher dimensions.

Let us consider the question as to the number of distinct pairs $`[\alpha _1,\alpha _2]`$ there are for given $`n`$ for which $`L=<\alpha _1,\alpha _2>`$ is a pandiagonal LS. To simplify let us assume in the following that $`n`$ is a prime $`p`$. For a given prime $`p5`$, consider the $`p3`$, $`[\alpha _1,\alpha _2]`$, pairs

$$[\alpha _1,\alpha _2]=[1,2],[1,3],[1,4],\mathrm{}[1,p2]$$
$`(2.7)`$
Obviously each of these satisfies (2.4) and no other pair for which $`\alpha _1=1`$ does so.

Theorem (2.3): The pandiagonal latin squares associated with the $`[\alpha _1,\alpha _2]`$ pairs in (2.7) are:

pairwise mutually orthogonal

any other pandiagonal LS can be obtained from one associated with $`<1,\mathrm{}>`$, $`\mathrm{}=2,3,\mathrm{}p2`$, in (2.7) by a permutation of symbols.

Proof:

Consider the two pandiagonal LS’s $`<1,\mathrm{}>,<1,m>;\mathrm{},m=2,3,\mathrm{}p2,`$ $`(\mathrm{}m)`$. According to (2.5) the associated determinant $`d_2`$ is

$$d_2=\left|\begin{array}{cc}1\hfill & \mathrm{}\hfill \\ 1\hfill & m.\hfill \end{array}\right|=m\mathrm{}$$

and obviously satisfies (2.6) since we assume $`m\mathrm{}`$.

Let $`x,y`$ be any two unequal positive integers $`ϵ\overline{S}=\{1,2,\mathrm{}p1\}`$ that satisfy each of (2.4) so that $`<x,y>`$ is a pandiagonal LS. If $`x^1`$ is the inverse of $`x(modp)`$ then making use of (2.3) and Remark (2.2), we find that $`<1,x^1y>`$ is also pandiagonal and is obtained from $`<x,y>`$ by a permutation of its symbols. Finally since $`x^1y(modp)`$ must be one of $`(2,3,\mathrm{}p2)`$ in (2.7) it follows from Remarks 2.1 and 2.2 that $`<x,y>`$ can be obtained from $`<1,x^1y>`$ by a permutation of its symbols.

For example, for $`p=7`$, successive multiplication $`(modp)`$ by use of $`k=2,3,\mathrm{}6`$, leads to $`<1,2>`$ $``$ $`<2,4>,<3,6>,<4,1>,<5,3>,<6,5>`$. Similarly $`<1,3>`$ $``$ $`<2,6>`$, $`<3,2>`$, $`<4,5>`$, $`<5,1>`$, $`<6,4>`$; and $`<1,4>`$ $``$ $`<2,1>`$, $`<3,5>`$, $`<4,2>`$, $`<5,6>`$, $`<6,3>`$.

Remark 2.5: If we allow for permutation of symbols then all pandiagonal LS’s in (2.3) can be obtained, for given $`p`$ by use only of the $`[\alpha _1,\alpha _2]`$ pairs in (2.7).

Finally, as shown for pandiagonal magic squares, by Ball and Coxeter \[3, p. 203\] and by Martin Gardner  for a $`5\times 5`$ pandiagonal magic square, we find that pandiagonal LS’s have an analogous unusual property. If in an $`n\times n`$ pandiagonal LS we move the left hand column so that it becomes the right hand column (or vice versa) or similarly move the top row to the bottom, the resultant array is again a pandiagonal LS. This has the consequence that if we “tile” the plane with a given $`n\times n`$ pandiagonal LS we can outline any $`n\times n`$ square on this infinite pattern and obtain a pandiagonal LS.

The underlying result is contained in:

Theorem 2.4: If the left column of a pandiagonal latin square is moved so it becomes the right column, the resulting LS is a pandiagonal LS obtained from the original one by a permutation of symbols.

Proof: Let $`L=<\alpha _1,\alpha _2>`$ be the original pandiagonal LS and define a second one $`L^{}`$ with elements given by

$$L_{ij}^{}L_{ij}+\alpha _1(\alpha _1+1)i+\alpha _2j(modn)$$

It follows from (2.3) that $`L^{}`$ is simply $`L`$ with its elements permuted. Further we have

$$L_{oj}^{}L_{1j};L_{1j}^{}L_{2j};\mathrm{}L_{n1,j}^{}L_{oj};(modn);j=0,1,\mathrm{}n1$$

which shows that $`L^{}`$ is simply $`L`$ with its left column moved so it becomes the right column.

For example, if the left hand column of the LS in Figure IIc, which incidentally is simply $`<1,3>`$, (with $`i`$ increasing to the right and j increasing upwards) is moved to the right side, the original LS, but with permuted symbols, in IId results! Of considerable interest perhaps is that a form of this property, as will be shown below, has an analogue in higher dimensions.

B. Pandiagonal Latin Cubes

We define a pandiagonal latin cube as an $`n\times n\times n`$ cube each of whose $`n^3`$ lattice points contains one of the members of $`S=\{0,1,2,\mathrm{}n1\}`$ and in a way so that each of its $`3n+6`$ constituent squares is a pandiagonal latin square. Recall in this connection that an $`n\times n\times n`$ cube has $`3n`$ squares parallel to a cube face plus 6 “diagonal” squares which contain its 4 body diagonals. Note that, as here defined, in a pandiagonal latin cube each row, column, file, diagonal and broken diagonal of each of its squares, i.e., its ROWS, contain each element of $`S`$ once and only once.

In an obvious generalization of (2.3) to three dimensions, we define an $`n\times n\times n`$ array C by the formula

$$C_{ijk}\alpha _1i+\alpha _2j+\alpha _3k(modn)$$
$`(2.8)`$
where each element of $`C_{ijk}`$ $`ϵS,\{0,1,\mathrm{}n1\}`$. Here i,j,k are integer variables each running over $`\{0,1,\mathrm{}n1\}`$, and $`\alpha _1,\alpha _2,\alpha _3`$ are elements of $`\overline{S}=\{1,2,\mathrm{}n1\}`$. We shall use the notation $`C=<\alpha _1,\alpha _2,\alpha _3>`$ when we wish to focus on the dependence of $`C`$ on $`\alpha _1,\alpha _2,\alpha _3`$. We now establish the following:

Theorem (2.5): The cube defined in (2.8) is a pandiagonal latin cube—in that each of its constituent 3n + 6 squares is pandiagonal—provided $`\alpha _1,\alpha _2,\alpha _3`$ satisfy the constraints:

$$(\alpha _{\mathrm{}},n)=1;\mathrm{}=1,2,3$$
$`(2.9a)`$
$$(\alpha _{\mathrm{}}+\alpha _{\mathrm{}^{}},n)=1;\mathrm{},\mathrm{}^{}=1,2,3;\mathrm{}\mathrm{}^{}$$
$`(2.9b)`$
$$(\alpha _{\mathrm{}}\alpha _{\mathrm{}^{}},n)=1;\mathrm{},\mathrm{}^{}=1,2,3;\mathrm{}\mathrm{}^{}$$
$`(2.9c)`$
$$(A,n)=1;A=\alpha _1+\alpha _2+\alpha _3$$
$`(2.9d)`$
$$(A2\alpha _{\mathrm{}},n)=1;\mathrm{}=1,2,3$$
$`(2.9e)`$
Proof: Consider first the $`3n`$ squares parallel to the faces of the cube. For fixed $`k`$, say, 0 $`kn1`$, consider the square of C<sub>ijk</sub> that is parallel to the $`ij`$ plane as $`i`$ and $`j`$ vary over the elements of S. Then C<sub>ijk</sub> has essentially the same structure as does $`L_{ij}`$ in (2.3) since the constant $`\alpha _3k`$ is of no consequence. Making use of the restrictions on the $`\alpha `$’s in (2.9 a,b,c) and comparing with those in (2.4 a, b, c) we conclude that for any fixed $`k`$, the square $`C_{ijk}`$ is a pandiagonal latin square. Repeating this argument for fixed $`i`$, with $`j`$ and $`k`$ variable and for fixed $`j`$ with $`i`$ and $`k`$ variable we conclude that all 3n squares in the cube defined by (2.8) that are parallel to a cube face are pandiagonal latin squares.

With regard to the six squares not parallel to a cube face we proceed as follows. For that diagonal square, with vertices at $`(i,j,k)=(0,0,0),(0,0,n1),(n1,n1,n1),(n1,n1,0)`$, we have $`i=j`$ so that (2.8) becomes

$$C_{iik}(\alpha _1+\alpha _2)i+\alpha _3k(modn).$$

But this is again of the form in (2.3) if we make the replacements in (2.4a,b,c) $`\alpha _1\alpha _1+\alpha _2;\alpha _2\alpha _3`$. The first of (2.4) is satisfied because of (2.9a,b) while (2.4b) and (2.4c) become respectively $`(\alpha _1+\alpha _2+\alpha _3,n)=1`$ and $`(\alpha _1+\alpha _2\alpha _3,n)=1`$ and these are the same constraints as in (2.9d) and (2.9e), respectively. Thus $`C_{iik}`$, the given square, is a pandiagonal latin square.

Similarly for the diagonal square perpendicular to $`C_{iik}`$ whose vertices have the coordinates $`(i,j,k)=(n1,0,0),(0,n1,0),(0,n1,n1),(n1,0,n1)`$, we have $`i+j=n1`$ and

$$C_{i,(n1i),k}(\alpha _1\alpha _2)i+\alpha _3k+\alpha _2(n1)(modn)$$

which on comparison with (2.3) and (2,4) with the replacements $`\alpha _1\alpha _1\alpha _2;\alpha _2\alpha _3`$ leads to the conditions $`(\alpha _1\alpha _2+\alpha _3,n)=1`$, and $`(\alpha _1\alpha _2\alpha _3,n)=1`$ The first of these is the same as (2.9e) with $`\mathrm{}`$ = 2 and the second is the same as (2.9e) with $`\mathrm{}`$ = 1 if we make use of the fact that (-x, y) = 1 is equivalent to (x,y) = 1.

A similar argument shows that for the constraints in (2.9) the remaining four “diagonal squares” of the cube are also pandiagonal latin squares. The theorem is thus established.

Just as for the 2-D case, we require a 3-D analogue of orthogonality of latin cubes. For our purposes, we shall say that three latin cubes are orthogonal if when they are superposed none of the $`n^3`$ ordered triplets of elements that result occurs more than once. There are other definitions of orthogonal latin cubes, \[9, ch. 3\], but for purposes of producing magic cubes this definition is essential.

We now establish:

Theorem (2.6): Consider the 3 pandiagonal latin cubes

$$C_{ijk}^{(q)}\alpha _{1q}i+\alpha _{2q}j+\alpha _{3q}k(modn);q=1,2,3$$
$`(2.10)`$
where $`\alpha _{pq}`$ $`(p,q=1,2,3)`$ are elements of $`\overline{S}`$ and for each value for $`q`$ satisfy the conditions in (2.9) and let $`d_3`$ be the determinant $`d_3=|\alpha _{pq}|`$. Then if $`d_3`$ is relatively prime to $`n`$, that is

$$(d_3,n)=1$$
$`(2.11)`$
then the three cubes $`C^{(q)}`$, $`q=1,2,3`$ are orthogonal.

Proof: Suppose on the contrary that for a given $`(i,j,k)`$ there existed an integer triplet $`(u,v,w)`$ each $`ϵS`$ not equal to $`(i,j,k)`$ for which $`C_{ijk}^qC_{uvw}^q`$ $`(modn)`$ for each $`q`$. Then, as for the 2-D case, we could express this in matrix notation

$$\left(\begin{array}{ccc}\hfill \alpha _{11}& \hfill \alpha _{21}& \hfill \alpha _{31}\\ \hfill \alpha _{12}& \hfill \alpha _{22}& \hfill \alpha _{32}\\ \hfill \alpha _{13}& \hfill \alpha _{23}& \hfill \alpha _{33}\end{array}\right)\left(\begin{array}{c}iu\hfill \\ jv\hfill \\ kw\hfill \end{array}\right)0modn$$
$`(2.12)`$
Since by (2.11), the determinant $`d_3`$ is relatively prime to $`n`$, so that in particular $`d_30`$ $`modn`$, we may multiply both sides of (2.12) by $`d_3`$ times the inverse of the matrix on the left. The result is

$$d_3\left(\begin{array}{c}iu\hfill \\ jv\hfill \\ kw\hfill \end{array}\right)0modn.$$

Finally since $`(d_3,n)=1`$ it follows that $`i=u;j=v;k=w`$ and the theorem is proved.

An empirical study of the constraints in (2.9) shows that no triplets $`[\alpha _1,\alpha _2,\alpha _3]`$ satisfying (2.9) exists for $`n=5,7,9`$. (The latter clearly since $`n`$ would be divisible by 3.) Such an analysis is most easily carried out by recognizing that—as in the 2-D case—without loss of generality we can take $`\alpha _1=1`$ and consider simply the triplets $`[1,\alpha _2,\alpha _3]`$. For $`n=7`$, it is easily confirmed that no values, such $`[1,2,3],[1,2,4],[1,2,5]`$ would satisfy all of (2.9). For $`n=11`$ however we find, among others, the possibilities $`[1,2,4],[1,2,5],[1,2,6],[1,2,7],[1,5,8]`$, $`[1,6,8]`$ as well as these with $`\alpha _2`$ and $`\alpha _3`$ interchanged. As an example of orthogonal pandiagonal cubes we note that for appropriate integers, $`\mathrm{},m,p`$ for the three LS’s $`<1,2,\mathrm{}>`$, $`<1,m,2>`$ and $`<1,2,p>`$, for $`n=11`$, $`d_3`$ is given by

$$d_3=\left|\begin{array}{ccc}1\hfill & 2\hfill & \mathrm{}\hfill \\ 1\hfill & m\hfill & 2\hfill \\ 1\hfill & 2\hfill & p\hfill \end{array}\right|=(m2)(p\mathrm{}),$$

and will for values of $`\mathrm{},m,p`$ with $`\mathrm{}p`$ and $`4\mathrm{},m,p9`$ lead to $`(d_3,n)=1`$. Figure III shows two planar sections through the cube $`<1,2,7>`$ for $`n=11:`$ (a) corresponds to the square $`k=2`$ in (2.8) and (b) to the diagonal square $`i=j`$.

Note, however, that for the triplet $`<1,2,\mathrm{}>,<1,2,m>,<1,2,p>`$, $`d_3=0`$ so that these three do not constitute orthogonal cubes as we have defined them even though each cube is itself pandiagonal.

As for the analogous 2-D case (Theorem 2.4) we can easily establish the fact that if we move, say, a face of a pandiagonal latin cube to its opposite

$$\begin{array}{c}\begin{array}{c}\left[\begin{array}{ccccccccccc}\hfill 3& \hfill 5& \hfill 7& \hfill 9& \hfill 0& \hfill 2& \hfill 4& \hfill 6& \hfill 8& \hfill 10& \hfill 1\\ \hfill 4& \hfill 6& \hfill 8& \hfill 10& \hfill 1& \hfill 3& \hfill 5& \hfill 7& \hfill 9& \hfill 0& \hfill 2\\ \hfill 5& \hfill 7& \hfill 9& \hfill 0& \hfill 2& \hfill 4& \hfill 6& \hfill 8& \hfill 10& \hfill 1& \hfill 3\\ \hfill 6& \hfill 8& \hfill 10& \hfill 1& \hfill 3& \hfill 5& \hfill 7& \hfill 9& \hfill 0& \hfill 2& \hfill 4\\ \hfill 7& \hfill 9& \hfill 0& \hfill 2& \hfill 4& \hfill 6& \hfill 8& \hfill 10& \hfill 1& \hfill 3& \hfill 5\\ \hfill 8& \hfill 10& \hfill 1& \hfill 3& \hfill 5& \hfill 7& \hfill 9& \hfill 0& \hfill 2& \hfill 4& \hfill 6\\ \hfill 9& \hfill 0& \hfill 2& \hfill 4& \hfill 6& \hfill 8& \hfill 10& \hfill 1& \hfill 3& \hfill 5& \hfill 7\\ \hfill 10& \hfill 1& \hfill 3& \hfill 5& \hfill 7& \hfill 9& \hfill 0& \hfill 2& \hfill 4& \hfill 6& \hfill 8\\ \hfill 0& \hfill 2& \hfill 4& \hfill 6& \hfill 8& \hfill 10& \hfill 1& \hfill 3& \hfill 5& \hfill 7& \hfill 9\\ \hfill 1& \hfill 3& \hfill 5& \hfill 7& \hfill 9& \hfill 0& \hfill 2& \hfill 4& \hfill 6& \hfill 8& \hfill 10\\ \hfill 2& \hfill 4& \hfill 6& \hfill 8& \hfill 10& \hfill 1& \hfill 3& \hfill 5& \hfill 7& \hfill 9& \hfill 0\end{array}\right]\hfill \\ (a)\hfill \end{array}\hfill \\ \begin{array}{c}\left[\begin{array}{ccccccccccc}\hfill 0& \hfill 7& \hfill 3& \hfill 10& \hfill 6& \hfill 2& \hfill 9& \hfill 5& \hfill 1& \hfill 8& \hfill 4\\ \hfill 3& \hfill 10& \hfill 6& \hfill 2& \hfill 9& \hfill 5& \hfill 1& \hfill 8& \hfill 4& \hfill 0& \hfill 7\\ \hfill 6& \hfill 2& \hfill 9& \hfill 5& \hfill 1& \hfill 8& \hfill 4& \hfill 0& \hfill 7& \hfill 3& \hfill 10\\ \hfill 9& \hfill 5& \hfill 1& \hfill 8& \hfill 4& \hfill 0& \hfill 7& \hfill 3& \hfill 10& \hfill 6& \hfill 2\\ \hfill 1& \hfill 8& \hfill 4& \hfill 0& \hfill 7& \hfill 3& \hfill 10& \hfill 6& \hfill 2& \hfill 9& \hfill 5\\ \hfill 4& \hfill 0& \hfill 7& \hfill 3& \hfill 10& \hfill 6& \hfill 2& \hfill 9& \hfill 5& \hfill 1& \hfill 8\\ \hfill 7& \hfill 3& \hfill 10& \hfill 6& \hfill 2& \hfill 9& \hfill 5& \hfill 1& \hfill 8& \hfill 4& \hfill 0\\ \hfill 10& \hfill 6& \hfill 2& \hfill 9& \hfill 5& \hfill 1& \hfill 8& \hfill 4& \hfill 0& \hfill 7& \hfill 3\\ \hfill 2& \hfill 9& \hfill 5& \hfill 1& \hfill 8& \hfill 4& \hfill 0& \hfill 7& \hfill 3& \hfill 10& \hfill 6\\ \hfill 5& \hfill 1& \hfill 8& \hfill 4& \hfill 0& \hfill 7& \hfill 3& \hfill 10& \hfill 6& \hfill 2& \hfill 9\\ \hfill 8& \hfill 4& \hfill 0& \hfill 7& \hfill 3& \hfill 10& \hfill 6& \hfill 2& \hfill 9& \hfill 5& \hfill 1\end{array}\right]\hfill \\ (b)\hfill \end{array}\hfill \end{array}$$

Figure III

side, the resultant cube remains a pandiagonal latin cube. To see this, consider $`C=<\alpha _1,\alpha _2,\alpha _3>`$ in (2.8) and define $`C^{}`$ with elements given by $`C_{ijk}^{}=C_{ijk}+\alpha _1`$ so that

$$C_{ijk}^{}\alpha _1(i+1)+\alpha _2j+\alpha _3k(modn).$$

Obviously,

$$C_{ojk}^{}=C_{1jk};C_{1jk}^{}C_{2jk};\mathrm{}C_{n1,jk}^{}C_{o,jk}(modn);0j,kn1$$

Thus $`C^{}`$ is obtained from $`C`$ by transporting its $`i=0`$ face to $`i=n1`$. This implies that if we “tile” all of 3-D space with a given pandiagonal latin cube any $`n\times n\times n`$ cube selected out of this infinite array will be a pandiagonal latin cube but with its elements permuted in some way.

C. Four dimensional pandiagonal latin cubes

A 4-D pandiagonal latin cube is an arrangement of the integers $`\{0,1,2,\mathrm{}n1\}`$ among the cells of an $`n\times n\times n\times n`$ cube in a way so that each of its $`4n+12`$ constituent 3-D cubes is a pandiagonal latin cube.

In an obvious generalization of (2.8), we define an $`n\times n\times n\times n`$ array of integers by the formula

$$H_{ijk\mathrm{}}\alpha _1i+\alpha _2j+\alpha _3k+\alpha _4\mathrm{}(modn)$$
$`(2.13)`$
where each of $`H_{ijk\mathrm{}}`$ is $`ϵ`$ $`S`$ = $`\{0,1,2,\mathrm{},n1\}`$ and where $`i,j,k`$, $`\mathrm{}`$ are integer variables each with the range $`0,1,2,\mathrm{},n1`$. The four quantities $`\alpha _1,\alpha _2,\alpha _3,\alpha _4`$ are integer parameters $`ϵ\overline{S}=\{1,2,\mathrm{}n1\}`$. As before we use the notation $`<\alpha _1,\alpha _2,\alpha _3,\alpha _4>`$ when the dependence of $`H`$ on the $`\alpha `$’s is of interest. We now establish the following:

Theorem (2.7): The hypercube defined in (2.13) is a 4-D pandiagonal latin cube provided the four integer parameters $`\alpha _1,\alpha _2,\alpha _3,\alpha _4`$ satisfy the constraints

$$\begin{array}{ccc}(\alpha _m,n)=1\hfill & ;m=1,2,3,4\hfill & \hfill (2.14a)\\ (\alpha _m+\alpha _m^{},n)=1\hfill & ;m,m^{}=1,2,3,4;mm^{}\hfill & \hfill (2.14b)\\ (\alpha _m\alpha _m^{},n)=1\hfill & ;m,m^{}=1,2,3,4;mm^{}\hfill & \hfill (2.14c)\\ (B\alpha _m,n)=1\hfill & ;B=\alpha _1+\alpha _2+\alpha _3+\alpha _4;m=1,2,3,4\hfill & \hfill (2.14d)\\ (B,n)=1\hfill & ;\hfill & \hfill (2.14e)\\ (B2\alpha _m,n)=1\hfill & ;m=1,2,3,4\hfill & \hfill (2.14f)\\ (B\alpha _m^{}2\alpha _m,n)=1\hfill & ;m,m^{}=1,2,3,4;mm^{}\hfill & \hfill (2.14g)\\ (B2\alpha _m^{}2\alpha _m,n)=1\hfill & ;m,m^{}=1,2,3,4;mm^{}\hfill & \hfill (2.14h)\end{array}$$

Proof: Consider first the 3-D cube that results for fixed $`\mathrm{}`$ from (2.13) as $`i,j,k`$ vary over S. It will be a pandiagonal latin cube provided the constraints in (2.9a)-(2.9e) are satisfied. Now (2.9a)-(2.9e) are contained within (by appropriate choice of the subscripts) (2.14a), (2.14b), (2.14c), (2.14e) and (2.14g), respectively. Thus since the added constant $`\alpha _4\mathrm{}`$ plays no role the cube is a pandiagonal latin cube. Similarly, the cubes that result for fixed $`i`$, as $`j,k,\mathrm{}`$ vary, and for fixed $`j`$ as $`i,k,\mathrm{}`$ vary and for fixed $`k`$ as $`i,j,\mathrm{}`$ vary over S are all pandiagonal latin cubes. There are altogether $`4n`$ pandiagonal 3-D latin cubes of this type.

Similarly the “diagonal” 3-D cube that results for $`i=j`$, as $`i,k,\mathrm{}`$ vary over S has the form

$$H_{iik\mathrm{}}(\alpha _1+\alpha _2)i+\alpha _3k+\alpha _4\mathrm{}(modn)$$

On comparison with (2.8) and (2.9a) - (2.9e) we see that this is also a pandiagonal latin cube since we assumed (2.14a), (2.14b) (2.14c), (2.14d), (2.14e), (2.14g) to be satisfied . And similarly for the other five pairs: $`i=k,i=\mathrm{}`$, $`j=k,j=\mathrm{}`$, $`k=\mathrm{}`$.

Finally for the “diagonal” latin 3-D cube that results for $`i+j=n1`$, as $`i,k`$ and $`\mathrm{}`$ vary over S

$$H_{i,n1i,k,\mathrm{}}(\alpha _1\alpha _2)i+\alpha _3k+\alpha _4\mathrm{}\alpha _2(n1)(modn)$$

which on comparison with (2.9a) - (2.9c) is a pandiagonal latin cube by virtue of (2.14a), (2.14b), (2.14c), (2.14f), (2.14g), (2.14h). And similarly for the remaining five cubes $`i+k=n1,i+\mathrm{}=n1,j+k=n1,j+\mathrm{}=n1`$ and $`k+\mathrm{}=n1`$. This proves the theorem.

Turning to the question of orthogonality, we define four, 4-D pandiagonal cubes to be orthogonal, if when they are superposed no two of the $`n^4`$ ordered quartets of elements that result are the same. As in the lower dimensional cases there are other definitions of this orthogonality , but for our purposes this one is essential.

By analogy to Theorems (2.2) and (2.6) we have:

Theorem (2.8): Consider the four, 4-D pandiagonal latin cubes

$$H_{ijk\mathrm{}}^{(q)}\alpha _{1q}i+\alpha _{2q}j+d_{3q}k+\alpha _{4q}\mathrm{}(modn);q=1,2,3,4$$
$`(2.15)`$
where $`\alpha _{pq}(p,q=1,2,3,4)`$ are elements of $`\overline{S}`$ which for each value of $`q`$ satisfy the conditions in (2.14) and the variables $`i,j,k,\mathrm{}`$ range over $`S=\{0,1,2,\mathrm{}n1\}`$. Further, let $`d_4=|\alpha _{pq}|`$ be the determinant of the $`\alpha `$’s. Then if $`d_4`$ is relatively prime to $`n`$

$$(d_4,n)=1$$
$`(2.16)`$
then the four 4-D cubes in (2.15) are orthogonal.

Proof: Suppose to the contrary there existed a quartet of integers $`(u,v,w,x)`$ each $`ϵ`$ $`S`$, not equal to any given quartet $`(i,j,k,\mathrm{})`$ for which $`H_{ijk\mathrm{}}^{(q)}H_{uvwx}^{(q)}`$ $`(modn)`$ for each $`q`$. Then just as in deriving (2.12) we would find

$$\left(\begin{array}{cccc}\hfill \alpha _{11}& \hfill \alpha _{21}& \hfill \alpha _{31}& \hfill \alpha _{41}\\ \hfill \alpha _{12}& \hfill \alpha _{22}& \hfill \alpha _{32}& \hfill \alpha _{42}\\ \hfill \alpha _{13}& \hfill \alpha _{23}& \hfill \alpha _{33}& \hfill \alpha _{43}\\ \hfill \alpha _{14}& \hfill \alpha _{24}& \hfill \alpha _{34}& \hfill \alpha _{44}\end{array}\right)\left(\begin{array}{c}iu\hfill \\ jv\hfill \\ kw\hfill \\ \mathrm{}x\hfill \end{array}\right)0modn$$
$`(2.17)`$
and conclude following essentially the same steps as before that $`i,j,k,\mathrm{}`$ must be equal to $`u,v,w,x`$ respectively. Thus concluding the proof of orthogonality.

It is easily confirmed by enumerating the various possibilities that only for n $``$ 17 is it possible to find integers $`[\alpha _1,\alpha _2,\alpha _3,\alpha _4]`$ among $`\{1,2,\mathrm{}n1\}`$ for which (2.14) can be satisfied. In particular for $`n=17`$, possible 4-D latin hypercubes are given by $`<1,2,4,8>`$, $`<1,2,4,9>`$, $`<1,2,13,8>`$, $`<1,2,13,9>`$ as is readily confirmed. Thus a possible form for the determinant of the matrix in (2.17) is

$$d_4=\left|\begin{array}{ccccccccccc}& 1\hfill & 2\hfill & 4& 8\hfill & & & & & & \\ & 1\hfill & 2\hfill & 4& 9\hfill & & & & & & \\ & 1\hfill & 2\hfill & 8& 4\hfill & & & & & & \\ & 1\hfill & 4\hfill & 9& 2\hfill & & & & & & \end{array}\right|=8$$

and since (-8, 17) = 1, the four 4-D pandiagonal latin cubes $`<1,2,4,8>,<1,2,4,9>,<1,2,8,4>,<1,4,9,2>`$ constitute an orthogonal set of such hypercubes, for $`n=17`$.

Figure IV, shows a planar section through the hypercube $`<1,2,4,9>`$ corresponding to $`i=2,j+k=16`$ in (2.13).

$$\begin{array}{c}\left[\begin{array}{ccccccccccccccccc}\hfill 0& \hfill 9& \hfill 1& \hfill 10& \hfill 2& \hfill 11& \hfill 3& \hfill 12& \hfill 4& \hfill 13& \hfill 5& \hfill 14& \hfill 6& \hfill 15& \hfill 7& \hfill 16& \hfill 8\\ \hfill 2& \hfill 11& \hfill 3& \hfill 12& \hfill 4& \hfill 13& \hfill 5& \hfill 14& \hfill 6& \hfill 15& \hfill 7& \hfill 16& \hfill 8& \hfill 0& \hfill 9& \hfill 1& \hfill 10\\ \hfill 4& \hfill 13& \hfill 5& \hfill 14& \hfill 6& \hfill 15& \hfill 7& \hfill 16& \hfill 8& \hfill 0& \hfill 9& \hfill 1& \hfill 10& \hfill 2& \hfill 11& \hfill 3& \hfill 12\\ \hfill 6& \hfill 15& \hfill 7& \hfill 16& \hfill 8& \hfill 0& \hfill 9& \hfill 1& \hfill 10& \hfill 2& \hfill 11& \hfill 3& \hfill 12& \hfill 4& \hfill 13& \hfill 5& \hfill 14\\ \hfill 8& \hfill 0& \hfill 9& \hfill 1& \hfill 10& \hfill 2& \hfill 11& \hfill 3& \hfill 12& \hfill 4& \hfill 13& \hfill 5& \hfill 14& \hfill 6& \hfill 15& \hfill 7& \hfill 16\\ \hfill 10& \hfill 2& \hfill 11& \hfill 3& \hfill 12& \hfill 4& \hfill 13& \hfill 5& \hfill 14& \hfill 6& \hfill 15& \hfill 7& \hfill 16& \hfill 8& \hfill 0& \hfill 9& \hfill 1\\ \hfill 12& \hfill 4& \hfill 13& \hfill 5& \hfill 14& \hfill 6& \hfill 15& \hfill 7& \hfill 16& \hfill 8& \hfill 0& \hfill 9& \hfill 1& \hfill 10& \hfill 2& \hfill 11& \hfill 3\\ \hfill 14& \hfill 6& \hfill 15& \hfill 7& \hfill 16& \hfill 8& \hfill 0& \hfill 9& \hfill 1& \hfill 10& \hfill 2& \hfill 11& \hfill 3& \hfill 12& \hfill 4& \hfill 13& \hfill 5\\ \hfill 16& \hfill 8& \hfill 0& \hfill 9& \hfill 1& \hfill 10& \hfill 2& \hfill 11& \hfill 3& \hfill 12& \hfill 4& \hfill 13& \hfill 5& \hfill 14& \hfill 6& \hfill 15& \hfill 7\\ \hfill 1& \hfill 10& \hfill 2& \hfill 11& \hfill 3& \hfill 12& \hfill 4& \hfill 13& \hfill 5& \hfill 14& \hfill 6& \hfill 15& \hfill 7& \hfill 16& \hfill 8& \hfill 0& \hfill 9\\ \hfill 3& \hfill 12& \hfill 4& \hfill 13& \hfill 5& \hfill 14& \hfill 6& \hfill 15& \hfill 7& \hfill 16& \hfill 8& \hfill 0& \hfill 9& \hfill 1& \hfill 10& \hfill 2& \hfill 11\\ \hfill 5& \hfill 14& \hfill 6& \hfill 15& \hfill 7& \hfill 16& \hfill 8& \hfill 0& \hfill 9& \hfill 1& \hfill 10& \hfill 2& \hfill 11& \hfill 3& \hfill 12& \hfill 4& \hfill 13\\ \hfill 7& \hfill 16& \hfill 8& \hfill 0& \hfill 9& \hfill 1& \hfill 10& \hfill 2& \hfill 11& \hfill 3& \hfill 12& \hfill 4& \hfill 13& \hfill 5& \hfill 14& \hfill 6& \hfill 15\\ \hfill 9& \hfill 1& \hfill 10& \hfill 2& \hfill 11& \hfill 3& \hfill 12& \hfill 4& \hfill 13& \hfill 5& \hfill 14& \hfill 6& \hfill 15& \hfill 7& \hfill 16& \hfill 8& \hfill 0\\ \hfill 11& \hfill 3& \hfill 12& \hfill 4& \hfill 13& \hfill 5& \hfill 14& \hfill 6& \hfill 15& \hfill 7& \hfill 16& \hfill 8& \hfill 0& \hfill 9& \hfill 1& \hfill 10& \hfill 2\\ \hfill 13& \hfill 5& \hfill 14& \hfill 6& \hfill 15& \hfill 7& \hfill 16& \hfill 8& \hfill 0& \hfill 9& \hfill 1& \hfill 10& \hfill 2& \hfill 11& \hfill 3& \hfill 12& \hfill 4\\ \hfill 15& \hfill 7& \hfill 16& \hfill 8& \hfill 0& \hfill 9& \hfill 1& \hfill 10& \hfill 2& \hfill 11& \hfill 3& \hfill 12& \hfill 4& \hfill 13& \hfill 5& \hfill 14& \hfill 6\end{array}\right]\hfill \end{array}$$

Figure IV

Just as for the 2 and 3 dimensional cases, if $`H=<\alpha _1,\alpha _2,\alpha _3,\alpha _4>`$ is a 4-D pandiagonal latin cube, then so is $`H^{}`$ with elements given by

$$H_{ijk\mathrm{}}^{}=H_{ijk\mathrm{}}+\alpha _1(\alpha _1+1)i+\alpha _2j+\alpha _3k+\alpha _4\mathrm{}(modn).$$

Again we can imagine “tiling” all of 4-D space with $`H`$ and be assured that any $`n\times n\times n\times n`$ subcube in this space will be a 4-D pandiagonal latin cube, with its elements a permutation of the original elements of $`H`$.

III. Magic, Pandiagonal Squares and Cubes in Two and Three Dimensions

With the results of the preceding section available, it is now straightforward , to generate magic pandiagonal cubes in three and four dimensions. To set the stage we first illustrate the matter in two dimensions.

A. Magic Pandiagonal Squares

Recall that a magic pandiagonal square of order $`n`$ is an arrangement, without repeats, of the integers $`(0,1,2,\mathrm{}n^21)`$ on the $`n^2`$ lattice points of an $`n\times n`$ array so that the sum of the elements in each of the $`n`$rows, $`n`$columns, $`n`$diagonals (including the $`n2`$ broken diagonals) has the same value. (As above let us use the generic ROW to represent any one of these rows, columns, diagonals, etc. ) This common sum of the ROWS has been given in (1.1).

Theorem 3.1 \[9, p. 178\] Let $`L^{(1)}=<\alpha _1,\alpha _2>`$, $`L^{(2)}=<\beta _1,\beta _2>`$ be two orthogonal pandiagonal latin squares that separately satisfy the conditions of theorems (2.1) and (2.2) and define an $`n\times n`$ array $`M^{(2)}`$ with elements $`M_{ij}^{(2)}`$ by

$$M^{(2)}=nL^{(1)}+L^{(2)}$$
$`(3.1)`$
Then $`M^{(2)}`$ is a magic pandiagonal square.

Proof: Since the elements of $`L^{(1)}`$ and $`L^{(2)}`$ range over $`S=\{0,1,2,\mathrm{}n1\}`$ it follows from (3.1) that each of the elements $`M_{ij}^{(2)}`$ must be one of the integers $`0,1,2,\mathrm{}n^21`$. Further, since the sum of the elements in each ROW of the pandiagonal $`L^{(1)}`$ and $`L^{(2)}`$, is $`0+1+2+\mathrm{}+n1`$ $`=n(n1)/2`$ it follows that the sum of each of the $`4n`$ ROWS of $`M^{(2)}`$ is

$$n[n(n1)/2]+n(n1)/2=n(n^21)/2=\sigma _2$$

with $`\sigma _2`$ defined in (1.1). Finally since $`L^{(1)}`$ and $`L^{(2)}`$ are orthogonal, no two elements of $`M^{(2)}`$ can be the same and the theorem is established.

Since $`L^{(1)}`$ and $`L^{(2)}`$ are pandiagonal LS’s for which we know $`n5`$, only for $`n5`$ can $`M^{(2)}`$ in (3.1) be a magic pandiagonal square. Further if we tile the plane with $`M^{(2)}`$ any $`n\times n`$ subsquare will also be a pandiagonal magic square since $`L^{(1)}`$ and $`L^{(2)}`$ have this same tiling property . (See the discussion at the end of IIA.)

Since the elements of the orthogonal pandiagonal latin squares $`L^{(1)}`$ and $`L^{(2)}`$ may be permuted among themselves without changing their essential properties, it follows that from a given pandiagonal LS we may generate $`n`$! different versions. This leads to the number $`N_2`$ of pandiagonal magic squares obtainable by our method to be

$$N_2=\mathrm{}(n)(n!)^2$$
$`(3.2)`$
where $`\mathrm{}(n)`$ is a low order polynomial in $`n`$. For $`n=5`$, since the only independent LS’s are $`<1,2>`$ and $`<1,3>`$ (Theorem 2.3), and since we may interchange their roles in (3.1) $`\mathrm{}(5)=2`$ and we obtain consistent with the results of Rosser and Walker

$$N_2=2880$$

Clearly because of the $`(n!)^2`$ factor $`N_2`$ rises very rapidly with $`n`$. For $`n=11`$ for example, $`N_2`$ is $``$ $`9.0\times 10^{16}`$.

B. Pandiagonal Magic Cubes

A magic pandiagonal $`n\times n\times n`$ cube is an arrangement, without repeats, of the integers $`0,1,2,\mathrm{}n^31`$ onto the $`n^3`$ lattice sites of a cube so that the sums along each ROW (i.e. along each of the $`n^2`$-rows, $`n^2`$-columns, $`n^2`$-files and the $`n(3n+6)`$ diagonals (including the broken diagonals) in each of its $`(3n+6)`$ squares are the same. The common sum along the ROWS is given by $`\sigma _3`$ in (1.2).

Theorem (3.2): Let $`C^{(q)}=<\alpha _{1q},\alpha _{2q},\alpha _{3q}>`$, $`q=1,2,3`$, be three orthogonal pandiagonal $`n\times n\times n`$ latin cubes that satisfy the conditions of theorems (2.5),(2.6) and define an $`n\times n\times n`$ array $`M^{(3)}`$ with elements $`M_{ijk}^{(3)}`$ by

$$M^{(3)}=n^2C^{(1)}+nC^{(2)}+C^{(3)}.$$
$`(3.3)`$
Then $`M^{(3)}`$ is a magic pandiagonal cube.

Proof: Since the elements of $`C^{(1)},C^{(2)}`$ and $`C^{(3)}`$ range over $`S=\{0,1,2,\mathrm{}n1\}`$ it follows from (3.3) that each element of $`M^{(3)}`$ must be one of the integers $`\{0,1,2,\mathrm{}n^31\}`$. Consider now any ROW of $`M^{(3)}`$. The sum of the elements of this ROW is, according to (3.3), given by

$$n^2[n(n1)/2]+n[n(n1)/2]+n(n1)/2=n(n^31)/2$$

and yields the value $`\sigma _3`$ in (1.2). Thus the sum of the elements in each ROW of the cube is the same and since $`C^{(1)}`$, $`C^{(2)}`$ and $`C^{(3)}`$ are orthogonal, it follows from (3.3) that no two elements of $`M^{(3)}`$ can be the same. The theorem is thus established.

According to the discussion in IIB, since the C’s are pandiagonal latin cubes, for which we know $`n11`$, $`M^{(3)}`$ will exist only for $`n11`$.

As for the two dimensional case, if we allow for the interchange in $`M^{(3)}`$, of $`\alpha _{1q},\alpha _{2q},\alpha _{3q}`$ $`(q=1,2,3)`$ with each other, and of permuting the symbols in each of $`C^{(1)}`$, $`C^{(2)}`$ and $`C^{(3)}`$ independently we can conclude that $`N_3`$, the number of cubes $`M^{(3)}`$ in (3.3) is given by

$$N_3=\mathrm{}_3(n)(n!)^3$$
$`(3.4)`$
where $`\mathrm{}_3`$ is an appropriate polynomial in $`n`$. For $`n=11,13`$ and 17, $`(n!)^3`$ assumes the approximate values $`6.4\times 10^{22}`$, $`2.4\times 10^{29}`$ and $`4.5\times 10^{43}`$ respectively. Because of this rapid rise of $`N_3`$ with $`n`$, we can anticipate that the factor $`\mathrm{}_3(n)`$ will not affect this variation, qualitatively.

As for the analogous 2-D case, it follows that by virtue of the “tiling” properties of pandiagonal latin cubes (see IIB), we can also tile 3-D space with any magic pandiagonal cube, and be assured that any $`n\times n\times n`$ cube selected out of this infinite array will also be a magic pandiagonal cube. It will differ from the original cube in that the elements of its underlying latin cubes will have been permuted.

Four Dimensional Pandiagonal Magic Cubes

By analogy to the above, we define a magic four dimensional pandiagonal cube, as an arrangement, without repeats, of the integers $`0,1,2,\mathrm{}n^41`$ among the $`n^4`$ lattice sites of an $`n\times n\times n\times n`$ cube so that the sum of the elements in each ROW of the 4-D cube, has the same value $`\sigma _4`$ as given in (1.3).

Theorem (3.3): Let $`Q^{(q)}=<\alpha _{1q},\alpha _{2q},\alpha _{3q},\alpha _{4q}>`$, $`q=1,2,3,4`$ be four orthogonal pandiagonal $`n\times n\times n\times n`$ latin cubes that satisfy the conditions of theorems (2.6) and (2.7) and define an $`n\times n\times n\times n`$ array $`M^{(4)}`$ with elements $`M_{ijk\mathrm{}}^{(4)}`$ by

$$M^{(4)}=n^3Q^{(1)}+n^2Q^{(2)}+nQ^{(3)}+Q^{(4)}$$
$`(3.5)`$
Then $`M^{(4)}`$ is a magic, pandiagonal four dimensional cube.

Proof: Since the elements of each of $`Q^{(1)},Q^{(2)},Q^{(3)}`$ and $`Q^{(4)}`$ range over $`S=\{0,1,2,\mathrm{}n1\}`$, it follows from (3.5) that the elements of $`M_{ijk\mathrm{}}^{(4)}`$ assume values from the set $`\{0,1,2,\mathrm{}n^41\}`$. Further, since the $`Q`$’s are pandiagonal it follows that the sum of the elements in any ROW of $`M^{(4)}`$ is given according to (3.5) by

$$n^3[n(n1)/2]+n^2[n(n1)/2]+n[n(n1)/2]+n(n1)/2=n(n^41)/2$$

and this is $`\sigma _4`$ in (1.3). Thus the sum of the elements of any ROW of $`M^{(4)}`$ is given by (1.3). Finally because of the assumed orthogonality of the $`Q`$’s, it follows that no two elements of $`M^{(4)}`$ are the same and $`M^{(4)}`$ is a magic pandiagonal 4-D cube without repeats. The theorem is established.

According to the discussion in IIC, since the $`Q`$’s are 4-D pandiagonal latin cubes, for which we found $`n17`$, it follows that $`M^{(4)}`$ will exist only for $`n17`$.

The number $`N_4`$ of different $`M^{(4)}`$’s in (3.5) can be estimated as above to vary for large values of $`n`$ as

$$N_4(n!)^4$$

so that for $`n=17,19`$ and 23, $`N_4`$ assumes the approximate and rapidly growing values of $`n`$ of $`1.6\times 10^{58}`$, $`2.19\times 10^{68}`$ and $`4.5\times 10^{89}`$ respectively.

It is also possible, using the tiling properties of the $`Q^{}s`$ to “tile” $`M^{(4)}`$ throughout four dimensional space and obtain a pandiagonal magic 4-D cube by selecting any $`n\times n\times n\times n`$ cube in this space. Such a cube will be the same as the original cube but with the elements of its underlying latin cubes permuted.

Evidentally, these arguments are extendable to dimensions higher than 4, but the resulting constraints on the $`\alpha `$ parameters, analogous to those in (2.9) in three dimensions and those in (2.14) in four, can be expected to become increasingly complex.

REFERENCES

1. Brian Alspach and Katherine Heinrich, “Perfect Magic Cubes of Order 4m”, The Fibonacci Quarterly, 19 (1981): 97-106.
2. W.S. Andrews, “Magic Squares and Cubes”, New York; Dover, 1960.
3. W. Rouse Ball and H.M.S. Coxeter “Mathematical Recreations Essays”, U. of Toronto Press, Cambridge, 1974; chapter VII.
4. J. Denes and A.D. Keedwell, Latin Squares and their Applications; New York Academic Press, 1974.
5. Maurice Kraitchik, “Mathematical Recreations”, New York, Dover, 1942.
6. Leonardi Euler, “De quatratis magicis”, Opera Omni, Série 1, 7 (1923) 441-457; and “Recherches sur un novelle espéce de quarres magiques,” Opera Omni, série 1, 7, (1923), 291-392.
7. Martin Gardner, “Mathematical Games”, Scientific American, 234, 1 (January 1976), 118-123.
8. Martin Gardner, “Mathematical Games”, Scientific American, 234, 2 (February 1976) 122-127.
9. Charles F. Laywine and Gary L. Mullen, “Discrete Mathematics Using Latin Squares”, New York, Wiley 1998.
10. Paul. C. Pasles, “The Lost Squares of Dr. Franklin”: Amer. Math. Monthly, 108, (2001) 489-511.
11. Barkeley Rosser and R.J. Walker, “The Algebraic Theory of Diabolic Magic Squares”, Duke Math. J., 5, (1939): 705-728.
12. Bayard E. Wynne, “Perfect Magic Cubes of Order Seven” J. Recreational Mathematics, 8:4, 1976, 285-293.