Convergence of energy stable finite-difference schemes with interfaces

Contents lists available atScienceDirect

Journal

of

Computational

Physics

www.elsevier.com/locate/jcp

Short note

Convergence

of

energy

stable

finite-difference

schemes

with

interfaces

Magnus Svärd

a

,

∗

,

Jan Nordström

b

,

c

a_{Dept.ofMathematics,UniversityofBergen,P.O.Box7803,5020Bergen,Norway}

b_{ComputationalMathematics,DepartmentofMathematics,LinköpingUniversity,SE-58183Linköping,Sweden}

c_{DepartmentofMathematicsandAppliedMathematics,UniversityofJohannesburg,P.O.Box524,AucklandPark2006,SouthAfrica}

a

r

t

i

c

l

e

i

n

f

o

a

b

s

t

r

a

c

t

Articlehistory: Received6May2020

Receivedinrevisedform19November2020 Accepted20November2020

Availableonline1December2020

Keywords: Convergencerate Finitedifference Stability Well-posedness Interface

We extendtheconvergenceresultsinSvärdandNordström(2019)[7] forsingle-domain energy-stable high-orderfinite differenceschemes,to includedomainssplitintoseveral grid blocks. The analysis also demonstrates that reflective boundary conditions enjoy the sameconvergenceproperties.Finally, webrieflyindicatethattheseresults (andthe previousonesin[7])alsoholdinmultipledimensions.

©2020TheAuthor(s).PublishedbyElsevierInc.Thisisanopenaccessarticleunderthe CCBYlicense(http://creativecommons.org/licenses/by/4.0/).

1. Background

The convergencerateoffinite-differenceschemes thatare closedattheboundarieswithstencils oflower order accu-racy1 thantheinteriorstencilsisalong-standingproblemthathasbeentreatedextensivelyintheliteratureine.g.[3,4,1,5]. Recentlytheproblemwas addressed in[7],whereitwasshownthatenergy-stable schemes,satisfyingsome natural con-straints,automatically“gain”asmanyordersontheboundaryasthehighestspatialderivative.

Thecasewhenthedomainissplitinto(atleast)twocomputationaldomainsthatareconjoinedataninterface,wasnot treatedin[7].Theinterfacecaseisfundamentallydifferentfromtheboundarycase.Inthetheoryforboundaryconditions, the semi-discretisation ofthe initial-boundary value problemis Laplacetransformed. This allowsa decomposition ofthe modesintothose that decayintothedomain fromtheboundary,andshould be suppliedwithboundary conditions,and those thatdonot. Theboundednessofthe boundarydataandtheenergystability ofthescheme,enableda proofofthe convergenceresult.

Thisproof,however,doesnotimmediatelyencompassinterfaces,sincetheytakedatafromtheconjoiningdomain,which mayormaynotbesufficientlybounded.Moreover,thetreatmentataninterfaceneednotabidetotheminimalnumberof boundaryconditionsrequirement,sincedataataninterfacecanbeconsideredexact.Likewise,theprevious theorydidnot consider reflectiveboundaries,i.e., those thatfeedout-goingwaves backinto thedomainvia thein-going characteristics, sincesufficientboundednessofthesolutionthatisfedbackhastobeestablished.

*

Correspondingauthor.

E-mailaddresses:Magnus.Svard@uib.no(M. Svärd),jan.nordstrom@liu.no(J. Nordström).

1 _{Convergencerateistheratebywhichthesolutionerrordecreases withgridrefinements.Theorderofaccuracyisgivenbythetruncationerrorand} canbedifferentatdifferentpoints.Itrepresentshowaccuratelytheequation(notthesolution)isapproximated.

https://doi.org/10.1016/j.jcp.2020.110020

0021-9991/©2020TheAuthor(s).PublishedbyElsevierInc.ThisisanopenaccessarticleundertheCCBYlicense (http://creativecommons.org/licenses/by/4.0/).

Here,weaddressbothinterfacesandreflectiveboundaryconditionsandwediscusspossiblegeneralisationstomultiple dimensions.Specifically,werequire

(1) thatthesemi-discreteschemeisenergystable,hasapolynomialnullspace(possiblyempty),isnullspaceconsistent,and nullspaceinvariant.(Alltheassumptionsusedin[7].)

(2) that thenumber ofinterface conditions,as seenfrom oneside, is consistent withtheminimal number ofboundary conditions.Thatis,datafromacrosstheinterface canbereplaced bysmooth andboundeddata,andresultina well-posedproblem.

Remark.Thesecondrequirementimmediatelyridsoneoftheproblemsdiscussedabove.Itprohibitsover-specificationattheinterface. Furthermore,itisthemostrobustchoicethatismostoftenusedinpractice.

Theprogrammeforanalysinginterfaceaccuracyistoviewtheinterfaceasanexternalboundaryand,withtheaidofthe requirementsabove,theaprioriestimatesandprevioustheory,infertheoptimalconvergencerates.

Wewilldemonstratethisforafewexamples.Thegeneralisationsareobviousandwillonlybediscussedbriefly. 2. Thehyperboliccase

Weconsidertheadvectionequationontherealline,

ut

+

aux

=

0

,

a

>

0

.

(1)

We splitthedomainintotwo: x

<

0 and x

>

0 withaninterface atx

=

0 anddiscretise eachwithaSummation-by-Parts (SBP) scheme (see [6]). (We assume appropriate trail-off conditions to the left andright and denote the left and right domain as

L

,R.) The gridis equidistantwith indexingxi

=

ih, i

=

0

,

±

1

,

±

2

,

...

,such that x0 representsthe interface in bothdomains.

v_tL

+

a P−1Q vL

=

0 (2)

v_tR

+

a P−1Q vR

= −

P−1a

(

v₀R

−

v₀L

)

e0 where vL

_{= (...}

_vL

−1

,

vL0

)

T and vR

= (

v0R

,

vR1

,

...)

T are the unknowns to the left and right and e0

= (

1

,

0

,

...

0

)

T. The SBP operator has the following properties: P

=

PT

_>

_{0, Q}

₊

_QT

₌

_B

₌

_diag

₍

₋

₁

_,

₀

_,

_...,

₀

_,

₁

₎

_{and P}−1_{Q v is a first-derivative} approximation.(See[7] formoreinformation.)Toadherewiththesecondrequirementabove,thereisonlyapenaltyterm totherightinagreementwiththecharacteristicdirection.

Remark.Generally,onecanpenalisebothequationsatx0and,withaproperscaling,stabilitycanbeproven.Althoughstableand

convergent,sucha“two-way”couplingrenderstheleftproblemill-posedviewedasastand-aloneproblem.Hence,Requirement(2)is notsatisfied.

Weapplytheenergymethodinasomewhatunconventionalway.Multiplythefirstequationin(2) byvL_{P toobtain}





vL



2



t

+

a

(

v L 0

)

2

=

0

.

The left problemis well-posed asa stand aloneproblemand we obtainthe bound vL₀

∈

L2

(

0

,

T

)

. Next,we multiply the secondequationin(2) by

(

vR

₎

T_{P ,}





vR



2



t

−

a

(

v R 0

)

2

= −

2av0R

(

v0R

−

vL0

).

Byrecastingtheboundaryterms,weobtain





vR



2



t

+

a

((

v R 0

−

v0L

)

2

=

a

(

v0L

)

2

.

This problemis bounded, andwell-posed,since vL₀

∈

L2

(

0

,

T

)

. Hence,the theory of[7] appliesto both theleft andright problemsindividuallyandweconcludethat(2) willconvergewithoptimalrates.

Thisprocedurecanbegeneralisedtosymmetricsystems,ut

+

Aux

=

0 discretisedby

v_tL

+ (

A

⊗

P−1Q

)

vL

= +(

A−

⊗

P−1

)((

vL₀

−

vR₀

)

⊗

e₀L

)

(3)

v_tR

+ (

A

⊗

P−1Q

)

vR

= −(

A+

⊗

P−1

)((

vR₀

−

v₀L

)

⊗

e₀R

)

whereeL,R _{arebothoneatx}

₌

_{0 and zeroelsewhere.(FormoreinformationonKroneckerproducts, see[6].)Thedivision}

A

=

A+

+

A− refers to thesigns ofthe eigenvalues such thatonly in-going characteristicsare penalisedon each side of the interface. As in the scalar case, the outgoing characteristics on each side serve as data for the other side and are independentlyboundedinL2

₍

₀

_,

_T

₎

_{suchthatoptimalratesareobtained.}

2.1. Reflectiveboundaryconditions

Another wayof looking at (1) is obtained by transforming the left domain, x

<

0, to

ξ >

0. Thatis, by x

= −ξ

and

∂x

= −∂

ξ,suchthat,

u_tL

−

au_ξL

=

0

,

ξ >

0

.

(4)

Relabelling

ξ

asx,uL

=

u1,uR

=

u2,weobtainthesystem,



u1 u2



t

+



−

a 0 0 a

 

u1 u2



x

=

0

,

x

>

0

withtheboundarycondition, u1

(

0

,

t

)

=

u2

(

0

,

t

)

.(Inparticularwe note thattheproblemiswell-posed ifu2

(

0

,

t

)

=

g

(

t

)

∈

L2

(

0

,

T

)

.) The previously derived apriori estimate guarantees that u1

(

0

,

t

)

∈

L2

(

0

,

T

)

andwe canview theproblemasa well-posedproblemwithareflective externalboundarythatconsequentlyenjoyoptimalconvergencerates.

3. Theparaboliccase

Weconsiderthesamesetupfortheadvection-diffusionequation,ut

+

aux

=



uxx,withaninterfaceatx

=

0 andtrail-off conditionstotheleftandright.TheSBP second-derivativeoperatorisgivenas:D2

=

P−1

(

−

A

+

B S

)

where A is symmet-ric positive semi-definite and S holds afirst-derivative stencilsat theboundary points. Following [2], theapproximation becomes,

P vL

+

a Q vL

=



(

−

A

+

B S

)

vL

+

σ

1e0

(

v0L

−

v0R

)

+

σ

2

e

0

((

SvL

)

0

− (

SvR

)

0

)

(5)

P vR

+

a Q vR

=



(

−

A

+

B S

)

vR

+

σ

3e0

(

vL0

−

v00

)

+

σ

4

e

0

((

SvR

)

0

− (

SvL

)

0

).

In[2] thestablechoicesofthe

σ

’sarederived.Furthermore,these

σ

’shavethefollowingdependenciesofh:

σ

1,3

∼

O(

h−1

)

and

σ

2,4

∼

O(

1

)

.

Considertheleftproblem.Itisanapproximationoftheadvection-diffusionequationwithboundarycondition

α

LuL

+

h

β

LuL_x

=

α

LuR

+

h

β

Lu_xR

=

f

(

uR

,

uR_x

).

(6)

Asufficientconditionforthisstand-aloneproblemtobewell-posedisuR

∈

L2

(

0

,

T

)

andu_xR

∈

L2

(

0

,

T

)

.(Ormoreprecisely,

huR_x

∈

L2

(

0

,

T

)

.)Furthermore,wesimultaneouslydemandthemirroredconditionsontherightside.

Remark.NotethattheinterfaceconditionisanapproximationofuL

(

0

,

t

)

=

uR

(

0

,

t

)

andnotaRobin-typeconditionsince

σ

1,3are

O(

h−1

₎

_and

_σ

2,4are

O(

1

)

.Furthermore,parabolicityrequiresaboundaryconditiononallboundariesimplyingthataninterface

conditiononbothsidesoftheinterfacedoesnotimplyover-specificationofthestand-aloneproblems,inaccordancewithRequirement (2).

Next, we investigatethe stability of the stand-aloneproblems. Thatis, we want thenumerical approximationsof (6) to be bounded in L2

(

0

,

T

)

.By perusing the stability proof in[2], itis clearthat

(

SvL

)

0 and

(

SvR

)

0 are indeed bounded in L2

₍

₀

_,

_T

₎

_{. (See p. 358-9. The eigenvalues associated with}

₍

_SvL

₎

0

+ (

SvR

)

0 and

(

SvL

)

0

− (

SvR

)

0 are non-zero and non-vanishing.) However, a bound on vR

0 or vL0,which isneeded forstand-alonewell-posedness of thetwo problems,is not obtainedfromtheinterface stability.(Thecorresponding eigenvalueiszeroonp.359.)Torecoverthenecessaryestimates on vR₀ orvL₀,weproceedviatheenergyestimateof(5) whichgivesthebounds

sup t∈{0,T}



v L

_{ ≤}

_C

_,

T



0

(

vL

)

AvLdt

≤

C

,

(7) sup t∈{0,T}



v R

_{ ≤}

_C

_,

T



0

(

vR

)

AvRdt

≤

C

.

Sincetheestimatesobtainedwiththe A-matrixabovecorrespondtoL2_{-estimatesofthegradients,weintendtouseSobolev} embeddingtoboundvL _andvR _inL2

₍

₀

_,

_T

_;

_L∞

_(L

D

=

h−1

⎛

⎜

⎜

⎜

⎜

⎜

⎝

−

1 1 0

. . .

0

−

1 1 0

. . .

. .

_.

−

1 1 0 0

⎞

⎟

⎟

⎟

⎟

⎟

⎠

,

such that Dv is an approximation of the derivative in the domain and



Dv



2_h

=

hvTDTDv bounds the L2-norm of the derivativeofv (interpretedasapiecewiselinearfunction).

Lemma3.1.Foranygridfunctionv,vT_Av

_>

_c

_

_Dv

_

2

hforsomec

>

0.

Proof. By assumptionofnullspaceconsistency, A has onezeroeigenvalue associatedwithits rowssummingtozero.(By symmetryitscolumnsthereforealsosumtozero.)Furthermore,itscomponentsare

O(

h−1

₎

_{andweassumethat A isN}

_×

_N. ACholeskyfactorisationisgivenasUTU

=

A.Since A issemi-definite,itwillresultinonezeroonthediagonalofthe upper-triangularmatrixU .Furthermore,thebandedstructureofA ispreservedinU .Since A1

=

0,wemusthaveU 1

=

0. Hence,theelementinthelowerrightcornerhastobe0.

SincetherowsofU sumtozero,wecanwriteU

=

R D,whereR isan N

×

N matrixwiththelastrowbeingzero.(The lastcolumnisalsotakentobezerosinceitdoesnotcontributeanywayduetothebandstructure.)Furthermore,noother rowof R sumto zero,since (most) rowsof A

=

DTRTR D aresecond derivative approximationsandthe two derivatives areaccountedforby DT _{andD.Hence, R}T_{R issymmetricpositivesemi-definitewithonezeroeigenvalueassociatedwith} the last rowandcolumn. Moreover, we denotethe upper-left

(

N

−

1

)

× (

N

−

1

)

submatrix of RTR as M, which isthen symmetricpositivedefiniteandhaselementsof

O(

h

)

.Furthermore,since R isbanded, M

≥

chI forsome c

>

0 whereI is

the

(

N

−

1

)

2_{-identitymatrix.Consequently,}

₍

_Dv

₎

T_RT_R

₍

_Dv

₎

_≥

_c

_

_Dv

_

2

h,sincethezerointhelast rowof RTR coincidewith thezeroinDv andM definesanorm.



By(7) andLemma3.1,wehaveDvL,R

_∈

_L2

₍

₀

_,

_T

_;

_L2

_(L

,R

))

which,togetherwiththeL2boundsin(7) onvL,R andSobolev embedding,impliesvL_i,R

∈

L2

₍

₀

_,

_T

_;

_L∞

_(L

,R

))

foralli.Withtheseaprioriestimates,andsincewealreadyknowfromthe stabilityproofthat

(

SvL,R

₎

0

∈

L2

(

0

,

T

)

,therightproblemboundsthe“data” f in(6) fortheleftproblem(andviceversafor therightproblem).

Remark.Above, we used thebound on

(

SvL,R

₎

0 obtained from the stabilityproof. However, since onlyh

(

SvL,R

)

0

∈

L2

(

0

,

T

)

wasneeded, we couldalsohaveobtainedthatdirectlyfrom theL2

(

0

,

T

;

L2

(L

,R

))

bounds on DvL,R,whichgives

√

hDvL,R

∈

L2

₍

₀

_,

_T

_;

_L∞

_(L

,R

))

.Seeingthat

(

SvL,R

)

0arelinearcombinationsof

(

DvL,R

)i

,theL2

(

0

,

T

)

boundson

√

h

(

SvL

₎

0and

√

h

(

SvR

₎

0 follow.

Weendthissectionwithafewremarks.ThesamereasoningappliestoparabolicsystemssincetheytooyieldL2

₍

₀

_,

_T

₎

boundsonboththevariablesandtheirgradients,whichistheonlyrequirementneededtoapplythetheoryin[7]. Further-more,by mirroringone domainandrewritinga parabolicequation(scalar orsystem)such thattheinterface turnsintoa reflectiveboundary,theconvergencerateofthelatterfollowsimmediately.

4. Multi-dimensions

Apossible problemwiththeprevious procedureis that L∞ isnot embedded in H1 inmultipledimensions. However, in multipledimensionsthe boundednessrequirements are suitablymilder.Consider the 2D advection-diffusion equation,

ut

+

aux

+

buy

=



(

uxx

+

uy y),withperiodicityinthey-direction(andN points).

Using von Neumann analysisin the y-direction reduces the problemto k

=

1

...

N 1-D problems,one foreach mode,

ˆ

vL_k,R inthey-direction.They-derivativeisreducedtoalow-ordertermandcanbeignored.The H1 embeddingimpliesan

L2

(

L∞

)

boundonv

ˆ

k,forallk

=

1

...

N.EachmodeconvergeswithoptimalrateandbyParseval’srelation,sodoesvL,R. Thisreasoningshowsthataboundon

₀T

_∂

L,R

(

u

L,R

₎

2_{dsdt issufficientintheperiodic(intheextradimensions)case,} andgivenbythe H1 embedding.Theoryforthenon-periodiccaseislacking,butitappearsreasonablethatthesamebound wouldbesufficient.

CRediTauthorshipcontributionstatement

MagnusSvärd:Writingoriginaldraft.JanNordström:Reviewingandediting. Declarationofcompetinginterest

The authors declare that they haveno known competingfinancial interests or personal relationships that could have appearedtoinfluencetheworkreportedinthispaper.

References

[1]S.Abarbanel,A.Ditkowski,B.Gustafsson,Onerrorboundsoffinitedifferenceapproximationstopartialdifferentialequations-temporalbehaviorand rateofconvergence,J.Sci.Comput.15(2000)79–116.

[2]M.H.Carpenter,J.Nordström,D.Gottlieb,Astableandconservativeinterfacetreatmentofarbitraryspatialaccuracy,J.Comput.Phys.148(1999).

[3]B.Gustafsson,Theconvergencerate fordifferenceapproximationstomixedinitialboundary valueproblems,Math.Comput.29 (130)(Apr.1975) 396–406.

[4]B.Gustafsson,Theconvergenceratefordifferenceapproximationstogeneralmixedinitialboundaryvalueproblems,SIAMJ.Numer.Anal.18 (2)(Apr. 1981)179–190.

[5]M.Svärd,J.Nordström,Ontheorderofaccuracyfordifferenceapproximationsofinitial-boundaryvalueproblems,J.Comput.Phys.218 (1)(October 2006)333–352.

[6]M.Svärd,J.Nordström,Reviewofsummation-by-partsschemesforinitial–boundary-valueproblems,J.Comput.Phys.268(2014)17–38.