A new approach to derive transparent boundary conditions

(TBCs) for wave, Schr¨odinger and drift-diffusion

equations is presented. It relies on the pole condition approach and distinguishes physical reasonable and unreasonable solutions by the location of the singularities of the spatial Laplace transform U of the exterior solution. By the condition that U is analytic in some region TBCs are established. To realize the pole condition numerically, a Möbius transform is used to map the region of analyticity to the unit disc. There the Laplace transform is expanded in a power series. The equations coupling the coefficients of the power series with the interior provide the TBC. Numerical result for the damped wave equation show that the error introduced by truncating the power series decays exponentially in the number of coefficients

Ruprecht, D.

Schaedle, A.

Schmidt, F.

English

Repository: Freie Universität Berlin (FU), Math Department (fu_mi_publications)

TRANSPARENT BOUNDARY CONDITIONS – THE POLE CONDITION APPROACHA. Scha¨dle◦,∗, D. Ruprecht◦, F. Schmidt◦◦Zuse Institute, Berlin, Germany.∗Email: schaedle@zib.deAbstractA new approach to derive transparent boundary con-ditions (TBCs) for wave, Schro¨dinger and drift-diffusionequations is presented. It relies on the pole condition ap-proach and distinguishes physical reasonable and unrea-sonable solutions by the location of the singularities of thespatial Laplace transform U of the exterior solution. Bythe condition that U is analytic in some region TBCs areestablished. To realize the pole condition numerically, aMo¨bius transform is used to map the region of analyticityto the unit disc. There the Laplace transform is expandedin a power series. The equations coupling the coefficientsof the power series with the interior provide the TBC. Nu-merical result for the damped wave equation show that theerror introduced by truncating the power series decays ex-ponentially in the number of coefficients.IntroductionTransparent boundary conditions are a key ingredientfor the simulation of wave propagation on unbounded do-mains. In this talk work in progress is presented.Prototypes of the governing equations under consider-ation are the wave, drift-diffusion and Schro¨dinger equa-tions on the real line for t > 0 given by∂ttu = ∂xxu− k2u, (1)∂tu = ∂xxu + 2d∂xu, (2)i∂tu = ∂xxu− k2u. (3)All of these have to complemented by appropriate initialvalues. To treat (1) - (3) simultaneously the symbol p(∂t)is introduced. Hence the generic equation isp(∂t)u = ∂xxu + 2d∂xu− k2u. (4)For the procedure to derive exact non-local transparentboundary conditions we refer to the recent review arti-cles [1]. The pole condition approach is an alternativeand as we hope to show a more flexible way of deriv-ing transparent boundary conditions. Almost immedi-ately the pole condition approach yields an algorithm toimplement approximate local transparent boundary con-ditions. The pole condition for time-harmonic problemsis studied in [2], where it is shown that it coincides withthe Sommerfeld radiation condition.Alternative derivation of TBCsSuppose we are only interested in the solution u re-stricted to the interval [−a, a]. Furthermore suppose thatthe initial value(s) are compactly supported in [−a, a]. Totruncate the computational domain TBCs are needed. Theexact TBCs are in general convolution in time, i.e. theyare non-local.Variational formulationMultiplying (4) by a test function and integrating overthe real line yields∫p(∂t)uv dx =∫−∂xu∂xv + 2d∂xuv − k2uv dx. (5)As test functions we chose v(x) = e−s(x−a) for x > aand v(x) = es(x+a) for x < −a with a complex parame-ter s with <s > 0. The integral over the real line is splitinto three parts: an integral from −∞ to −a, from −a toa and from a to ∞. DefiningU (r)(t, s) :=∫∞0u(t, x + a)e−sx,which is the Laplace transform of the solution u in theright exterior, and similar U (l) one obtains after somesimple manipulations∫ a−ap(∂t)uv dx + p(∂t)U(r) + p(∂t)U(l) =∫ a−a−∂xu∂xv + 2d∂xuv − k2uv dx+s(sU (l) − u−a)− 2d(s(U(l) − u−a))− k2U (l)+s(sU (r) − ua) + 2d(s(U(r) − ua))− k2U (r),where u±a are the boundary values of u at the left andright boundary.Pole ConditionConsider the equation for the right exterior only, sup-pose for the moment that u is given on [−a, a] and setu′ :=∫ a−ap(∂t)uv dx + ∂xu∂xv − 2d∂xuv + k2uv dxthen the equation for U (r) is given bys(sU (r)−ua)+2d(sU(r)−ua)−k2U (r)−p(∂t)U(r) = u′ .Taking a Laplace transform in time with dual variable ω,p(∂t) corresponds to a multiplication with p(ω) = ω,p(ω) = iω or p(ω) = ω2 depending on the type of equa-tion. Solving for U (r)(s) one obtainsU (r)(s) = (s2 + 2ds− k2 − p(ω))−1(u′ + sua + 2dua).Clearly U (r)(s) is analytic in s except for two poles (ormore generally for several singularities). If s− and s+ arethe roots of (s2 + 2ds − k2 − p(ω)) one can write byCauchy’s integral formulaU (r)(s) =12pi∮γ−(σ2 + 2dσ − k2 − p(ω))−1σ − sdσ+12pi∮γ+(σ2 + 2dσ − k2 − p(ω))−1σ − sdσwhere γ± are paths enclosing s±. In this simple settingthis is equivalent to a partial-fraction decomposition.U (r)(s) =r+(s, u′, ua)s+ − s+r−(s, u′, ua)s− − s,with r± = 1/2(ua ± (u′ + 2dua)/√p(ω) + k2 + d2).Transforming back to space domain we have the corre-spondence1s− − s↔ es−x and 1s+ − s↔ es+x.Suppose that we can identify es+x as an incoming waveor a exponentially increasing solution. Thus depend-ing on the location of the poles s± we can now distin-guish incoming/exponentially increasing waves from out-going/exponentially decreasing waves. So we are in theposition to formulate TBCs as a condition on U (r)(s).The pole condition says: A wave is outgoing if U (r)(s)is an analytic function in the half plane E of possiblelocations of s+. This is equivalent to the condition thatr+ = 0, which yields the classical transparent boundarycondition. But we are not required to form the expressionfor r+.Pole Condition in Hardy spaceHow to handle the pole condition numerically? An-alytic functions can be expanded into power serieses,which convergence in some ball, yet the pole conditionis a condition set on a complex half plane. The Mo¨biustransform is a conformal transformation that transformsa half plane to the unit ball. The Mo¨bius transform isthus the key ingredient to make our algorithm fly. Lets 7→ s˜ = M(s) be the Mo¨bius transform that maps thehalf plane E to the unit ball. We can now reformulate thepole condition: A wave is outgoing if U (r)(s˜) is analyticin the unit ball. ExpandingU (r)(s˜) =∞∑`=0a`s˜` (6)one has to deduce equations for the a`. Then simply trun-cating the series expansion by setting a` = 0 for ` > Lan algorithm is obtained, to realize TBCs.The details are as follows. The Mo¨bius transforms 7→ s˜ = M(s) :=s + s0s− s0maps the half plane {z : <(−z/s0) < 0} onto the unitdisk. (e.g. for positive real s0 the left half plane is mappedonto the unit ball; the imaginary axis is mapped to the unitcircle; −s0 is mapped to 0; and 0 is mapped to −1.) Theinverse is again a Mo¨bius transforms˜ 7→ s = M−1(s˜) := s0s˜ + 1s˜− 1.Space discretizationFor the sake of clearness we consider the case d = 0only. Space discretization is done using third order finiteelements resulting in the standard local mass and stiffnessmatrices, that are assemble to a global system. At theright boundary (and similar for the left boundary) we usethe special exp-element as test functionvs(x) ={e−s(x−a) x ≥ ax−(a−h)ha− h ≤ x ≤ aand obtainp(ω)U (r) + p(ω)u(0)a = u(2)a − k2u(0)a +s0s˜ + 1s˜− 1(s0s˜ + 1s˜− 1U (r) − ua)− k2U (r),(7)where u(0)a and u(2)a are the boundary contributionsu(0)a =∑j∫ aa−hujφjvs dx ; u(2)a =∑j∫ aa−hujφ′jvs dxSetting u′a = (p(ω) + k2)u(0)a − u(2), multiplying (7) by(s˜− 1)2 and rearranging terms yields(s20(s˜ + 1)2 − (s˜− 1)2(p(ω) + k2))U (r)= (s˜− 1)2u′a + s0(s˜2 − 1)ua.Inserting the power series (6), sorting for powers of s˜ andcomparing coefficients yields equations for the aj:(s20 − p− k2)a0 = u′a − s0ua, (8)2(s20 + p + k2)a0 +(s20 − p− k2)a1 = −2u′a, (9)(s20 − p− k2)a0 + 2(s20 + p + k2)a1+(s20 − p− k2)a2 = u′a + s0ua,(10)(s20 − p− k2)a`−1 + 2(s20 + p + k2)a`+(s20 − p− k2)a`+1 = 0, ` = 1, . . . , L(11)with aL+1 = 0. Similar equations hold for the left bound-ary. Transforming back to time-domain equations (8)to (11) yield a system of ordinary differential equationsfor the coefficients aj for j = 1, . . . , L.Take a closer look at (8). If one would choose s0 to betime or ω dependent, s0 =√p(ω) + k2 then (8) is thewell-known exact non-local TBC; equations (8) to (11)decouple and all a` vanish for ` ≥ 2. Choosing s0 to beconstant gives local approximate TBCs.In case of the wave equation (i.e. p(ω) = ω2) choosings0 = ω gives local TBCs. In case k = 0 this choice givesthe exact TBCs.Numerical resultsThe numerical results for the wave equation (1) inte-grated from t = 0, . . . , 20 with an extremely small step-size of ∆t = 10−4 using the trapezoidal rule are shownbelow. The computational domain is [−5, 5], k = 5, theinitial value is a Gaußian u(x, 0) = exp(−x2) and theinitial velocity is set to zero. Space discretization is doneby third order finite elements on an equidistant grid with∆x = 0.002. The reference solution is calculated on adomain [−15, 15]; this way the dominating error compo-nent should be the truncation error in the power seriesrepresentation. Figure 1 shows the evolution of the errorin energy norm for different L. Figure 2 shows the errorin energy-norm vs. the number of coefficients L in thepower series.Extensions and future workThe concepct is easily extended to systemsMp(∂t)u = A∂xxu + 2D∂xu−Ku.with matrices M , A, D and K . These type of systemsarise for example for two dimensional problems on a strip{(x, y), |y| < b,∞ < x < ∞} after a discretisation ofthe y component. The extension to general two or threedimensional problems is currently under investigation.0 5 10 15 2010−1010−810−610−410−2t error  L = 5L = 9L = 13L = 17L = 21L = 25Figure 1: Evolution of the error for different L.5 7 9 11 13 15 17 19 21 23 2510−810−610−410−2Lerror   error at t = 4.96 error at t = 10 error at t = 15 error at t = 20Figure 2: Error vs. L at different t.References[1] T. Hagstrom, “Radiation boundary conditions fornumerical simulation of waves”, Acta Numerica,vol. 8, pp. 47–106, 1999.[2] T. Hohage, F. Schmidt, and L. Zschiedrich, “SolvingTime-Harmonic Scattering Problems Based on thePole Condition I: Theory”, SIAM J. Math. Anal.,vol. 35, pp. 183–210, 2003.

Transparent boundary conditions - the pole condition approach

http://publications.imp.fu-berlin.de/795/1/SchaedleEtAl-PoleCondApproach-Waves2007.pdf

Transparent boundary conditions - the pole condition approach

Abstract

Similar works

Full text

Available Versions

Repository: Freie Universität Berlin (FU), Math Department (fu_mi_publications)