A NEW FINITE DIFFERENCE METHOD FOR THE HELMHOLTZ E

ANEWFINITEDIFFERENCEMETHODFORTHEHELMHOLTZEQUATIONUSINGSYMBOLICCOMPUTATIONLARRYA.LAMBE,RICHARDLUCZAK,ANDJOHNW.NEHRBASSAbstract.AnewnitedierencemethodfortheHelmholtzequationispre-sented.Themethodinvolvesreplacingthestandard\weightsinthecentraldierencequotients(Sects.2.1,2.2,and2.3)byweightsthatareoptimalinasensethatwillbeexplainedintheSects.justmentioned.Thecalculationoftheoptimalweightsinvolvessomecomplicatedanderrorpronemanipu-lationsofintegralformulasthatisbestdoneusingcomputeraidedsymboliccomputation(SC).Inaddition,wediscusstheimportantproblemofinterpo-lationinvolvingmeshesthathavebeenrenedincertainsubregions.AnalyticformulaearederivedusingSCfortheseinterpolationschemes.OurresultsarediscussedinSect.5.SomehintsaboutthecomputermethodsweusedtoaccomplishtheseresultsaregivenintheAppendix.Moreinformationisavailableandaccesstothatinformationisreferenced.WhilewedonotwanttomakeSCthefocusofthiswork,wealsodonotwanttounderestimateitsvalue.ArmedwithrobustandecientSClibraries,aresearchercancomfortablyandconvenientlyexperimentwithideasthatheorshemightnotexamineotherwise.1.IntroductionAstandardcomputationaltoolforapproximatingsolutionstosystemsofpartialdierentialequationswithboundaryconditionsisthenitedierencemethod[1][2][3].InthecaseoftheHelmholtzequationr2u=2u;(1)itwasshowninChapt.fourof[4]thatnumericalerrorsthatcanoccurinthecentraldierencequotientscouldbecorrectedwithoutincreasingtheorderbyoptimallyadjusting\weightsforthesequotients.ThismaterialwillbereviewedinSect.2.TheproblemofaccuratelyinterpolatingwhenameshisrenedinsomesubregionofinterestforagiveninitialmeshfortheHelmholtzequationwasalsoaddressedinChapt.veof[4]andthiswillbereviewedinSect.4.Thecalculationsneces-saryforthisinvolvedworkingoutsomerathercomplicatedintegralformulasandmanipulatingcomplexalgebraicexpressions,allofwhichcanbetediousanderrorprone.Inthispaper,itwillbeshownhowSCcanbeusedtorelievethetediumandeliminatetheinevitabletypographicalerrorsinvolvedinhandcalculations.2.OptimizingWeightsinaFiniteDifferenceMethodfortheHelmholtzEquation2.1.Dimensionone.Considerafunctionuofonevariable.Classicalnitedif-ferenceschemesarederivedundertheassumptionthatucanbeexpandedinterms12LARRYA.LAMBE,RICHARDLUCZAK,ANDJOHNW.NEHRBASSofaTaylorseries.Onethenhasu(x+h)=u(x)+u0(x)h1!+u00(x)h22!+:::;u(xh)=u(x)u0(x)h1!+u00(x)h22!:::(2)fromwhichtheapproximationsu0(x)u(x+h)u(xh)2h;u00(x)u(x+h)2u(x)+u(xh)h2(3)areeasilyderived.Theideain[4]istoreplacethecoecienttwointheapproxi-mationabovebyanewcoecient!whichminimizesu00u(x+h)!u(x)+u(xh)h2:(4)Toderivethisnewweight,consideru(x+h)+u(xh)ingeneral.Wehaveu(x+h)+u(xh)=2u(x)+u(2)(x)h22!+u(4)(x)h44!+:::(5)fromequations(2)above.Byiteratingtherelationu00=2u,itfollowsthatu(2n)=(1)n2nu:(6)SubstitutingthesevaluesintoEq.(5)givesu(x+h)+u(xh)=2cos(h)u:(7)Thus,theequation2u=u00(x)=u(x+h)!u(x)+u(xh)h2(8)hasanexactsolutioninthiscase,viz.!=2cos(h)+(h)2:(9)andwehavean\adjustedweightforanewnitedierencescheme.Remark2.1.Inthiscase,itisactuallywellknownthattheexactsolutiontotheonedimensionalHelmholtzequationisgivenbyekxj+ekxj(j2=1)wheretheconstantsandaredeterminedbytheboundaryvaluesandtheadjustedweightabovecanbedirectlycalculatedfromthisasisdonein[4].2.2.Dimensiontwo.Eq.(1)indimensiontwoisuxx+uyy=2u(x;y):(10)Theclassiccentraldierenceschemegivestheapproximationuxx+uyyu(x+h;y)+u(xh;y)4u(x;y)+u(x;y+h)+u(x;yh)h2(11)andwewanttoreplacethecoecientfourabovebyanoptimal(inasensetobemadeprecise)weight!.Unfortunately,themethodofthelastsectiondoesnotdirectlygeneralizeduetotheexistenceofcrossderivativeterms.However,anargumentwasgivenin[4]whichgivesasatisfactoryanswerthatisoptimalinthesensethatwerecallhere.ANEWFINITEDIFFERENCEMETHODFORTHEHELMHOLTZEQUATIONUSINGSC3ConsiderEq.(10)intheabsenceofboundaries.Itisknownthattheplanewavesf(x;y)=ej(xcos()+ysin())(12)aresolutions(j2=1).Astraightforwardcalculationgivesf(x+h;y)+f(xh;y)+f(x;y+h)+f(y;xh)=2(cos(hcos())+cos(hsin()))f(x;y):(13)Thus,inthiscase,wewouldliketond!suchthat2(cos(hcos())+cos(hsin()))f!fh2(14)isascloseto2faspossible.Equivalently,wewanttondanoptimal!suchthat2(cos(hcos())+cos(hsin()))(!2h2)0:(15)Eq.(15)doesnothaveasolutionthatisindependentoftheangle,however,itisreasonabletotrytominimizetheaverageoverallanglesandhopethatthereisauniquesolution.Inotherwords,weseekasolution!totheequationZ20(2(cos(hcos())+cos(hsin()))(!2h2))d=0:(16)Equivalently,2!=2Z20(cos(hcos())+cos(hsin()))d+22h2;(17)andthereisauniquesolution!=1Z20(cos(hcos())+cos(hsin()))d+2h2:(18)Infact,thereisananalyticexpressionforthisintegralintermsoftheBesselfunctionoftherstkind.WehaveZ20(cos(hcos())+cos(hsin()))d=4J0(h):(19)Thisfollowsfromtheclassicformula[5](wherez=x+yj)J0(z)=1Z0cos(zsin())d=1Z0cos(zcos())d(20)andtheeasilyprovenformulasZ0cos(zsin())d=Z2cos(zsin())d;Z0cos(zcos())d=Z2cos(zcos())d:(21)Thus,indimensiontwo,wehavetheadjustedweight!=4J0(h)+(h)2:(22)Notethat,inthiscase,theintegralinvolvediseasyenought