text:
# Decoupling and the radiatively-corrected MSSM Higgs sector This work was supported in part by the U.S. Department of Energy under the grant DE-FG03-92ER40689.

## 1 Introduction

Suppose a Higgs boson is discovered at the LHC, and its properties are observed to coincide (within experimental error) to those of the Standard Model (SM) Higgs boson. Moreover, suppose that evidence for supersymmetry is found, which suggests that in the minimal version of the model (MSSM), the observed Higgs boson is the lightest state of a two-Higgs-doublet model (2HDM). Finally, imagine that no evidence for the heavier Higgs states at the LHC is found. With precision measurements at a future high energy $`e^+e^{}`$ Linear Collider (LC), can one determine the mass scale of the heavier Higgs states?

At the LHC, more than one Higgs scalar of the MSSM will often be observed if $`\mathrm{tan}\beta 1`$ (due to enhanced Higgs couplings to down-type fermions). But, there is also a substantial region of moderate $`\mathrm{tan}\beta `$ in which only the lightest CP-even Higgs scalar ($`h`$) is observed . If the properties of $`h`$ approximate those of the SM Higgs boson ($`h_{\mathrm{SM}}`$), a program of precision Higgs measurements at the LC will play a critical role in elucidating the physics of Higgs bosons.

## 2 Decoupling Limit of the 2HDM

Given a non-minimal Higgs sector, the decoupling limit corresponds to the parameter regime in which all but one CP-even neutral Higgs scalar are significantly heavier than the $`Z`$. The properties of the lightest CP-even Higgs boson are nearly indistinguishable from those of the Standard Model (SM) Higgs boson. The decoupling limit is very general and exists in many multi-Higgs models. The MSSM Higgs sector provides a well-motivated example for the decoupling limit and is the main focus of this work.

First, consider the general 2HDM. The most general scalar potential is given by:

$`𝒱=m_{11}^2\mathrm{\Phi }_1^{}\mathrm{\Phi }_1+m_{22}^2\mathrm{\Phi }_2^{}\mathrm{\Phi }_2[m_{12}^2\mathrm{\Phi }_1^{}\mathrm{\Phi }_2+\mathrm{h}.\mathrm{c}.]`$
$`+\frac{1}{2}\lambda _1(\mathrm{\Phi }_1^{}\mathrm{\Phi }_1)^2+\frac{1}{2}\lambda _2(\mathrm{\Phi }_2^{}\mathrm{\Phi }_2)^2+\lambda _3(\mathrm{\Phi }_1^{}\mathrm{\Phi }_1)(\mathrm{\Phi }_2^{}\mathrm{\Phi }_2)`$
$`+\lambda _4(\mathrm{\Phi }_1^{}\mathrm{\Phi }_2)(\mathrm{\Phi }_2^{}\mathrm{\Phi }_1)+\frac{1}{2}[\lambda _5(\mathrm{\Phi }_1^{}\mathrm{\Phi }_2)^2+\mathrm{h}.\mathrm{c}.]`$
$`+\{[\lambda _6(\mathrm{\Phi }_1^{}\mathrm{\Phi }_1)+\lambda _7(\mathrm{\Phi }_2^{}\mathrm{\Phi }_2)]\mathrm{\Phi }_1^{}\mathrm{\Phi }_2+\mathrm{h}.\mathrm{c}.\}.`$ (1)
For simplicity, we assume that there is no explicit (or spontaneous) CP violation. Since the ground state must preserve U(1)EM, the scalar vacuum expectation values are $`\mathrm{\Phi }_i^0v_i/\sqrt{2}`$, with $`\mathrm{tan}\beta v_2/v_1`$ and $`v^2v_1^2+v_2^2=(246\mathrm{GeV})^2`$.

Diagonalizing the CP-even Higgs squared-mass matrix yields two CP-even scalar eigenstates:

$`h`$ $`=`$ $`(\sqrt{2}\mathrm{Re}\mathrm{\Phi }_1^0v_1)s_\alpha +(\sqrt{2}\mathrm{Re}\mathrm{\Phi }_2^0v_2)c_\alpha ,`$
$`H`$ $`=`$ $`(\sqrt{2}\mathrm{Re}\mathrm{\Phi }_1^0v_1)c_\alpha +(\sqrt{2}\mathrm{Re}\mathrm{\Phi }_2^0v_2)s_\alpha ,`$ (2)
where $`c_\alpha \mathrm{cos}\alpha `$, $`s_\alpha \mathrm{sin}\alpha `$. The other Higgs scalars of the model include a CP-odd state, $`A`$, and a charged scalar pair, $`H^\pm `$. The decoupling limit is defined as the limit of $`m_Am_Z`$, assuming $`\lambda _i<𝒪(1)`$. One can show that this limit corresponds to taking $`\beta \alpha \pi /2`$; i.e., $`\mathrm{cos}(\beta \alpha )0`$. In the approach to the decoupling limit, one finds :

$`m_A^2v^2\left[{\displaystyle \frac{\widehat{\lambda }}{c_{\beta \alpha }}}+\lambda _A\frac{3}{2}\widehat{\lambda }c_{\beta \alpha }\right],`$ (3)
$`m_h^2v^2(\lambda \widehat{\lambda }c_{\beta \alpha }),`$ (4)
$`m_H^2m_A^2+(\lambda \lambda _A+\widehat{\lambda }c_{\beta \alpha })v^2,`$ (5)
$`m_{H^\pm }^2=m_A^2+\frac{1}{2}(\lambda _5\lambda _4)v^2.`$ (6)
where $`c_{\beta \alpha }\mathrm{cos}(\beta \alpha )`$ and

$`\lambda \lambda _1c_\beta ^4+\lambda _2s_\beta ^4+\frac{1}{2}(\lambda _3+\lambda _4+\lambda _5)s_{2\beta }^2`$
$`+2s_{2\beta }(\lambda _6c_\beta ^2+\lambda _7s_\beta ^2),`$ (7)
$`\widehat{\lambda }\frac{1}{2}s_{2\beta }\left[\lambda _1c_\beta ^2\lambda _2s_\beta ^2(\lambda _3+\lambda _4+\lambda _5)c_{2\beta }\right]`$
$`\lambda _6c_\beta c_{3\beta }\lambda _7s_\beta s_{3\beta },`$ (8)
$`\lambda _Ac_{2\beta }(\lambda _1c_\beta ^2\lambda _2s_\beta ^2)+(\lambda _3+\lambda _4)s_{2\beta }^2`$
$`\lambda _5c_{2\beta }^2+2\lambda _6c_\beta s_{3\beta }2\lambda _7s_\beta c_{3\beta }.`$ (9)
In particular, eqs. (3)–(5) yield

$$\mathrm{cos}(\beta \alpha )\frac{\widehat{\lambda }v^2}{m_A^2\lambda _Av^2}\frac{\widehat{\lambda }v^2}{m_H^2m_h^2}.$$
(10)
It follow that: (i) $`m_h𝒪(m_Z)`$, (ii) $`m_Hm_Am_{H^\pm }`$, up to corrections of $`𝒪(m_Z^2/m_A)`$, and (iii) $`\mathrm{cos}(\beta \alpha )𝒪(m_Z^2/m_A^2)`$.

The couplings of the Higgs bosons to vector bosons, fermions and scalars typically depend on $`\alpha `$ and $`\beta `$. By examining the tree-level couplings of the lightest CP-even scalar $`h`$, one notes that in the limit of $`c_{\beta \alpha }=0`$, the couplings of $`h`$ reduce to the corresponding SM Higgs couplings. That is, in the approach to the decoupling limit, the properties of $`h`$ are nearly indistinguishable from those of the SM Higgs boson, whereas all the other Higgs states are significantly heavier. Thus, the effective low energy theory below the mass scale of $`𝒪(m_A)`$ is the Standard Model with one Higgs doublet.

The Higgs sector of the minimal supersymmetric extension of the Standard Model (MSSM) is a two-Higgs doublet model, with interactions constrained by supersymmetry . It particular,

$`\lambda v^2=m_Z^2c_{2\beta }^2,\lambda _4v^2=2m_W^2,`$
$`\widehat{\lambda }v^2=m_Z^2s_{2\beta }c_{2\beta },\lambda _Av^2=m_Z^2c_{4\beta },`$ (11)
and $`\lambda _5=\lambda _6=\lambda _7=0`$. Moreover, the following tree-level result can be derived:

$$\mathrm{cos}^2(\beta \alpha )=\frac{m_h^2(m_Z^2m_h^2)}{m_A^2(m_H^2m_h^2)}.$$
(12)
As expected, $`c_{\beta \alpha }0`$ in the decoupling limit where $`m_h𝒪(m_Z)`$ and $`m_Am_Z`$. Moreover, when $`m_Am_Z`$, eqs. (3)–(6) yield

$`m_h^2m_Z^2c_{2\beta }^2,m_H^2m_A^2+m_Z^2s_{2\beta }^2,`$
$`m_{H^\pm }^2=m_A^2+m_W^2,c_{\beta \alpha }{\displaystyle \frac{m_Z^4\mathrm{sin}^24\beta }{4m_A^4}}.`$ (13)
## 3 A SM-like Higgs boson without decoupling

It is possible for the theory to exhibit a SM-like Higgs boson without decoupling. For example, assume that $`\mathrm{tan}\beta 1`$. Two cases arise in the 2HDM where $`|c_{\beta \alpha }|1`$: (i) $`m_A^2\lambda _iv^2\mathrm{tan}\beta `$, with $`\lambda _i<𝒪(1)`$, and (ii) $`|\widehat{\lambda }|1`$ with $`m_A`$ arbitrary. (If $`\mathrm{tan}\beta 1`$, replace $`\mathrm{tan}\beta `$ with $`\mathrm{cot}\beta `$ above.) In the MSSM, the $`hb\overline{b}`$ coupling normalized to its SM value is given by:

$$\frac{\mathrm{sin}\alpha }{\mathrm{cos}\beta }=\mathrm{sin}(\beta \alpha )\mathrm{tan}\beta \mathrm{cos}(\beta \alpha ).$$
(14)
In case (i), we have $`|\mathrm{tan}\beta c_{\beta \alpha }|1`$ even if $`\mathrm{tan}\beta 1`$ \[see eq. (10)\], which implies that the $`hb\overline{b}`$ coupling is SM-like (corresponding to the decoupling limit). Case (ii) does not correspond to decoupling if $`m_A^2<𝒪(v^2)`$. If $`\mathrm{tan}\beta 1`$, it is possible to have $`|\mathrm{tan}\beta c_{\beta \alpha }|𝒪(1)`$ even if $`|c_{\beta \alpha }|1`$, in which case the $`hb\overline{b}`$ coupling deviates from its SM value. Nevertheless, for the Higgs couplings to $`t\overline{t}`$, vector bosons and scalars, $`h`$ is SM-like.

In the MSSM at tree-level, $`|c_{\beta \alpha }|1`$ is possible only in the decoupling regime \[see eq. (12)\], corresponding to case (i) above. However, one-loop effects mediated by supersymmetric particles can generate significant modifications to the tree-level MSSM Higgs sector. For example, the tree-level upper bound, $`m_hm_Z`$, can be significantly raised . Allowing for maximal mixing in the top squark sector and supersymmetric mass parameters of order $`M_S1`$ TeV, one finds a radiatively-corrected Higgs mass bound of $`m_h<135`$ GeV . We shall demonstrate in section 4 that for moderate values of $`m_A`$, radiative corrections to $`\widehat{\lambda }`$ can result in $`|\widehat{\lambda }|1`$ in certain regions of the MSSM parameter space. This would correspond to case (ii) above, and the $`h`$ of the MSSM could indeed exhibit SM-like properties outside the domain of the decoupling limit.

## 4 Radiatively-corrected MSSM Higgs couplings

In order to study the decoupling properties of the MSSM Higgs sector, it is crucial to examine the Higgs couplings, including the most significant loop-corrections. The leading contributions to the radiatively-corrected Higgs couplings arise in two ways. First, the radiative corrections to the CP-even Higgs squared-mass matrix results in a shift of the CP-even Higgs mixing angle $`\alpha `$ from its tree-level value. That is, the dominant Higgs propagator corrections can to a good approximation be absorbed into an effective (“radiatively-corrected”) mixing angle $`\alpha `$ . In this approximation, we can write:

$$^2\left(\begin{array}{cc}_{11}^2& _{12}^2\\ \\ _{12}^2& _{22}^2\end{array}\right)=_0^2+\delta ^2,$$
(15)
where the tree-level contribution is denoted by $`_0^2`$ and $`\delta ^2`$ is the contribution from the radiative corrections. Then, $`c_{\beta \alpha }`$ is given by

$$c_{\beta \alpha }=\frac{(_{11}^2_{22}^2)\mathrm{sin}2\beta 2_{12}^2\mathrm{cos}2\beta }{2(m_H^2m_h^2)\mathrm{sin}(\beta \alpha )}.$$
(16)
Inserting the tree-level values for $`_0^2`$, one can rewrite eq. (16) as

$`c_{\beta \alpha }={\displaystyle \frac{m_Z^2\mathrm{sin}4\beta }{2(m_H^2m_h^2)\mathrm{sin}(\beta \alpha )}}`$
$`+{\displaystyle \frac{(\delta _{11}^2\delta _{22}^2)\mathrm{sin}2\beta 2\delta _{12}^2\mathrm{cos}2\beta }{2(m_H^2m_h^2)\mathrm{sin}(\beta \alpha )}}.`$ (17)
Using tree-level Higgs couplings with $`\alpha `$ replaced by its effective one-loop value provides a useful first approximation to the radiatively-corrected Higgs couplings.

Second, contributions from the one-loop vertex corrections to tree-level Higgs-fermion couplings can modify these couplings in a significant way, especially in the limit of large $`\mathrm{tan}\beta `$. When radiative corrections are included, all possible dimension-four Higgs-fermion couplings are generated. In particular, the effects of higher dimension operators can be ignored if $`M_Sm_Z`$, which we henceforth assume. These results can be summarized by an effective Lagrangian that describes the coupling of the neutral Higgs bosons to the third generation quarks:

$`_{\mathrm{eff}}`$ $`=`$ $`(h_b+\delta h_b)\overline{b}_Rb_L\mathrm{\Phi }_1^0+(h_t+\delta h_t)\overline{t}_Rt_L\mathrm{\Phi }_2^0`$ (18)
$`+`$ $`\mathrm{\Delta }h_t\overline{t}_Rt_L\mathrm{\Phi }_1^0+\mathrm{\Delta }h_b\overline{b}_Rb_L\mathrm{\Phi }_2^0+\mathrm{h}.\mathrm{c}.,`$
resulting in a modification of the tree-level relation between $`h_t`$ \[$`h_b`$\] and $`m_t`$ \[$`m_b`$\] as follows :

$`m_b`$ $`=`$ $`{\displaystyle \frac{h_bv}{\sqrt{2}}}\mathrm{cos}\beta \left(1+{\displaystyle \frac{\delta h_b}{h_b}}+{\displaystyle \frac{\mathrm{\Delta }h_b\mathrm{tan}\beta }{h_b}}\right)`$ (19)
$`{\displaystyle \frac{h_bv}{\sqrt{2}}}\mathrm{cos}\beta (1+\mathrm{\Delta }_b),`$
$`m_t`$ $`=`$ $`{\displaystyle \frac{h_tv}{\sqrt{2}}}\mathrm{sin}\beta \left(1+{\displaystyle \frac{\delta h_t}{h_t}}+{\displaystyle \frac{\mathrm{\Delta }h_t\mathrm{cot}\beta }{h_t}}\right)`$ (20)
$`{\displaystyle \frac{h_tv}{\sqrt{2}}}\mathrm{sin}\beta (1+\mathrm{\Delta }_t).`$
The dominant contributions to $`\mathrm{\Delta }_b`$ are $`\mathrm{tan}\beta `$-enhanced, with $`\mathrm{\Delta }_b(\mathrm{\Delta }h_b/h_b)\mathrm{tan}\beta `$; whereas for $`\mathrm{tan}\beta 1`$, $`\delta h_b/h_b`$ provides a small correction to $`\mathrm{\Delta }_b`$. \[In the same limit, $`\mathrm{\Delta }_t\delta h_t/h_t`$, with the additional contribution of $`(\mathrm{\Delta }h_t/h_t)\mathrm{cot}\beta `$ providing a small correction.\]

From eq. (18) we can obtain the couplings of the physical neutral Higgs bosons to third generation quarks. The resulting couplings of $`h`$ to $`b\overline{b}`$ and $`t\overline{t}`$ pairs are given by:

$`g_{ht\overline{t}}`$ $`=`$ $`{\displaystyle \frac{m_t}{v}}{\displaystyle \frac{\mathrm{cos}\alpha }{\mathrm{sin}\beta }}\left[1{\displaystyle \frac{1}{1+\mathrm{\Delta }_t}}{\displaystyle \frac{\mathrm{\Delta }h_t}{h_t}}(\mathrm{cot}\beta +\mathrm{tan}\alpha )\right]`$
$`g_{hb\overline{b}}`$ $`=`$ $`{\displaystyle \frac{m_b}{v}}{\displaystyle \frac{\mathrm{sin}\alpha }{\mathrm{cos}\beta }}[1+{\displaystyle \frac{1}{1+\mathrm{\Delta }_b}}`$ (21)
$`\times ({\displaystyle \frac{\delta h_b}{h_b}}\mathrm{\Delta }_b)(1+\mathrm{cot}\alpha \mathrm{cot}\beta )].`$
We now turn to the decoupling limit. First consider the implications for the radiatively-corrected value of $`c_{\beta \alpha }`$. Since $`\delta _{ij}^2𝒪(m_Z^2)`$, and $`m_H^2m_h^2=m_A^2+𝒪(m_Z^2)`$, one finds

$$\mathrm{cos}(\beta \alpha )=c\left[\frac{m_Z^2\mathrm{sin}4\beta }{2m_A^2}+𝒪\left(\frac{m_Z^4}{m_A^4}\right)\right],$$
(22)
in the limit of $`m_Am_Z`$, where

$$c1+\frac{\delta _{11}^2\delta _{22}^2}{2m_Z^2\mathrm{cos}2\beta }\frac{\delta _{12}^2}{m_Z^2\mathrm{sin}2\beta }.$$
(23)
Equivalently, the radiative corrections have modified the tree-level definition of $`\widehat{\lambda }`$:

$$\widehat{\lambda }v^2=cm_Z^2\mathrm{sin}2\beta \mathrm{cos}2\beta .$$
(24)
Eq. (22) exhibits the expected decoupling behavior for $`m_Am_Z`$. However, eqs. (4) and (22) exhibit another way in which $`\mathrm{cos}(\beta \alpha )=0`$ can be achieved—simply choose the MSSM parameters (which govern the Higgs mass radiative corrections) such that the right hand side of eq. (4) vanishes. That is,

$$\mathrm{sin}2\beta =\frac{2\delta _{12}^2\mathrm{tan}2\beta \left(\delta _{11}^2\delta _{22}^2\right)}{2m_Z^2}.$$
(25)
This condition is equivalent to $`c=0`$, which implies $`\widehat{\lambda }=0`$ \[eq. (24)\] and thus corresponds to case (ii) in section 3. That is, although $`\widehat{\lambda }0`$ at tree-level, the one-loop radiative corrections to $`\widehat{\lambda }`$ can approximately cancel its tree-level value, resulting in $`|\widehat{\lambda }|1`$. (Note that the one-loop corrections arise from the exchange of supersymmetric particles, whose contributions can be enhanced for certain MSSM parameter choices. One can show that the two-loop corrections are subdominant, so that the approximation scheme is under control.) In particular, eq. (25) is independent of the value of $`m_A`$. Typically, eq. (25) yields a solution at large $`\mathrm{tan}\beta `$. That is, by approximating $`\mathrm{tan}2\beta \mathrm{sin}2\beta 2/\mathrm{tan}\beta `$, one can determine the value of $`\beta `$ at which $`\widehat{\lambda }0`$ :

$$\mathrm{tan}\beta \frac{2m_Z^2\delta _{11}^2+\delta _{22}^2}{\delta _{12}^2}.$$
(26)
Hence, there exists a value of $`\mathrm{tan}\beta `$ (which depends on the choice of MSSM parameters) where $`\mathrm{cos}(\beta \alpha )0`$ independently of the value of $`m_A`$. If $`m_A`$ is not much larger than $`m_Z`$, then $`h`$ is a SM-like Higgs boson outside the decoupling regime.

Finally, note that for $`m_Am_Z`$,

$$\mathrm{cot}\alpha =\mathrm{tan}\beta \left[1+\frac{2m_Z^2}{m_A^2}c_{2\beta }\right]+𝒪\left(\frac{m_Z^4}{m_A^4}\right).$$
(27)
Applying this result to eq. (21), it follows that in the decoupling limit, $`g_{hq\overline{q}}=g_{h_{\mathrm{SM}}q\overline{q}}=m_q/v`$. Away from the decoupling limit, the Higgs couplings to down-type fermions can deviate significantly from their tree-level values due to enhanced radiative corrections at large $`\mathrm{tan}\beta `$ \[where $`\mathrm{\Delta }_b𝒪(1)`$\]. In particular, because $`\mathrm{\Delta }_b\mathrm{tan}\beta `$, the leading one-loop radiative correction to $`g_{hb\overline{b}}`$ is of $`𝒪(m_Z^2\mathrm{tan}\beta /m_A^2)`$, which decouples only when $`m_A^2m_Z^2\mathrm{tan}\beta `$ (this behavior was called delayed decoupling in ).

## 5 Implications for precision Higgs measurements at the LC

As noted in section 1, a program of precision Higgs measurements at the LC may be critical in determining whether the properties of the lightest CP-even Higgs boson differ from those of $`h_{\mathrm{SM}}`$. In particular, recent simulations of Higgs branching ratio measurements  suggest that the Higgs couplings to vector bosons and the third generation fermions can be determined with an accuracy in the range of 1–3% at the LC. In the exact decoupling limit (of infinitely large $`m_A`$), $`h=h_{\mathrm{SM}}`$. However, for finite values of $`m_A`$, the fractional deviations of the couplings of $`h`$ relative to those of $`h_{\mathrm{SM}}`$ scale as $`m_Z^2/m_A^2`$. Thus, if precision measurements reveal a significant deviation from SM expectations, one could in principle derive a constraint (e.g., upper and lower bounds) on the heavy Higgs masses of the model.

In the MSSM, this constraint is sensitive to the supersymmetric parameters that control the radiative corrections to the Higgs couplings. This is illustrated in fig. 1, where the constraints on $`m_A`$ are derived for two different sets of MSSM parameter choices . Here, a simulation of a global fit of measured $`hbb`$, $`h\tau \tau `$ and $`hgg`$ couplings is made (based on the anticipated experimental accuracies given in ) and $`\chi ^2`$ contours are plotted indicating the constraints in the $`m_A`$$`\mathrm{tan}\beta `$ plane, assuming that a deviation from SM Higgs boson couplings is seen. In the maximal mixing scenario shown in fig. 1(a), the constraints on $`m_A`$ are significant and rather insensitive to the value of $`\mathrm{tan}\beta `$. However in some cases, as shown in fig. 1(b), a region of $`\mathrm{tan}\beta `$ may yield almost no constraint on $`m_A`$. This corresponds to the value of $`\mathrm{tan}\beta `$ given by eq. (26), and is a result of $`\widehat{\lambda }0`$ generated by radiative corrections \[$`c0`$ in eq. (24)\]. Thus, one cannot extract a fully model-independent upper bound on the value of $`m_A`$ beyond the kinematical limit that would be obtained if direct $`A`$ production were not observed at the LC.


format: text

text:
# 1 Introduction

## 1 Introduction

We will in this letter study the zero modes of tensile BPS strings in an interacting (2,0) theory in five space plus one time dimension. By zero modes we will mean modes that do not change the energy. We will thus not consider excited strings here. We show that it is possible to have such strings in a theory as solitons. The model we consider thus contains no fundamental strings. The only fields we consider are the massless ones which are in the (2,0) tensor multiplet, and which includes a self-dual tensor field $`H^+`$ plus five scalars $`\varphi ^a`$, $`a=1,\mathrm{},5`$. The strings being tensile means that we study the theory away from the origin of some M-theory moduli space, corresponding to separated M5-branes. We thus break the conformal invariance by giving non-zero expectation value to one of the five scalar fields. We then can expand the quantum fields about such classical solitons in the spirit of . This should presumably work if we consider energies much lower than $`D^{\frac{1}{2}}`$ as the work done in  may indicate. Here $`D`$ is the expectation value of the fifth scalar field at infinity. We can not make an expansion in a small coupling constant since we do not have any adjustable coupling constant in these six-dimensional theories.

The structure of this paper is as follows. In section 2 we present the model. In section 3 we show how to obtain charged solitonic BPS strings. In sections 3 we show that these strings can have fermionic zero modes, thus obtaining spinning strings. In section 4 we examine the bosonic zero modes of BPS strings and compute the metric on the moduli space. In section 5 we make contact with four dimensional results.

## 2 The model

We will throughout this letter assume that we have a flat six-dimensional Minkowski metric $`\eta _{MN}=diag(1,1,1,1,1,1)`$. We will take as our six-dimensional action the following

$$\frac{1}{2\pi }_W(\frac{1}{2}HHd\varphi ^ad\varphi _a+HA(\varphi ^a))+fermions,$$
(1)
where we will express the three-form potential $`A`$ (coming from eleven dimensional supergravity, having field strength $`F=dA`$) in terms of the five scalars living on the brane. That the term $`HA`$ should be there in the eleven dimensional context has been shown in . In  this term was used to show that an M2 brane can end on an M5 brane. We will assume that we have two parallel M5 branes, and we will factorize out the ‘center of mass’ tensor multiplet, so that we get the (spontaneously broken) $`A_1`$ theory described by one massless tensor multiplet (plus massive degrees of freedom, which we do not know how to describe properly).

We will now examine the effects that the second five brane $`W^{}`$ have on $`W`$. Small fluctuations of $`W^{}`$ can presumably be neglected when the two M5 branes are widely separated. But if there is a M2 brane stretching between the two M5 branes, then we can not trust such a model when being close to the boundary of such a M2 brane, which will be our selfdual string in the M5 brane. But far away from such strings we should be able to solve $`A`$ to a good approximation from

$$dF=2\pi \delta _W^{}$$
(2)
with $`\delta _W^{}`$ being the Poincare dual of a flat five-brane world volume $`W^{}`$. It has the solution $`F=2\pi \frac{1}{V_4}\mathrm{\Omega }_4`$ where $`\mathrm{\Omega }_4`$ and $`V_4`$ are the volume form and volume respectively of the unit four-sphere that surrounds $`W^{}`$. If we let

$`\varphi ^5`$ $`=`$ $`r\mathrm{cos}\theta `$ (3)
$`\varphi ^4`$ $`=`$ $`r\mathrm{sin}\theta \mathrm{cos}\varphi `$ (4)
$`\varphi ^3`$ $`=`$ $`r\mathrm{sin}\theta \mathrm{sin}\varphi \mathrm{cos}\psi `$ (5)
$`\varphi ^2`$ $`=`$ $`r\mathrm{sin}\theta \mathrm{sin}\varphi \mathrm{sin}\psi \mathrm{cos}\chi `$ (6)
$`\varphi ^1`$ $`=`$ $`r\mathrm{sin}\theta \mathrm{sin}\varphi \mathrm{sin}\psi \mathrm{sin}\chi `$ (7)
denote the five transverse coordinates to $`W^{}`$, then

$`\mathrm{\Omega }_4`$ $`=`$ $`{\displaystyle \frac{1}{4!}}ϵ_{abcde}{\displaystyle \frac{1}{r^5}}\varphi ^ad\varphi ^bd\varphi ^cd\varphi ^dd\varphi ^e`$ (8)
$`=`$ $`\mathrm{sin}^3\theta \mathrm{sin}^2\phi \mathrm{sin}\psi d\theta d\phi d\psi d\chi `$ (9)
$`=`$ $`\mathrm{sin}^3\theta d\theta \mathrm{\Omega }_3.`$ (10)
Here $`\mathrm{\Omega }_3`$ is the volume form on the equator $`\varphi ^5=0`$. Locally $`F=dA`$. We may for instance choose $`A=\frac{2\pi }{V_4}a(\theta )\mathrm{\Omega }_3`$ for some function $`a(\theta )`$ such that $`da=\mathrm{sin}^3\theta d\theta `$. Such solutions can only be obtained locally by excluding either the south pole ($`\theta =\pi `$) or the north pole ($`\theta =0`$). They are given respectively by

$$A_\pm =\frac{2\pi }{V_4}\left(\pm \frac{2}{3}\mathrm{cos}\theta +\frac{1}{3}\mathrm{cos}^3\theta \right)\mathrm{\Omega }_3.$$
(11)
In Cartesian coordinates (with $`A,B,\mathrm{}=1,2,3,4`$) we get

$`A_\pm `$ $`=`$ $`{\displaystyle \frac{\pi (\varphi ^5\pm 2r)}{3V_4r^3(\varphi ^5\pm r)^2}}\epsilon _{ABCD}\varphi ^Ad\varphi ^Bd\varphi ^Cd\varphi ^D,`$ (12)
so we see that $`A_{}`$ ($`A_+`$) is singular along the positive (negative) $`\varphi ^5`$-axis.

In order for the anti-selfdual piece of $`H`$ to decouple we should define the field strength as

$$H=dB+A(\varphi )$$
(13)
where $`B`$ is the two-form gauge potential. This definition implies the Bianchi identity

$$dH=F(\varphi ).$$
(14)
The equations of motion one derives from this action are

$`dH`$ $`=`$ $`{\displaystyle \frac{2\pi }{4!V_4}}ϵ_{abcde}{\displaystyle \frac{1}{r^5}}\varphi ^ad\varphi ^bd\varphi ^cd\varphi ^dd\varphi ^e`$ (15)
$`dd\varphi _a`$ $`=`$ $`{\displaystyle \frac{\pi }{3!V_4}}ϵ_{abcde}{\displaystyle \frac{1}{r^5}}(H+H)\varphi ^bd\varphi ^c\varphi ^dd\varphi ^e.`$ (16)
These equations have also been considered in . We see that only $`H^+=\frac{1}{2}(H+H)`$ couples to the current $`F(\varphi )`$. The equations for the five scalars can be derived as follows. The variation of the Lagrangian is

$$\frac{1}{2\pi }[2\delta \varphi ^add\varphi ^a\delta A(\varphi )(H+H)].$$
(17)
We choose the gauge corresponding to $`A_{}`$ and define $`\rho =\sqrt{\varphi ^A\varphi _A}`$ ($`A=1,2,3,4`$) and compute the variation with respect to $`\varphi ^5`$,

$$\delta A(\varphi )=\frac{2\pi }{V_4}\delta (a(\theta )\mathrm{\Omega }_3)=\frac{2\pi }{V_4}\mathrm{sin}^3\theta \delta \theta \mathrm{\Omega }_3=\frac{2\pi }{V_4}\mathrm{sin}^2\theta \delta (\mathrm{cos}\theta )\mathrm{\Omega }_3=\frac{2\pi }{V_4}\delta \varphi ^5\frac{1}{r^5}\rho ^4\mathrm{\Omega }_3.$$
(18)
Hence the corresponding equation of motion is

$$2dd\varphi ^5=\frac{2\pi }{V_4}(H+H)\frac{1}{r^5}\rho ^4\mathrm{\Omega }_3.$$
(19)
where $`\rho ^4\mathrm{\Omega }_3=\frac{1}{3!}ϵ_{5bcde}\varphi ^bd\varphi ^cd\varphi ^dd\varphi ^e`$. The equation of motion for any other $`\varphi ^a`$ is most easily derived by choosing the gauge in which $`A`$ is singular along the positive $`\varphi ^a`$-axis. The equations of motion are of course gauge invariant, so this is a legitimate way to proceed.

We notice that the term $`HA(\varphi )`$ in the Lagrangian depends linearly on $`\dot{B}`$. This implies that this term vanishes under a Legendre transformation. If we choose the gauge $`A_{0IJ}=B_{0I}=0`$ ($`I,J=1,2,\mathrm{},5`$) then the conjugate momenta of $`B_{IJ}`$ become

$$2\pi E^{IJ}=\frac{1}{2}H^{0IJ}$$
(20)
Letting $`H^\pm =\frac{1}{2}(1\pm )H`$, the bosonic Hamiltonian in the rest frame becomes

$`2\pi M`$ $`=`$ $`{\displaystyle }(H^+_5H^++H^{}_5H^{}+d\varphi ^a_5d\varphi ^a).`$ (21)
Rewriting the mass as

$$2\pi M=d^5x\left(\right|H^+\pm _5(d\varphi ^5dx_5)|^2+(_5\varphi ^5)^2+|H^{}|^2+|d\varphi ^A|^2)2H^+d\varphi ^5dx_5,$$
(22)
we obtain a lower bound on the mass given by

$$2\pi M2\left|H^+d\varphi ^5dx_5\right|.$$
(23)
This is not the most general lower bound that we can have. The situation that we have described is the bound that corresponds to a single selfdual string in the $`x^5`$-direction. In a supersymmetric theory we know that this bound should be given through a central charge in the supersymmetry algebra. Although the supersymmetry generators may have to be corrected by complicated nonlinear terms in order to be able to close the supersymmetry algebra on the new shell that we get when adding nonlinear interaction terms to the action, we were able to deduce that the central charge one gets still should be given by the simple expression that one derives in the free theory, namely $`Z_\mu ^a=2H^+d\varphi ^adx_\mu `$. In the bound above this central charge has rank one (viewed as a matrix in the indices $`\mu `$ and a). In general it could have higher rank, corresponding to several non-parallel strings. Now if we choose the central charge as $`2QD𝑑x^5`$ with $`D>0`$ we get the following Bogomolnyi equations which have also been analysed in ,

$`H^+`$ $`=`$ $`sign(Q)_4d\varphi ^5`$ (24)
$`d\varphi ^A`$ $`=`$ $`0`$ (25)
$`_5\varphi ^5`$ $`=`$ $`0`$ (26)
$`H^{}`$ $`=`$ $`0.`$ (27)
Of course we can always rotate the coordinate system so that the rank one central charge will be aligned in the $`x^5`$-direction.

## 3 A classical BPS string solution

We will look for a BPS solution where the scalar field $`\varphi ^5`$ is spherically symmetric, i.e. depend only on $`R=\sqrt{x^Ax_A}`$, ($`A=1,2,3,4`$). All fields are assumed to be translationally invariant in the $`x^0,x^5`$ directions (which we will associate to the directions of a string world-sheet). We will assume the asymptotic behaviour

$`\varphi ^5(\mathrm{})`$ $`=`$ $`D>0`$ (28)
$`\varphi ^A(\mathrm{})`$ $`=`$ $`0`$ (29)
$`\varphi ^A(0)`$ $`=`$ $`0.`$ (30)
Of course it now follows from the Bogomolnyi equations that $`\varphi ^A=0`$ everywhere. But in that case we can not find any non-trivial solution. We therefore first look for a charged solution, which is not necessarily BPS. The fact that the solution is charged and thereby topologically non-trivial prevents it from becoming trivial. So then we can take the limit that $`\varphi ^A`$ tends to zero to obtain a solution that is BPS.

We may integrate up the equation of motion $`dH=F`$ if $`A`$ is globally defined on $`W`$ by using Stokes theorem. If the map $`x^A\frac{\varphi ^A}{\rho }`$, $`\rho =\sqrt{\varphi ^A\varphi ^A}`$, has winding number $`N`$ then we get

$$_{S^3}H=\frac{2\pi }{V_4}\left(\frac{2}{3}\mathrm{cos}\theta +\frac{1}{3}\mathrm{cos}^3\theta \right)NV_3.$$
(31)
We have thus assume that the point $`(\rho ,\varphi ^5)=(0,D)`$ never is reached for any finite $`R`$ and that $`\varphi ^5<0`$ at $`R=0`$ where $`\rho =0`$. This means that we may use the potential $`A_{}`$ which is singular on the positive $`\varphi ^5`$-axis. This potential will then be non-singular all over the M five-brane world-volume. Furthermore, if we compute the charge $`Q`$ by integrating $`H`$ over a three-sphere with infinite radius $`R`$ we get $`Q=_{S_{\mathrm{}}^3}H=2\pi N`$. We can now consider a limiting process which will not alter this topological twist to obtain a BPS configuration. We may for instance take the square integrable function $`\rho =ϵRe^R`$ and let $`ϵ0`$. Using the Bogomolnyi equation which relates $`H`$ and $`d\varphi ^5`$ we get,

$$R^3\varphi _{}^{5}{}_{}{}^{}(R)=|Q|\frac{1+sign(\varphi ^5)}{2V_3}+Ordo(ϵ^4).$$
(32)
In the limit $`ϵ0`$ we then get the BPS solution

$`\varphi ^5`$ $`=`$ $`\left(D{\displaystyle \frac{|Q|}{2V_3R^2}}\right)\mathrm{\Theta }\left(R\sqrt{{\displaystyle \frac{|Q|}{2V_3D}}}\right)`$ (33)
$`\varphi ^A`$ $`=`$ $`0.`$ (34)
Here $`\mathrm{\Theta }`$ is the Heaviside step function. This is thus the exact solution in our model. But our model is not realistic when $`\varphi ^5`$ is close to zero. The fact that charge quantization comes about from the topological twist of the scalars around the string should however be true no matter how many further interaction terms we add to our action to make it more realistic. One argument for why one should think so, apart from the obvious one that we got the correct Dirac charge quantization in this way, is that the charge should be a topological invariant (since it can not jump under a continuous deformation) and we may deform our BPS solution continuously in such a way that the two M5 branes never come close to each other. We can also deform $`W^{}`$ away from its true stable configuration into the flat brane that we have based our model on without changing the topology. In that situation we may rely on our model. That the charge corresponds to a winding number of the scalars has been anticipated earlier , , but perhaps not been shown this explicitly.

We may insert this solution into $`dH`$ to compute $`F`$. We then find that

$$dH=Q\delta \left(R\sqrt{\frac{|Q|}{2V_3D}}\right)dR\frac{\mathrm{\Omega }_3}{V_3}$$
(35)
If we have two strings far away from each other then we should be able to superimpose two source terms like that above, with one such source term centered at the each string. So although the equations we started with were non-linear, when restricting to BPS solutions the equations become linear and the superposition principle holds, at least when the strings are widely separated.

## 4 The fermionic zero modes

The equation of motion for the fermions can be obtained by making a supersymmetry transformation of the bosonic equations. Imposing the supersymmetry variations (which should be correct up to linear order in the fields)

$`\delta B_{MN}`$ $`=`$ $`ϵ\gamma _{MN}\psi `$ (36)
$`\delta \varphi ^a`$ $`=`$ $`ϵ\sigma ^a\psi `$ (37)
(where $`\gamma ^M`$ and $`\sigma ^a`$ are gamma-matrices, $`ϵ`$ a constant anticommuting symplectic Majorana spinor parameter and $`\psi `$ the symplectic Majorana spinor field which is in the $`(2,0)`$ tensor multiplet) we get from the supersymmetry variation of the bosonic equations of motion that

$$\gamma ^M_M\psi =\frac{\pi }{75V_4}ϵ_{abcde}\frac{1}{r^5}\varphi ^a_M\varphi ^b_N\varphi ^c_P\varphi ^d\gamma ^{MNP}\sigma ^e\psi .$$
(38)
The right-hand side will be a source term similar to (35) for BPS strings. If we consider a soliton in which all fermionic fields are put equal to zero, then of course the supersymmetry variations of the bosonic fields vanish. The supersymmetry variation of the fermions is (to linear order)

$$\delta \psi =\left(\frac{1}{24}\gamma ^{MNP}H_{MNP}^++\frac{1}{2}\sigma _a\gamma ^M_M\varphi ^a\right)ϵ.$$
(39)
When inserting our classical bosonic soliton solution into this variation (thus following the same semi-classical arguments as in ) we obtain four (complex) fermionic zero modes (corresponding to the eight real supersymmetries which are broken by our string soliton),

$$\psi _0=\frac{1}{\sqrt{|Q|DV_3𝑑x^5}}\gamma ^{05I}H_{05I}u$$
(40)
($`I=1,2,3,4`$) where $`u`$ is a constant spinor that survives the projection $`P=\frac{1}{2}(1+\sigma _5\gamma ^0\gamma ^5)`$. This is a projector onto the broken supersymmetries. There are four complex solutions $`u`$ that survives this projection. These are static solutions of the equation of motion and hence zero modes, but in general there could exist more fermionic zero modes. It seems plausible that for winding numbers $`\pm 1`$ of the scalars around the string, i.e. for the charges $`Q=\pm 2\pi `$, these are all the fermionic zero modes. Expanding $`\psi =\psi _0`$ \+ (non-zero modes) in an orthonormal basis, the equal-time canonical quantization relations $`\{\psi (x),\psi (y)^{}\}=\delta ^5(xy)/4`$ (together with the completeness and orthonormality relations of the basis) implies that $`\{u,u^{}\}=1/4`$. We may write $`u=a_iu_i`$ where $`a_i`$ are operators obeying $`\{a_i,a_j^{}\}=\delta _{ij}/4`$ and $`u_i`$ ($`i=1,\mathrm{},4`$) are c-number valued constant complex orthonormal spinors.

We may insert our classical static soliton solution into the supercharges (they are given in  in terms of the fields), and find that the supercharges are equal to $`4\sqrt{|Q|D𝑑x^5}Pu`$, so the supersymmetry algebra will be obeyed. Since the Hamiltonian commutes with the supercharges we can construct new string states which have the same energy by acting on the soliton state with vanishing $`SO(4)`$-spin 111The string breaks the Lorentz group $`SO(1,5)`$ down to $`SO(4)`$ by any one of the following 16 creation operators $`1,a_i^{},a_i^{}a_j^{},a_i^{}a_j^{}a_k^{},a_1^{}a_2^{}a_3^{}a_4^{}`$ to obtain the short BPS multiplet of spinning string states that was constructed in .

## 5 The moduli space of two bosonic strings

The minimal tension (we assume that we have translational invariance in the $`x^0`$ and $`x^5`$ directions) in a system with fields $`H_T`$ and $`\varphi _T`$ is given by

$$2\pi T=2\left|_{R^4}H_Td\varphi _T\right|$$
(41)
and this is so only if $`H_T=\pm d\varphi _T`$ for some sign choice. Now if these fields are produced by a system of two strings, each of which give rise to the fields $`H`$ and $`\stackrel{~}{H}`$ where $`\stackrel{~}{H}(x)=cH(xR)`$, c being some constant, is the field produced by the second string separated from the first string by a distance $`\stackrel{}{R}`$, then since the superposition principle holds outside the strings where the source terms vanish, we get the total field as $`H_T=H+\stackrel{~}{H}`$. Charge quantization requires $`c`$ to be a rational number. Similarly $`d\stackrel{~}{\varphi }=c^{}d\varphi `$ for some constant $`c^{}`$. The two strings are BPS separately only if

$`H`$ $`=`$ $`sign(Q)d\varphi `$ (42)
$`\stackrel{~}{H}`$ $`=`$ $`sign(\stackrel{~}{Q})d\stackrel{~}{\varphi }.`$ (43)
Hence $`c^{}=sign(Q)sign(\stackrel{~}{Q})c`$. The whole system will then in general not be BPS since now we get

$$H_T=sign(Q)d\varphi +sign(\stackrel{~}{Q})d\stackrel{~}{\varphi }$$
(44)
and this is equal to $`sign(\stackrel{~}{\stackrel{~}{Q}})(d\varphi +d\stackrel{~}{\varphi })`$ for some $`\stackrel{~}{\stackrel{~}{Q}}`$ (which will be equal to $`\stackrel{~}{Q}+Q`$), as is necessary for the configuration to be BPS, if and only if $`sign(Q)=sign(\stackrel{~}{Q})`$.

The tension of one of the strings we define as follows. We first notice that

$$H_Td\varphi _T=Hd\varphi _T+\stackrel{~}{H}d\varphi _T$$
(45)
and that we may compute these integrals by performing the integrations as follows,

$$Hd\varphi _T=𝑑\varphi _T_{\varphi _T=const}H$$
(46)
The level curve $`\varphi _T=const`$ will for sufficiently small $`\varphi _T`$ consist of disconnected components, each of which encloses just one string. But only that closed contour that surrounds the string that is producing the field $`H`$ will give a contribution to the integral $`_{\varphi _T=const}H`$ and this will have the value $`Q`$. This value is independent of $`\varphi _T`$. So all in all we get the result

$$H_Td\varphi _T=\underset{i}{}Q_iD$$
(47)
where the sum runs over all strings, having charges $`Q_i`$ respectively, and $`D=lim_R\mathrm{}\varphi _T`$. It is now natural to interpret $`2|Q_i|D`$ as the tension of string number $`i`$ in this configuration. (The factor $`2`$ is just a convention.)

Furthermore, the contributions $`D_i=lim\varphi `$ from each string separately must add up to $`D`$. From $`c^{}=c`$ it follows that $`D_i=\frac{|Q_i|}{|_jQ_j|}D`$.

As a consequence of the zero force property of BPS configurations we have a moduli space for two parallel BPS strings. Since we found four fermionic zero modes we expect of a supersymmetric theory to find exactly four bosonic zero modes as well (and not more, if we ignore center of mass motion). We therefore expect the asymptotic moduli space to be $`R^4`$ where a point in $`R^4`$ corresponds to the relative separation of the two strings. It might not be so interesting to scatter infinitely long and therefore infinitely heavy strings against each other. But if we compactify the direction along which they are stretching on a circle then their masses become finite, and it could be of interest to compute the metric on this moduli space. In order to do that we first notice that the bosonic fields produced by a string, outside that region where the fields are zero, happens to be, in our model, exactly reproducible by an action like222If this action also can be used for fluctuating strings is not clear, but we will not consider fluctuating strings in this letter. As the strings appears to look like singular objects in this action, one might wonder if these solitonic strings really could be fundamental objects. But as the strings are self-dual nothing essential changes under electric-magnetic duality. There is no way to escape the magnetic coupling. In Yang-Mills theory that would have implied that such objects would have to be solitonic. But of course this is not a Yang-Mills theory, so we should not exclude any possibilities.

$`{\displaystyle }{\displaystyle \frac{1}{2}}HH|Q|{\displaystyle _\mathrm{\Sigma }}B`$ (48)
$`{\displaystyle }d\varphi d\varphi 2|Q|{\displaystyle _\mathrm{\Sigma }}d^2\sigma \sqrt{g}\varphi `$ (49)
where $`\mathrm{\Sigma }`$ is the string world-sheet oriented such that $`_{R^4}\delta _\mathrm{\Sigma }=sign(Q)`$, and $`H=dB+|Q|\delta _D`$ with $`D`$ being a Dirac brane extending from $`\mathrm{\Sigma }`$ out to infinity. The induced metric on $`\mathrm{\Sigma }`$ is denoted $`g_{ab}`$ and $`g=det(g_{ab})`$. The peculiar normalization of the kinetic term for the scalars follows from our conventions for the supersymmetry variations . In this action there is implicitly a coupling term to the magnetic charge as well, coming from333Our definition of Poincare dual is $`\omega \delta _D=_D\omega `$. Then we must define $`D=\mathrm{\Sigma }`$ if we want $`d\delta _D=\delta _\mathrm{\Sigma }`$

$$\frac{1}{2}HH=\frac{1}{2}dBdB|Q|\delta _DH=\frac{1}{2}dBdB+|Q|_DH.$$
(50)
If we define the dual potential $`\stackrel{~}{B}`$ as $`d\stackrel{~}{B}=dB`$, then the last term can be written as

$$|Q|_\mathrm{\Sigma }\stackrel{~}{B}.$$
(51)
Due to self-duality we can on-shell simply choose $`B=\stackrel{~}{B}`$ when we are away from the Dirac-brane $`D`$, and we may simply collect the coupling terms to a term that is twice the electric coupling term. It is now important to note that this action, in which we have used the self-duality of $`H`$, then can be used only to derive the equations of motion for the strings which move in a given on-shell $`B`$-field. It can then not be used to derive the equations for the $`B`$-field.

As we chose the orientations according to the sign of the charges, this is exactly the same as saying that the whole system is BPS if and only if the two strings are parallel. But $`Q`$ and $`\stackrel{~}{Q}`$ need not be equal. It suffices that they have same sign. But if $`Q=2\pi N`$ for some $`N\pm 1`$, then we think that we can view this as $`N`$ copies of a string with $`N=\pm 1`$ which are on top of each other, and which then could be separated from each other without any energy cost. So we will in the sequel assume that $`N=\pm 1`$.

We now set out to compute the metric on the asymptotic part of the moduli space which is $`R^4`$ by following the same method as in . The fields produced by a string can be expressed in terms of a Greens function $`J(x)`$ as follows

$`B_{MN}(x)`$ $`=`$ $`|Q|{\displaystyle _\mathrm{\Sigma }}𝑑X_MdX_NJ(xX)+|Q|ϵ_{MNPQRS}{\displaystyle _D}𝑑Y^PdY^QdY^R{\displaystyle \frac{}{x^S}}J(xY)`$ (52)
$`\varphi (x)`$ $`=`$ $`|Q|{\displaystyle _\mathrm{\Sigma }}d^2\sigma \sqrt{g}J(xX(\sigma )).`$ (53)
where $`J(x)`$ satisfies $`^M_MJ(x)=\delta ^6(x)`$. The retarded solution to this equation is

$$J(x^0,r)=\frac{1}{3V_3}\left[\frac{1}{r^3}\delta (x^0r)+\frac{1}{r^2}\delta ^{}(x^0r)\right]$$
(54)
where $`r=\sqrt{R^2+(x^5)^2}`$ and $`R=\sqrt{(x^1)^2+(x^2)^2+(x^3)^2+(x^4)^2}`$. The potential above will be singular on the Dirac-brane $`D`$.

If string number one moves with a velocity $`v_1`$ which is much less than one (the speed of light), and we parameterize $`\mathrm{\Sigma }_1`$ as

$`X^0`$ $`=`$ $`\tau `$ (55)
$`X^I`$ $`=`$ $`v_{1}^{}{}_{}{}^{I}\tau `$ (56)
$`X^5`$ $`=`$ $`\sigma `$ (57)
and the associated Dirac brane $`D`$ as

$`Y^0`$ $`=`$ $`\tau `$ (58)
$`Y^{1,2,3}`$ $`=`$ $`v_{1}^{}{}_{}{}^{1,2,3}\tau `$ (59)
$`Y^4`$ $`=`$ $`\stackrel{~}{\sigma }+v_{1}^{}{}_{}{}^{4}\tau `$ (60)
$`Y^5`$ $`=`$ $`\sigma ,`$ (61)
then the scalar field produced by $`\mathrm{\Sigma }_1`$ becomes

$$\varphi =D_1\sqrt{1\stackrel{}{v_1}^2}\frac{|Q_1|}{2V_3S^2}+Ordo((v_1)^3)$$
(62)
where $`S=R(1+Ordo(v_1))`$. Similarly, the retarded potential becomes

$`B_{05}`$ $`=`$ $`{\displaystyle \frac{Q_1}{2V_3S^2}}+Ordo((v_1)^3)`$ (63)
$`B_{IJ}`$ $`=`$ $`\omega _{IJ}+Ordo(v_1)`$ (64)
$`B_{0J}`$ $`=`$ $`v^I\omega _{IJ}+Ordo((v_1)^2)`$ (65)
$`B_{I5}`$ $`=`$ $`v_{1}^{}{}_{I}{}^{}B_{05}.`$ (66)
Here $`\omega _{MN}`$ is the potential generated by a static magnetically charged string with Dirac-brane singularity along $`D`$. We assume that we can neglect the accelaration. It is then exactly the same quantity $`S`$ that enters in all these expressions. This is a crucial observation.

The expressions above which involve $`\omega `$ turn out to be irrelevant for us since these components are not in the direction of $`\mathrm{\Sigma }`$. (This is the essential difference as compared to the four dimensional case, which there gave an interpretation of electric charge as a fourth moduli . Here we are lucky that we do not get any contribution from the magnetic components as we do not want any more moduli than the four that we already have found.)

Notice that we can make periodic $`x^5`$ and still get the same retarded solutions as above because that solution is constant in the $`x^5`$ direction and hence satisfy both the correct equations of motion as well as any periodic boundary conditions that we might wish to impose in the $`x^5`$-direction.

Inserting these fields into the action for a second string parallel to the first (which implies that their charges have the same sign) and at a distance $`\stackrel{}{R}`$ from it, we get (now we can replace $`S`$ by $`R`$ as we see that the $`Ordo(v_1)`$-corrections cancel, so that the effect of replacing $`S`$ by $`R`$ will only give corrections of $`Ordo(v_{1}^{}{}_{}{}^{3})`$)

$`|Q_2|\left(2{\displaystyle d^2\sigma \sqrt{g}\varphi }2{\displaystyle B}\right)`$ (67)
$`=`$ $`2|Q_2|{\displaystyle d^2\sigma \left[D_1\sqrt{1\stackrel{}{v_2}^2}+\frac{|Q_1|}{R^2}\left(\sqrt{1\stackrel{}{v}_{1}^{}{}_{}{}^{2}}\sqrt{1\stackrel{}{v}_{2}^{}{}_{}{}^{2}}+1+\stackrel{}{v}_1\stackrel{}{v}_2\right)\right]}.`$ (68)
Upon adding the kinetic term of $`\mathrm{\Sigma }_2`$, which has $`D_2=\frac{|Q_2|}{|Q_1+Q_2|}D`$, the center of mass motion decouples and the action governing the relative motion is given by

$$d^2\sigma \left(\frac{Q_1Q_2}{2|Q_1+Q_2|}D\frac{Q_1Q_2}{\stackrel{}{R}\stackrel{}{R}}\right)\dot{\stackrel{}{R}}\dot{\stackrel{}{R}}$$
(69)
The strings will move in geodesics on this $`R^4`$ with metric

$$(ds)^2=\left(\frac{Q_1Q_2}{2|Q_1+Q_2|}D\frac{Q_1Q_2}{\stackrel{}{R}^2}\right)dR^idR^i$$
(70)
where $`\stackrel{}{R}`$ is the relative separation between the strings. This is a valid approximation only when $`|\stackrel{}{R}|>>D^{\frac{1}{2}}`$.

### 5.1 Reduction to four dimensions

The moduli space of these strings should be connected to the moduli space in four dimensions . At least we should be able to reduce to the moduli space for two monopoles that are far away from each other by following the route in . We can get the four-dimensional action that was considered in  by compactifying $`x^4`$ and $`x^5`$ on a small $`T^2`$. We thus make the identifications $`x^4x^4+qR`$ and $`x^5x^5+gR`$. We consider two parallel strings on $`T^2`$ which are aligned along the ‘diagonal’ (that is, they are straight, and wind around both the cycles arbitrarily many turns). We parameterize the string world-sheets $`\mathrm{\Sigma }_i`$ $`(i=1,2)`$ as follows

$`X^0`$ $`=`$ $`\tau `$ (71)
$`X^i`$ $`=`$ $`v^i\tau `$ (72)
$`X^4`$ $`=`$ $`nRq\sigma +v^4\tau `$ (73)
$`X^5`$ $`=`$ $`mRg\sigma ,`$ (74)
$`0\sigma 1`$ and $`n,m`$ being integers. (Thus they are winding numbers of the string around the a- and b-cycles). We will call such a string an (n,m)-string. The four-dimensional electric and magnetic potentials should now be defined as

$`A_\mu `$ $`=`$ $`RB_{\mu 4}`$ (75)
$`\stackrel{~}{A}_\mu `$ $`=`$ $`RB_{\mu 5}`$ (76)
($`R`$ is assumed to have dimension length, and we choose $`q`$ and $`g`$ such that $`qg=2\pi \mathrm{}`$). We will only care about the coupling terms to the string (or particle upon dimensional reduction) and will not worry about getting a the Lagrangian whose variation with respect to $`B`$ yields the correct equation of motion for the gauge field. Instead we will assume that the gauge field already satisfies the equation of motion. That is, we assume that $`H`$ is self-dual. That will in turn relate $`A`$ and $`\stackrel{~}{A}`$ as each others dual potentials.

If we let $`\mathrm{\Sigma }=\gamma \times S^1`$ where $`\gamma `$ thus describes the trajectory of $`X^\mu (\tau )`$ for $`\mu =0,1,2,3`$ and $`S^1`$ describes the trajectory of the coordinates $`X^{4,5}(\tau ,\sigma )`$ around the $`T^2`$, then the six-dimensional coupling of a string to an on-shell gauge field reduces to

$$_\mathrm{\Sigma }B=_\gamma \left\{nqA+mg\stackrel{~}{A}+v^4nmqgR^2B_{45}\right\}.$$
(77)
The last term is nothing that we recognize, but it will turn out to be irrelevant when we take the limit $`R0`$.

If we use the notation $`\dot{X}^\mu =\frac{X^\mu }{\tau }`$, $`X_{}^{}{}_{}{}^{\mu }=\frac{X^\mu }{\sigma }`$, $`(\dot{X})^2=\dot{X}_\mu \dot{X}^\mu `$, and let $`a,b,\mathrm{}=0,1`$ correspond to $`\tau `$ and $`\sigma `$ respectively, then the induced metric on the world-sheet becomes

$$g_{ab}=\left(\begin{array}{cc}(\dot{X})^2+(\dot{X}^4)^2& \dot{X}^4X_{}^{}{}_{}{}^{4}\\ \dot{X}^4X_{}^{}{}_{}{}^{4}& (X_{}^{}{}_{}{}^{4})^2+(X_{}^{}{}_{}{}^{5})^2\end{array}\right).$$
(78)
Then

$$g:=det(g_{ab})=\left((X_{}^{}{}_{}{}^{4})^2+(X_{}^{}{}_{}{}^{5})^2\right)(\dot{X})^2+\left(\dot{X}^4X_{}^{}{}_{}{}^{5}\right)^2.$$
(79)
We can now reduce the world-sheet action as follows

$$_\mathrm{\Sigma }\sqrt{g}\varphi =_\gamma 𝑑\tau \varphi \sqrt{\left((X_{}^{}{}_{}{}^{4})^2+(X_{}^{}{}_{}{}^{5})^2\right)(\dot{X})^2\left(\dot{X}^4X_{}^{}{}_{}{}^{5}\right)^2}$$
(80)
Classically this action is equivalent to the action

$$\frac{1}{2}𝑑\tau \left\{\eta ^1\left((\dot{X})^2+\frac{(X_{}^{}{}_{}{}^{5})^2}{(X_{}^{}{}_{}{}^{4})^2+(X_{}^{}{}_{}{}^{5})^2}(\dot{X}^4)^2\right)\eta \varphi ^2\left((X_{}^{}{}_{}{}^{4})^2+(X_{}^{}{}_{}{}^{5})^2\right)\right\}.$$
(81)
where we have introduced an auxiliary tetrad field $`\eta `$ on the world-line. We may eliminate $`\dot{X}^4`$ in favour of the conserved conjugate momentum $`P=\frac{L}{\dot{X}^4}=\eta ^1\frac{(X_{}^{}{}_{}{}^{5})^2}{(X_{}^{}{}_{}{}^{4})^2+(X_{}^{}{}_{}{}^{5})^2}\dot{X}_4`$ by making a Legendre transformation with respect to $`X^4`$. We then get the action

$$𝑑\tau \left\{L(X^\mu ,\dot{X}^\mu ,P)P\dot{X}^4\right\}=\frac{1}{2}𝑑\tau \left\{\eta ^1(\dot{X})^2\eta \left(\frac{(X_{}^{}{}_{}{}^{4})^2+(X_{}^{}{}_{}{}^{5})^2}{(X_{}^{}{}_{}{}^{5})^2}P^2+\varphi ^2\left((X_{}^{}{}_{}{}^{4})^2+(X_{}^{}{}_{}{}^{5})^2\right)\right)\right\}$$
(82)
and this is classically equivalent to the particle action used in  if we identify $`P`$ with the rest mass of the dyon in the absence of the scalar field. We now see that we should define the four-dimensional scalar fields as $`\mathrm{\Phi }^a=R\varphi ^a`$ in order to get the action in . The sixth scalar should then be defined as $`\mathrm{\Phi }^6RB_{45}`$ in order to have the same mass-dimension on all the scalars. (The normalization constant can be determined from the (2,0) supersymmetry, and the requirement that $`\mathrm{\Phi }^𝒜`$ ($`𝒜=1,\mathrm{},6`$) should transform as a 6-vector under the R-symmetry $`SO(6)`$.) We now see that the term $`R^2B_{45}R\mathrm{\Phi }^6`$ becomes irrelevant compared to the other scalars, as $`R0`$.

Charge quantization is a quantum effect. In the classical limit the electric charge takes a continuum of values. From Dirac quantization of the charges it follows that if we take the limit $`\mathrm{}0`$ while we keep both the electric and magnetic charges $`nq`$ and $`mg`$ fixed, then we must take $`n,m\mathrm{}`$, and hence the charges will take a continuum of values. If we at the same time take the weak coupling limit $`q0`$ as we let $`\mathrm{}0`$ then we may $`m`$ fixed and finite. When taking any of these limits we may thus continuously rotate the strings relative to each other. Two strings of types $`(n,m)`$ and $`(n^{},m)`$ will correspond to two dyons with same magnetic charges, but with different electric charges.

We thus regain the action that was considered in . From this action the Taub-NUT moduli space was obtained there. The corresponding exact metric and moduli space was obtained in . Our six dimensional theory does not reproduce that exact moduli space in four dimensions. It would be very interesting if one could determine the exact moduli space metric of the six dimensional theory from knowing the asymptotic behaviour as well as the exact moduli space of the four dimensional theory.

Thus far almost everything we know about these theories have been deduced from string theory. But if we want to learn something new about string theory by studying these six-dimensional theories (which of course is our ultimate goal), then it seems to be more logical to see what can learned about them by using the fact that they should reduce to superYang-Mills theories in four dimensions.


text:
# The Dirichlet Problem on Quadratic Surfaces

## 1. Introduction

In this paper we present a fast, exact algorithm for solving Dirichlet problems with polynomial boundary functions on a quadratic surface in $`𝐑^n`$ $`(n2)`$. To illustrate the kind of Dirichlet problem we study, fix $`b=(b_1,\mathrm{},b_n)𝐑^n`$. For $`x=(x_1,\mathrm{},x_n)𝐑^n`$, we will write

$$bx^2=b_1^2x_1^2+\mathrm{}+b_n^2x_n^2.$$

Suppose we are given a polynomial $`p`$ on $`𝐑^n`$. We wish to find a harmonic polynomial that equals $`p`$ on the quadratic surface $`\{x𝐑^n:bx^2=1\}`$.

Even if all the $`b_j`$ are nonzero (so that our quadratic surface is bounded and is, in particular, an ellipsoid), computing a solution to this Dirichlet problem presents several difficulties. A standard means of expressing the solution to the Dirichlet problem for bounded domains involves the Green’s function and integration. However, the Green’s function of an ellipsoid does not have a known formula allowing for exact computations. An alternative approach avoids integration by employing a finite difference, finite element, or Galerkin-type scheme to approximate the solution, but again this procedure will not produce an exact solution.

If $`b_1=b_2=\mathrm{}=b_n=1`$, then our quadratic surface is the unit sphere. In this case, a fast algorithm for finding exact solutions is presented in . That algorithm involves differentiation but no integration. The basis of that algorithm is that any polynomial $`p`$ of degree $`m`$ on $`\text{R}^n`$ can be decomposed in the form

$$p=h+(x^21)f,$$

where $`h`$ is a harmonic polynomial of degree at most $`m`$ and $`f`$ is a polynomial of degree at most $`m2`$. Because $`h`$ is harmonic and equals $`p`$ on the unit sphere, it is the solution to our Dirichlet problem. The algorithm presented in  shows how the polynomials $`h`$ and $`f`$ in the decomposition above can be computed via differentiation from $`p`$. Unfortunately these techniques work only on spheres and so do not provide an algorithm for nonspherical ellipsoids or other quadratic surfaces.

In this paper we solve the Dirichlet problem discussed above, getting solutions for ellipsoids as well as for elliptic cylinders (for example, $`\{x𝐑^3:x_1^2+2x_2^2=1\}`$) and paraboloids (for example, $`\{x𝐑^3:x_3=x_1+x_2^2\}`$). We will begin by extending the decomposition above to a collection of quadratic surfaces. We then use this decomposition to produce a system of linear equations whose solution will give an exact solution to our Dirichlet problem. We will show how this system of linear equations has a structure allowing it to be reduced to smaller systems of linear equations, thus producing a fast algorithm. The algorithm requires differentiation of the boundary function, but no integration.

Before we turn to these matters, we need to present some background, much of which appears in . A multi-index is an $`n`$-tuple $`\alpha =(\alpha _1,\mathrm{},\alpha _n)`$ of nonnegative integers. The order of $`\alpha `$, denoted $`|\alpha |`$, is defined by

$$|\alpha |=\alpha _1+\mathrm{}+\alpha _n.$$

We let $`x^\alpha `$ denote the monomial $`x_{1}^{}{}_{}{}^{\alpha _1}\mathrm{}x_{n}^{}{}_{}{}^{\alpha _n}`$ and $`D^\alpha `$ denote the differential operator $`D_{1}^{}{}_{}{}^{\alpha _1}\mathrm{}D_{n}^{}{}_{}{}^{\alpha _n}`$, where $`D_j`$ denotes differentiation with respect to $`x_j`$. If $`q`$ is a polynomial on $`𝐑^n`$ given by $`q(x)=_\alpha c_\alpha x^\alpha `$, then $`q(D)`$ is the differential operator defined by $`q(D)=_\alpha c_\alpha D^\alpha `$. A polynomial is called homogeneous of degree $`m`$ if it is a linear combination of monomials of degree $`m`$.

Ernst Fischer proved that given a homogeneous polynomial $`q`$ on $`𝐑^n`$, every homogeneous polynomial $`p`$ of degree $`m`$ can be decomposed uniquely as $`p=h+qf`$, where $`h`$ is a homogeneous polynomial of degree $`m`$ satisfying $`q(D)h=0`$ and $`f`$ is a homogeneous polynomial of suitable degree.

In  the subject of more general decompositions is discussed. Given two polynomials $`g`$ and $`q`$ on $`𝐑^n`$, the relevant question is whether an arbitrary polynomial $`p`$ can be decomposed as $`p=h+qf`$, where $`h`$ and $`f`$ are polynomials, with $`h`$ satisfying $`g(D)h=0`$. Shapiro refers to a pair $`(g,q)`$ with this property as a generalized Fischer pair. He asks: Which $`(g,q)`$ form generalized Fischer pairs? Note that if $`g(x)=x^2`$, then $`g(D)`$ is the Laplacian and so this decomposition would require $`h`$ to be harmonic. We will provide examples of a robust class of quadratic polynomials $`q`$ that form a generalized Fischer pair with $`g(x)=x^2`$, and we give explicit examples of the decomposition via our algorithm.

Some of the surfaces that we consider are unbounded (for example, the elliptic cylinders and paraboloids mentioned above). Thus unique solutions to Dirichlet problems on these surfaces, even in the class of polynomials, are neither automatic nor expected. For example, the set of harmonic polynomials that vanish on the hyperplane $`\{x𝐑^n:x_n=0\}`$ is not trivial. However, we will show that for the quadratic surfaces we consider, polynomial solutions to the Dirichlet problem with polynomial boundary functions are unique.

The paper is organized as follows: In Section 2 we begin by presenting Fischer’s lemma and the corresponding decomposition theorem. We then extend these results to cover a wider class of generalized Fischer pairs. These generalized Fischer pairs are then used to solve the Dirichlet problem. We prove that the polynomial produced by this technique is the unique polynomial solving the Dirichlet problem, even when our quadratic surface is unbounded. In Section 3 we describe a fast algorithm for computing the solution to the Dirichlet problem promised by the results in Section 2. In Section 4 we present some examples, computed using our algorithm, of solutions to Dirichlet problems on ellipsoids. The Appendix contains a differentiation formula needed by our algorithm.

## 2. Fischer’s Lemma and the Dirichlet Problem

In this section, we state and prove our generalization of Fischer’s results. Then we prove the decomposition theorem that solves the Dirichlet problem. Then we show that even when our quadratic surfaces are unbounded, the solution given by our decomposition theorem is the only polynomial solution.

We begin by stating Fischer’s results, which are nicely restated and proved in ; also see . Fix an integer $`n2`$. We will always use $`m`$ to denote a nonnegative integer. Let $`𝒫_m`$ denote the vector space of polynomials (with real coefficients) of degree at most $`m`$ on $`𝐑^n`$. For convenience, we define $`𝒫_k`$ to be $`\{0\}`$ for $`k<0`$. As usual, $`\mathrm{\Delta }`$ denotes the Laplacian.

###### Lemma 2.1 (Fischer’s Lemma).

Suppose $`b=(b_1,\mathrm{},b_n)`$, where each $`b_j0`$. Define $`L:𝒫_m𝒫_m`$ by

$$L(f)=\mathrm{\Delta }\left((bx^21)f\right).$$

Then $`L`$ is a linear, degree-preserving, bijection of $`𝒫_m`$ onto itself.

Fischer’s Lemma leads to Fischer’s Decomposition Theorem, which gives a solution to the Dirichlet problem for ellipsoids.

###### Theorem 2.2 (Fischer’s Decomposition Theorem).

Suppose $`b=(b_1,\mathrm{},b_n)`$, where each $`b_j0`$. Let $`p𝒫_m`$. Then there exists a unique harmonic polynomial $`h𝒫_m`$ such that

$$p=h+(bx^21)f$$

for some $`f𝒫_{m2}`$.

Let $`E=\{x𝐑^n:bx^2<1\}`$. Then $`h`$ is the unique continuous function on $`\overline{E}`$ that is harmonic on $`E`$ and equals $`p`$ on the ellipsoid $`E`$.

We will need the following generalization of Corollary 5.3 of .

###### Lemma 2.3.

If $`b=(b_1,\mathrm{},b_n)`$, where each $`b_j0`$, then no nonzero polynomial multiple of $`bx^2`$ is harmonic.

###### Proof.

Suppose, to the contrary, that $`p`$ is a nonzero polynomial of degree $`m`$ such that $`bx^2p`$ is harmonic. Let $`E`$ be as in the theorem above, which states that there exists a harmonic polynomial $`h𝒫_m`$ such that $`h`$ equals $`p`$ on the ellipsoid $`E`$. But $`bx^2p`$ is also a harmonic polynomial that equals $`p`$ on $`E`$. Because $`bx^2p`$ has degree $`m+2`$, the polynomials $`bx^2p`$ and $`h`$ cannot be equal. But this contradicts the uniqueness of solutions to the Dirichlet problem on bounded domains. ∎

In order to generate new generalized Fischer pairs, let $`q`$ be a quadratic (degree $`2`$) polynomial on $`𝐑^n`$. We want to look at the map on $`𝒫_m`$ defined by

$$L(f)=\mathrm{\Delta }(qf).$$

Our goal is to identify choices of $`q`$ for which $`L`$ is a bijection of $`𝒫_m`$ onto $`𝒫_m`$. This leads us to the definition of a nonhyperbolic quadratic.

###### Definition 2.4.

A nonhyperbolic quadratic is a polynomial $`q`$ on $`𝐑^n`$ of the form

$$\underset{j=1}{\overset{n}{}}b_j^2x_j^2+\underset{j=1}{\overset{n}{}}c_jx_j+d,$$

where at least one $`b_j0`$.

Note that the next theorem, which gives the desired bijectivity, implies that no nonzero polynomial multiple of a nonhyperbolic quadratic is harmonic, generalizing Lemma 2.3. This result does not hold for arbitrary quadratic polynomials. In fact, even for a nonharmonic quadratic polynomial, a nonzero polynomial multiple might be harmonic. For example, $`x_1^23x_2^2`$ is not harmonic, but $`(x_1^23x_2^2)x_1`$ is harmonic.

###### Theorem 2.5.

Let $`q`$ be a nonhyperbolic quadratic. Define $`L:𝒫_m𝒫_m`$ by

$$L(f)=\mathrm{\Delta }(qf).$$

Then $`L`$ is a linear bijection of $`𝒫_m`$ onto $`𝒫_m`$.

###### Proof.

Clearly $`L`$ is a linear map of $`𝒫_m`$ into $`𝒫_m`$. Since $`𝒫_m`$ is finite dimensional, we need only show that $`L`$ is injective. So suppose that this is not true. Then there exists $`f𝒫_m`$, $`f0`$, such that $`L(f)=0`$. Without loss of generality, we may suppose that $`f`$ is of degree $`m`$ (otherwise, replace $`m`$ by a lower integer). We can write $`f=f_m+f_{m1}`$, where $`f_m`$ is homogeneous of degree $`m`$ and $`f_{m1}𝒫_{m1}`$. Write $`q=q_2+q_1`$, where $`q_2`$ is homogeneous degree $`2`$ (so $`q_2(x)=_{j=1}^nb_j^2x_j^2`$ in the notation above) and $`q_1𝒫_1`$. Because $`\mathrm{\Delta }(qf)=0`$, we know that

$$\mathrm{\Delta }(q_2f_m+q_2f_{m1}+q_1f)=0.$$

Thus

$$\mathrm{\Delta }(q_2f_m)=0,$$

because the other terms have lower degrees. We will now show that this implies that $`f_m=0`$, which is a contradiction.

Note that we have reduced our theorem to the case where $`q=q_2`$. Reordering the variables (if some of the $`b_j=0`$), we see that it suffices to prove our theorem in the case when $`q(x)=_{j=1}^rb_j^2x_j^2`$, where $`1rn`$ and $`b_1,\mathrm{},b_r`$ are all nonzero. To simplify notation, we will also replace $`f_m`$ in the previous paragraph with $`f`$. So again we have the assumption that $`\mathrm{\Delta }(qf)=0`$ and we want to prove that $`f=0`$, but now we have a special form for $`q`$.

If $`r=n`$, our desired conclusion that $`f=0`$ follows from Lemma 2.3. So suppose $`r<n`$. Let $`k`$ denote the degree of $`f`$ thought of as a polynomial in $`x_{r+1},\mathrm{},x_n`$ (temporarily think of $`x_1,\mathrm{},x_r`$ as constants to define $`k`$). Write

(2.6)
$$f=p+g,$$

where $`p`$ is the part of $`f`$ that is homogeneous of degree $`k`$ in the variables $`x_{r+1},\mathrm{},x_n`$ and $`g`$ is the remaining part of $`f`$, consisting of lower degree terms in the variables $`x_{r+1},\mathrm{},x_n`$.

Using the product formula for the Laplacian, which states that

(2.7)
$$\mathrm{\Delta }(qp)=p\mathrm{\Delta }q+q\mathrm{\Delta }p+2qp,$$

we obtain from (2.6) the equation

$$\mathrm{\Delta }(qf)=p\mathrm{\Delta }q+q\mathrm{\Delta }p+2qp+\mathrm{\Delta }(qg).$$

Let’s consider the degree, as a polynomial in $`x_{r+1},\mathrm{},x_n`$, of each term on the right side of this equation. Because $`\mathrm{\Delta }q`$ is a positive constant, this degree is $`k`$ for the first term. Because $`q`$ is independent of the variables $`x_{r+1},\mathrm{},x_n`$, the second and third terms have degree (as a function of $`x_{r+1},\mathrm{},x_n`$) less than $`k`$. Because the degree of $`g`$ (as a function of $`x_{r+1},\mathrm{},x_n`$) is less than $`k`$ and $`q`$ is independent of $`x_{r+1},\mathrm{},x_n`$, the degree of the fourth term is less than $`k`$.

Thus the only part of the right side of the equation above with degree $`k`$ (as a function of $`x_{r+1},\mathrm{},x_n`$) is the first term, $`p\mathrm{\Delta }q`$. The left side of the equation is $`0`$, so $`p\mathrm{\Delta }q=0`$. Hence $`p=0`$. But $`p`$ was the part of $`f`$ of highest degree in $`x_{r+1},\mathrm{},x_n`$. Hence $`f`$ is independent of $`x_{r+1},\mathrm{},x_n`$. Thus we can think of the equation $`\mathrm{\Delta }(qf)=0`$ as taking place in $`𝐑^r`$. Lemma 2.3 now implies that $`f=0`$, as desired. ∎

We now apply Theorem 2.5 to obtain the general decomposition theorem.

###### Theorem 2.8.

Suppose $`p𝒫_m`$. Let $`q`$ be a nonhyperbolic quadratic. Then there exists a unique harmonic polynomial $`h𝒫_m`$ such that

$$p=h+qf$$

for some $`f𝒫_{m2}`$.

###### Proof.

Note that $`\mathrm{\Delta }p𝒫_{m2}`$. Thus by Theorem 2.5, there exists $`f𝒫_{m2}`$ such that $`\mathrm{\Delta }(qf)=\mathrm{\Delta }p`$. Let $`h=pqf`$. Then $`h`$ is harmonic polynomial in $`𝒫_m`$ and $`p=h+qf`$, as desired.

To prove the uniqueness part of this theorem, suppose also that $`\stackrel{~}{h}`$ is a harmonic polynomial in $`𝒫_m`$ and that $`p=\stackrel{~}{h}+q\stackrel{~}{f}`$ for some $`\stackrel{~}{f}𝒫_{m2}`$. Then

$$h\stackrel{~}{h}=q(\stackrel{~}{f}f).$$

The left side of the equation above is harmonic, and hence $`\mathrm{\Delta }\left(q(\stackrel{~}{f}f)\right)=0`$. Theorem 2.5 now implies that $`\stackrel{~}{f}f=0`$, which implies that $`\stackrel{~}{h}=h`$, as desired. ∎

If $`p𝒫_m`$ and $`q`$ is a nonhyperbolic quadratic, we can consider the following Dirichlet problem: find a harmonic polynomial $`h𝒫_m`$ such that $`h`$ equals $`p`$ on the set $`\{x𝐑^n:q(x)=0\}`$. Clearly the $`h`$ produced by the theorem above solves this Dirichlet problem. Of course, $`\{x𝐑^n:q(x)=0\}`$ could be the empty set or a single point. Because $`q(x)\mathrm{}`$ as $`|x|\mathrm{}`$, the existence of a point in $`𝐑^n`$ where $`q`$ is negative is a convenient condition to ensure that $`\{x𝐑^n:q(x)=0\}`$ is a nondegenerate quadratic surface. For example, using the notation of Definition 2.4, $`\{x𝐑^n:q(x)=0\}`$ will be a nondegenerate quadratic surface if

$$d<\underset{\{j:b_j0\}}{}\frac{c_j^2}{4b_j^2}$$

or if $`c_j0`$ for some $`j`$ with $`b_j=0`$.

We now turn to the question of whether the polynomial $`h`$ produced by Theorem 2.8 is the unique polynomial solution to the Dirichlet problem discussed in the paragraph above. The following lemma will help us answer this uniqueness question. If we were working in $`𝐂^n`$ instead of $`𝐑^n`$, then Hilbert’s Nullstellensatz could be used to provide information about when a polynomial $`h`$ vanishing on the zero set of another polynomial $`q`$ is a polynomial multiple of $`q`$. A theorem called the Real Nullstellensatz (see, for example, , Chapter 3, Theorem 3.3) provides some information about polynomials on $`𝐑^n`$ vanishing on the zero set of another polynomial. However, we do not see how the Real Nullstellensatz can be used to prove the lemma below, so we have provided a proof without using such machinery.

###### Lemma 2.9.

Suppose $`q`$ is a nonhyperbolic quadratic that is negative at some point of $`𝐑^n`$. If $`h`$ is a polynomial on $`𝐑^n`$ such that $`h(x)=0`$ whenever $`q(x)=0`$, then $`h`$ is a polynomial multiple of $`q`$.

###### Proof.

We will prove this lemma by induction on the dimension $`n`$.

To get started, suppose $`n=1`$ and that $`q(x)=b^2x^2+cx+d`$, where $`b0`$, is a nonhyperbolic quadratic that is negative at some point of $`𝐑`$. Suppose $`h`$ is a polynomial on $`𝐑`$ such that $`h(x)=0`$ whenever $`q(x)=0`$. Because $`q`$ is negative at some point of $`𝐑`$ and $`q(x)\mathrm{}`$ as $`x\mathrm{}`$, we see that $`q`$ has precisely two distinct zeros. The polynomial $`h`$ vanishes on both these zeros, and thus $`h`$ is a polynomial multiple of the quadratic polynomial $`q`$, as desired.

Now suppose that the lemma holds in dimension $`n1`$. Let $`q`$ be a nonhyperbolic quadratic that is negative at some point of $`𝐑^n`$. Relabelling coordinates, if necessary, we can assume that

$$q(x)=\underset{j=1}{\overset{n}{}}b_j^2x_j^2+\underset{j=1}{\overset{n}{}}c_jx_j+d,$$

where $`b_j0`$ for some $`j\{1,\mathrm{},n1\}`$. Let $`y`$ denote a typical point of $`𝐑^{n1}`$ and let $`z`$ denote a typical point of $`𝐑`$; thus $`(y,z)`$ denotes a typical point of $`𝐑^n`$.

Suppose $`h`$ is a polynomial on $`𝐑^n`$ of degree $`m`$ such that $`h(x)=0`$ whenever $`q(x)=0`$. We need to show that $`h`$ is a polynomial multiple of $`q`$.

For $`z𝐑`$ such that $`q(y,z)`$ is negative for some $`y𝐑^{n1}`$, define a nonhyperbolic quadratic $`q_z`$ on $`𝐑^{n1}`$ by

$$q_z(y)=q(y,z),$$

and define a polynomial $`h_z`$ on $`𝐑^{n1}`$ of degree at most $`m`$ by

$$h_z(y)=h(y,z).$$

Then $`h_z(y)=0`$ whenever $`q_z(y)=0`$, and thus by our induction hypothesis there is a polynomial $`f_z`$ on $`𝐑^{n1}`$ such that

(2.10)
$$h_z(y)=f_z(y)q_z(y)$$

for all $`y𝐑^{n1}`$. Clearly $`f_z`$ has degree at most $`m2`$.

For $`y𝐑^{n1}`$ such that $`q(y,z)`$ is negative for some $`z𝐑`$, define a polynomial $`q^y`$ on $`𝐑`$ by

$$q^y(z)=q(y,z),$$

and define a polynomial $`h^y`$ on $`𝐑`$ of degree at most $`m`$ by

$$h^y(z)=h(y,z).$$

Then $`h^y(z)=0`$ whenever $`q^y(z)=0`$. We claim that there is a polynomial $`g^y`$ on $`𝐑`$ such that

(2.11)
$$h^y(z)=g^y(z)q^y(z)$$

for all $`z𝐑`$. If $`b_n0`$, then $`q^y`$ is a nonhyperbolic quadratic on $`𝐑`$ and the claim follows from the dimension $`1`$ case that was proved at the beginning of this proof. If $`b_n=0`$ but $`c_n0`$, then $`q^y`$ is a polynomial on $`𝐑`$ of degree 1 and the claim follows easily. Finally, if $`b_n=0`$ and $`c_n=0`$, then $`q^y`$ is a negative constant, in which case the claim is trivially true. In any case, we see that $`g^y`$ has degree at most $`m`$.

Let $`\mathrm{\Omega }=\{(y,z)𝐑^n:q(y,z)<0\}`$. Combining (2.10) and (2.11), we see that

$$h(y,z)=f_z(y)q(y,z)=g^y(z)q(y,z)$$

for all $`(y,z)\mathrm{\Omega }`$. Thus we can define a function $`p`$ on $`\mathrm{\Omega }`$ by

$$p(y,z)=f_z(y)=g^y(z).$$

Because

$$p(y,z)=\frac{h(y,z)}{q(y,z)}$$

for all $`(y,z)\mathrm{\Omega }`$, we see that $`p`$ is real-analytic on $`\mathrm{\Omega }`$.

Suppose $`\alpha =(\alpha _1,\mathrm{},\alpha _n)`$ is a multi-index of order greater than $`2m`$. Then either

$$\alpha _1+\mathrm{}+\alpha _{n1}>m1\text{or}\alpha _n>m+1.$$

This implies that $`D^\alpha p(y,z)=0`$ for all $`(y,z)\mathrm{\Omega }`$. Because all sufficiently high-order partial derivatives of $`p`$ equal $`0`$ on the open set $`\mathrm{\Omega }`$, we conclude that $`p`$ is a polynomial on $`\mathrm{\Omega }`$. Hence we can think of $`p`$ as a polynomial defined everywhere on $`𝐑^n`$.

Finally, because $`h(y,z)=p(y,z)q(y,z)`$ for all $`(y,z)\mathrm{\Omega }`$, and because polynomials that agree on a nonempty open subset of $`𝐑^n`$ must agree everywhere, we have $`h=pq`$. Thus $`h`$ is a polynomial multiple of $`q`$, as desired. ∎

Now we can combine the previous lemma and Theorem 2.5 to prove the desired uniqueness result. Of course in the ellipsoidal case (where each $`b_j0`$, in the notation of Definition 2.4) uniqueness follows easily from the boundedness of the surface, but we want to consider also elliptic cylinders and paraboloids. Note that the uniqueness result in the theorem below fails on some nondegenerate quadratic surfaces, so the hypothesis that $`q`$ is nonhyperbolic cannot be deleted. For example, on the quadratic surface defined by $`\{x𝐑^n:x_1^23x_2^21=0\}`$, any solution to a Dirichlet problem can be used to produce another solution by adding to it the harmonic polynomial $`(x_1^23x_2^21)x_1`$, which vanishes on the quadratic surface in question.

To obtain uniqueness results on half-spaces, even in the class of polynomial solutions, a growth condition on the solutions is needed (see ). However, the theorem below shows that we have unique polynomial solutions on our quadratic surfaces without the requirement of a growth condition.

###### Theorem 2.12.

Suppose $`q`$ is a nonhyperbolic quadratic that is negative at some point of $`𝐑^n`$. If $`p𝒫_m`$, then there is a unique harmonic polynomial $`h`$ that equals $`p`$ on $`\{x𝐑^n:q(x)=0\}`$. Furthermore,

$$h=pqf$$

for some $`f𝒫_{m2}`$.

###### Proof.

Take $`h`$ and $`f`$ as in Theorem 2.8. Then $`h`$ is a harmonic polynomial that equals $`p`$ on $`\{x𝐑^n:q(x)=0\}`$; furthermore $`f𝒫_{m2}`$ and $`h=pqf`$.

To prove uniqueness, suppose $`\stackrel{~}{h}`$ is also a harmonic polynomial that equals $`p`$ on $`\{x𝐑^n:q(x)=0\}`$. Then $`\stackrel{~}{h}h`$ equals $`0`$ on $`\{x𝐑^n:q(x)=0\}`$. By Lemma 2.9, $`\stackrel{~}{h}h`$ is a polynomial multiple of $`q`$. However, Theorem 2.5 implies that no nonzero polynomial multiple of $`q`$ is harmonic. Thus $`\stackrel{~}{h}h=0`$ and hence $`\stackrel{~}{h}=h`$, completing the proof of uniqueness. ∎

## 3. Algorithm

In this section, we will turn the results of the previous section into a computationally useful algorithm for solving the Dirichlet problem. Every polynomial is a sum of homogeneous polynomials, and thus to solve Dirichlet problems with polynomial boundary functions it suffices to do so for homogeneous polynomials.

Let $`_m`$ denote the vector space of polynomials on $`𝐑^n`$ that are homogeneuos of degree $`m`$ (the polynomial $`0`$ is homogeneous of every degree, so $`0_m`$). Suppose $`q`$ is a nonhyperbolic quadratic and $`p_{m+2}`$. We use the decomposition in Theorem 2.8 to write $`p=h+qf`$, where $`h`$ is a harmonic polynomial in $`𝒫_{m+2}`$ and $`f𝒫_m`$. Breaking each polynomial into its homogeneous components, we obtain

$$p=\underset{j=0}{\overset{m+2}{}}h_j+(q_2+q_1+q_0)\underset{j=0}{\overset{m}{}}f_j,$$

where each $`h_j`$ is harmonic. We can break the equation above into homogeneuous equations by degree to obtain the system

$`p`$ $`=h_{m+2}+q_2f_m,`$
$`q_1f_m`$ $`=h_{m+1}+q_2f_{m1},`$
$`q_0f_mq_1f_{m1}`$ $`=h_m+q_2f_{m2},`$
(3.1) $`\mathrm{}`$
$`q_0f_2q_1f_1`$ $`=h_2+q_2f_0,`$
$`q_0f_1q_1f_0`$ $`=h_1,`$
$`q_0f_0`$ $`=h_0.`$
Here $`p`$ and $`q_2,q_1,q_0`$ are known and we need to compute $`h_{m+2},\mathrm{},h_0`$ and $`f_m,\mathrm{},f_0`$. Our plan of attack is to use the first equation to find $`f_m`$, which will then give us both $`h_{m+2}`$ and the left side of the next equation. We can repeat the procedure with the second equation to find $`f_{m1}`$, which will give us both $`h_{m+1}`$ and the left side of the next equation. We continue this process until we have found $`h_j`$ for each $`j`$. Then $`h=_{j=0}^{m+2}h_j`$ is the harmonic function that agrees with $`p`$ on $`\{x𝐑^n:q(x)=0\}`$.

Thus, we turn our attention to equations of the form

(3.2)
$$p=h+q_2f$$

where $`p`$ is a known homogeneous polynomial, $`q_2`$ is the highest degree part of a known nonhyperbolic quadratic, $`h`$ is an unknown harmonic polynomial homogeneous of degree $`\mathrm{deg}(p)`$, and $`f`$ is an unknown polynomial homogeneous of degree $`\mathrm{deg}(p)2`$. Note that all the equations in the system above are of this form (except the last two equations, which are trivial to solve for $`h_1`$ and $`h_0`$). Hence an algorithm for finding the solution to (3.2) will give us an algorithm for solving our Dirichlet problem.

To solve equation (3.2), eliminate $`h`$ by taking the Laplacian of both sides, getting

(3.3)
$$\mathrm{\Delta }p=\mathrm{\Delta }(q_2f).$$

Now our problem has been reduced to finding $`f_m`$ satisfying the equation above, where $`p_{m+2}`$ and $`q_2(x)=b_1^2x_1^2+\mathrm{}+b_n^2x_n^2`$ are known, with at least one $`b_j0`$. Note that Theorem 2.5 implies that (3.3) has a unique solution $`f_m`$.

The algorithm we introduce for finding $`f`$ involves repeated differentiation of (3.3). Because $`f_m`$, we will know $`f`$ once we know the constants $`D^\alpha f`$ for every multi-index $`\alpha `$ of order $`m`$. We need one more piece of notation before finding the constants $`D^\alpha f`$. For $`1jn`$, let $`e_j`$ denote the multi-index whose $`j^{\text{th}}`$-coordinate equals $`1`$ and whose other coordinates equal $`0`$.

For $`\alpha `$ a multi-index of order $`m`$, apply the differential operator $`D^\alpha `$ to both sides of (3.3), getting

(3.4) $`D^\alpha (\mathrm{\Delta }p)`$ $`=D^\alpha \left(\mathrm{\Delta }(q_2f)\right)`$
$`=\mathrm{\Delta }\left(D^\alpha (q_2f)\right)`$
$`=\mathrm{\Delta }\left(q_2D^\alpha f+{\displaystyle \underset{j=1}{\overset{n}{}}}2\alpha _jb_j^2x_jD^{\alpha e_j}f+{\displaystyle \underset{j=1}{\overset{n}{}}}\alpha _j(\alpha _j1)b_j^2D^{\alpha 2e_j}f\right),`$
where the last equality comes from Proposition A.2 in our Appendix and the explicit form of $`q_2`$. Because $`D^\alpha f`$ is a constant, the Laplacian of the first term in parentheses above equals $`2b^2D^\alpha f`$. Because $`D^{\alpha e_j}f`$ has degree 1, the Laplacian of the second term in parentheses above can be easily computed using the product formula for the Laplacian (2.7), and the last equation becomes

(3.5)
$$D^\alpha (\mathrm{\Delta }p)=\left(2b^2+4\underset{j=1}{\overset{n}{}}\alpha _jb_j^2\right)D^\alpha f+\underset{j=1}{\overset{n}{}}\left[\alpha _j(\alpha _j1)b_j^2\underset{k=1}{\overset{n}{}}(D^{\alpha +2e_k2e_j}f)\right].$$

Note that multi-index $`\alpha +2e_k2e_j`$ in the equation above has order $`m`$. Thus as $`\alpha `$ ranges over all multi-indices of order $`m`$, (3.5) gives us a system of linear equations in the unknowns $`\{D^\alpha f\}`$. This system of equations can be solved using Gaussian elimination, giving us $`f`$ and thus solving our Dirichlet problem.

Any solution to the system of equations (3.5) in the unknowns $`\{D^\alpha f\}`$ gives a function $`f_m`$ satisfying (3.4) for all multi-indices $`\alpha `$ of order $`m`$, which implies that $`f`$ satisfies $`\mathrm{\Delta }p=\mathrm{\Delta }(q_2f)`$. However, we already know from Theorem 2.5 that there is a unique $`f_m`$ satisfying $`\mathrm{\Delta }p=\mathrm{\Delta }(q_2f)`$. Thus the system of equations (3.5) in the unknowns $`\{D^\alpha f\}`$ has a unique solution.

The number of operations needed to compute the solution to a system of $`s`$ linear equations in $`s`$ unknowns is on the order of $`(2/3)s^3`$. Thus we can expect that solving the system (3.5) will take on the order of

(3.6)
$$\frac{2(m+n1)!^3}{3m!^3(n1)!^3}$$

operations (the formula for the number of multi-indices of order $`m`$ can be found, for example, on page 78 of ). However, a careful look at the system (3.5) shows that we can do much better.

For a fixed multi-index $`\alpha `$, the multi-indices that appear in (3.5) all have the form $`\alpha +2e_k2e_j`$. This leads us to define an equivalence relation on the set of multi-indices of order $`m`$ by declaring that two multi-indices $`\alpha `$ and $`\beta `$ are equivalent if $`\alpha _i\beta _i(mod2)`$ for each $`i=1,\mathrm{},n`$. This equivalence relation breaks the set of multi-indices of order $`m`$ into equivalence classes, and the system of equations (3.5) breaks into corresponding systems of equations. Hence instead of solving one large system of equations, we can solve several smaller systems of equations, which leads to considerable computational savings, as we will soon see.

If $`m`$ is larger than $`n`$, then we have $`2^{n1}`$ equivalence classes, corresponding to a choice of even or odd entry in the first $`n1`$ coordinates of a multi-index (the parity of the last coordinate is forced by the condition that the coordinates add up to $`m`$). As $`m`$ increases, the number of elements in each of these $`2^{n1}`$ equivalence classes divided by the total number of multi-indices of order $`m`$ approaches $`1/2^{n1}`$. Thus instead of solving one large system of $`\left(\genfrac{}{}{0pt}{}{m+n1}{m}\right)`$ equations in the same number of variables, we can solve $`2^{n1}`$ systems of equations, each system containing on the order of $`\left(\genfrac{}{}{0pt}{}{m+n1}{m}\right)/2^{n1}`$ equations and variables. This computation requires on the order of

(3.7)
$$2^{n1}\frac{2(m+n1)!^3}{3(2^{n1})^3m!^3(n1)!^3}$$

operations.

The ratio of (3.6) to (3.7) is $`2^{2n2}`$. Thus breaking our system of equations into smaller systems of equations using our equivalence relation reduces the number of operations by a factor of $`2^{2n2}`$. For example, if $`n=6`$, then this technique should reduce computation time by a factor of over 1000. This savings is needed even for moderate size $`m`$ and $`n`$, because the system of equations $`(\text{3.5})`$ grows large rapidly.

Usually we can do even better than reducing computations by a factor of $`2^{2n2}`$. Suppose, for example, that $`p(x)=x_1^{20}x_2^7`$. Then $`\mathrm{\Delta }p=42x_1^{20}x_2^5+380x_1^{18}x_2^7`$. Note that the left side of (3.5) equals 0 for all multi-indices of order 25 except $`(20,5,0,\mathrm{},0)`$ and $`(18,7,0,\mathrm{},0)`$. In other words, in only one of our equivalence classes (the equivalence class consisting of those multi-indices of order 25 whose first coordinate is even, second coordinate is odd, and coordinates 3 through $`n`$ are even) is the left side of (3.5) anything other than 0. When the left side of every equation in an equivalence class is 0, there is no need to perform Gaussian elimination to solve the system of equations in that equivalence class, because obviously all the unknowns equal 0 (recall that the system (3.5) has a unique solution). Thus in the example at hand, instead of solving $`2^{n1}`$ smaller systems of equations, we need only solve one smaller system of equations.

As can be seen from the reasoning in the previous paragraph, if $`p`$ is a monomial then our computation time is reduced, through the use of our equivalence classes, by another factor of $`2^{n1}`$. Thus if $`p`$ is a monomial, then computation time through the use of equivalence classes is reduced by a factor of $`2^{3n3}`$. For example, if $`n=6`$, then this technique should reduce computation time by a factor of over 32,000.

The full strength of the reduction discussed above by a factor of $`2^{3n3}`$, as opposed to a reduction by a factor of $`2^{2n2}`$, holds only for ellipsoids and elliptic cylinders. For ellipsoids and elliptic cylinders we can assume (perhaps after a translation) that each $`c_j`$ in Definition 2.4 equals $`0`$, which gives $`q_1=0`$ in the system (3.1). However, if $`q_10`$ then multiplication by $`q_1`$, when solving the second and successive equations in the system (3.1), can lead to nonzero left sides in equations in the system (3.5) other than the equations in the equivalence class corresponding to the exponents in the monomial $`p`$.

## 4. Examples

The algorithm outlined in the previous section has been implemented by the authors in Mathematica and in MATLAB. The Mathematica version produces exact solutions, in considerably less time than we believed possible when we started this project. Even with simple boundary polynomials in low dimensions, solutions to Dirichlet problems on nonhyperbolic quadratic surfaces tend to involve fractions with large numerators and denominators, as we will see in the examples presented below. Our Mathematica software can use floating point arithmetic to produce only decimal approximations to the exact solutions, with even faster times than when working with exact rational arithmetic. Our MATLAB implementation of the algorithm works only in floating point arithmetic, again quickly producing decimal approximations to the exact solutions.

Our Mathematica implementation of the algorithm is available on the first author’s web site; our MATLAB implementation of the algorithm is available on the third author’s web site. The appropriate web addresses are listed at the end of this paper; within those web sites look for information about this paper to find the software. Although our software is available without charge, its use requires Mathematica or MATLAB. The examples presented below were generated by our Mathematica implementation of the algorithm.

Our first example will be in dimension $`3`$ with a boundary function of degree $`7`$. Here, as elsewhere throughout this section, all fractions are given in reduced form.

###### Example 4.1.

Suppose $`p(x_1,x_2,x_3)=x_1^4x_2^3`$ and

$$q(x_1,x_2,x_3)=2x_1^2+3x_2^2+4x_3^21.$$

Then the following function is harmonic on $`𝐑^3`$ and agrees with $`p`$ on the ellipsoid $`\{x𝐑^3:q(x)=0\}`$:

$$\begin{array}{cc}\hfill x_1^4& x_2^3+(2x_1^2+3x_2^2+4x_3^21)(\frac{97950}{20144813}x_2^5\frac{2524856930}{100139865423}x_1^2x_2\hfill \\ & \frac{3423451}{60434439}x_1^4x_2\frac{148091}{33379955141}x_2^3\frac{2306686}{20144813}x_1^2x_2^3\frac{701980831}{500699327115}x_2\hfill \\ & +\frac{32326712}{7703066571}x_2x_3^2+\frac{3712712}{60434439}x_1^2x_2x_3^2+\frac{53836}{20144813}x_2^3x_3^2\frac{236464}{60434439}x_2x_3^4).\hfill \end{array}$$

The function above obviously equals $`x_1^4x_2^3`$ on the ellipsoid in question. Thus to verify that the function above is indeed the solution to our Dirichlet problem, it is only necessary to verify that the Laplacian of the function above equals $`0`$; Mathematica can easily perform this calculation.

Note that, as expected from our discussion in the previous section, every exponent of each $`x_j`$ in the solution above has the same parity as in the boundary function. Specifically, in the solution above the exponents of $`x_1`$, $`x_3`$, and $`x_4`$ are even and the exponent of $`x_2`$ is odd, thus following the pattern of the boundary function $`x_1^4x_2^3`$.

Our Mathematica implementation of the algorithm can handle quadratic surfaces defined with symbols as well as concrete numbers. This capability is illustrated in our second example, which takes place in dimension $`4`$ with a boundary function of degree $`7`$.

###### Example 4.2.

Suppose $`p(x_1,x_2,x_3,x_4)=x_1^3x_2^2x_3x_4`$ and

$$q(x_1,x_2,x_3,x_4)=cx_1^2+3x_2^2+4x_3^2+5x_4^21.$$

Then the following function is harmonic on $`𝐑^4`$ and agrees with $`p`$ on the ellipsoid $`\{x𝐑^4:q(x)=0\}`$:

$$\begin{array}{cc}\hfill x_1^3& x_2^2x_3x_4\hfill \\ & +\frac{cx_1^2+3x_2^2+4x_3^2+5x_4^21}{788400+367920c+48712c^2+2520c^3+45c^4}(4(50+3c)(36+5c)x_1x_3^3x_4\hfill \\ & 12(2190+281c+9c^2)x_1^3x_3x_4\frac{(82800+21868c+1764c^2+45c^3)}{3(10+c)}x_1x_3x_4\hfill \\ & (50400+21118c+1845c^2+45c^3)x_1x_2^2x_3x_4+5(46+3c)(36+5c)x_1x_3x_4^3).\hfill \end{array}$$

We now return to dimension $`3`$ but increase the degree of the boundary function to $`10`$.

###### Example 4.3.

The harmonic polynomial that equals $`x_1^{10}`$ on the ellipsoid

$$\{x𝐑^3:2x_1^2+3x_2^2+4x_3^21=0\}$$

has value

$$\frac{500945213823452554440546462385400584789}{397263369506735959801289842040922215251461}$$

at the origin.

Of course, the value at the origin of the solution to this Dirichlet problem is simply the constant term in the polynomial that gives the solution. In the example above, we have given only the constant term because displaying the entire solution would require more than a page.

The example above illustrates the difference between the sphere and ellipsoids. The solution to the Dirichlet problem on the unit sphere in $`𝐑^3`$ with boundary function $`x_1^{10}`$ has value $`1/11`$ at the origin, and the largest integer appearing in the numerator or denominator of the coefficients of any term of the solution is $`46189`$. In contrast, if the ellipsoid $`\{x𝐑^3:2x_1^2+3x_2^2+4x_3^21=0\}`$ replaces the unit sphere, then the value of the solution at the origin is given by the fraction in the example above, and the coefficients of the other terms of the solution have similarly huge numerators and denominators.

We conclude this section by giving some data about the speed of our algorithm. The times given below are CPU times used by the Mathematica implementation of our software on a Windows desktop computer with a 1.7 GHz Intel Pentium 4 chip.

Examples 4.1 and 4.2 each took about 0.04 seconds. To produce Example 4.3, we found the solution to the Dirichlet problem with boundary function $`x_1^{10}`$ on the ellipsoid $`\{x𝐑^3:2x_1^2+3x_2^2+4x_3^21=0\}`$; this took less than 0.2 seconds. On the same ellipsoid, solving the Dirichlet problem with boundary function $`x_1^{20}`$ takes 55 seconds. Changing the boundary function to $`x_1^{25}`$ increases the solution time to 10 minutes, 5 seconds. Finally, changing the boundary function to $`x_1^{30}`$ (still on the same ellipsoid) increases the solution time to 2 hours, 29 minutes, 32 seconds.

The long solution times reported in the paragraph above arise because of the huge overhead associated with exact manipulation of the gigantic numerators and denominators that appear in the solutions to these high-degree problems. For example, most of the numerators and denominators are larger than $`10^{200}`$ for the Dirichlet problem above with boundary value $`x_1^{20}`$, rising to over $`10^{450}`$ for the Dirichlet problem above with boundary value $`x_1^{25}`$.

However, our algorithm is fast even for high-degree boundary functions if we use Mathematica’s floating point arithmetic instead of exact rational arithmetic. For example, if we ask the Mathematica implementation of our algorithm to use floating point arithmetic to solve the Dirichlet problem mentioned above with boundary function $`x_1^{30}`$, then the computation time is reduced from almost 2.5 hours to find the exact solution to just 8.9 seconds to find a very good decimal approximation of the solution.

## Appendix

This purpose of this appendix is to prove the differentiation formula given by Proposition A.2, which was used in the derivation of our algorithm in Section 3.

We begin with following lemma, which gives a formula for differentiating the product of a linear function and a polynomial.

###### Lemma A.1.

Suppose $`f`$ is a polynomial on $`𝐑^n`$ and $`\alpha `$ is a multi-index. Then

$$D^\alpha (gf)=gD^\alpha f+\underset{j=1}{\overset{n}{}}\alpha _j(D_jg)(D^{\alpha e_j}f)$$

for every $`g𝒫_1`$.

###### Proof.

Fix $`g𝒫_1`$. We will prove our desired result by induction on $`|\alpha |`$. To get started, note that the desired result is obviously true when $`|\alpha |=0`$.

Now suppose that $`|\alpha |>0`$ and that the desired result holds for all multi-indices of smaller order. Choose $`k`$ such that $`\alpha _k>0`$. Then

$`D^\alpha (gf)`$ $`=D^{\alpha e_k}D_k(gf)`$
$`=D^{\alpha e_k}\left(gD_kf+(D_kg)f\right)`$
$`=D^{\alpha e_k}(gD_kf)+(D_kg)(D^{\alpha e_k}f),`$
where the last equation holds because $`D_kg`$ is constant. Applying our induction hypothesis to evaluate $`D^{\alpha e_k}(gD_kf)`$ now gives
$`D^\alpha (gf)`$ $`=gD^\alpha f+{\displaystyle \underset{j=1}{\overset{n}{}}}\alpha _j(D_jg)(D^{\alpha e_j}f),`$
completing the proof. ∎

The next proposition, which gives a formula for differentiating the product of a quadratic polynomial and another polynomial, was used in deriving the system of equations (3.5).

###### Proposition A.2.

Suppose $`f`$ is a polynomial on $`𝐑^n`$ and $`\alpha `$ is a multi-index. Then

$$D^\alpha (qf)=qD^\alpha f+\underset{j=1}{\overset{n}{}}\alpha _j(D_jq)(D^{\alpha e_j}f)+\frac{1}{2}\underset{j=1}{\overset{n}{}}\alpha _j(\alpha _j1)(D_j^2q)(D^{\alpha 2e_j}f)$$

for every nonhyperbolic quadratic $`q`$.

###### Proof.

Fix a nonhyperbolic quadratic $`q`$. We will prove our desired result by induction on $`|\alpha |`$. To get started, note that the desired result is obviously true when $`|\alpha |=0`$.

Now suppose that $`|\alpha |>0`$ and that the desired result holds for all multi-indices of smaller order. Choose $`k`$ such that $`\alpha _k>0`$. Then

$`D^\alpha (qf)`$ $`=D^{\alpha e_k}D_k(qf)`$
(A.3) $`=D^{\alpha e_k}(qD_kf)+D^{\alpha e_k}\left((D_kq)f\right).`$
Using Lemma A.1 to evaluate the last term (applicable because $`D_kq𝒫_1`$) gives

(A.4)
$$D^{\alpha e_k}\left((D_kq)f\right)=(D_kq)(D^{\alpha e_k}f)+(\alpha _k1)(D_k^2q)(D^{\alpha 2e_k}f),$$

where we have used the fact that $`D_jD_kq=0`$ whenever $`jk`$. Using our induction hypothesis to evaluate the first term on the right side of (A.3) and using (A.4) to evaluate the second term on the right side of (A.3) now gives

$$D^\alpha (qf)=qD^\alpha f+\underset{j=1}{\overset{n}{}}\alpha _j(D_jq)(D^{\alpha e_j}f)+\frac{1}{2}\underset{j=1}{\overset{n}{}}\alpha _j(\alpha _j1)(D_j^2q)(D^{\alpha 2e_j}f),$$

completing the proof. ∎