1IntroductionThepurposeofthisworkistoconstructsmoothinterpolantsforfun - PDF document

1IntroductionThepurposeofthisworkistoconstructsmoothinterpolantsforfunctionsthatareonlyknownonthenodesofaregularrectangulargrid.Cubicsplineshavere-centlybeenusedinthecontextofparticleorvirtualparticletrackingincomplexelds,[1,2,3].Amoreelaboratediscussiononintegratingparticletrajectoriesininterpolatedelds,andtheadvantagesofusinghigherordersplines,istobepublishedelsewere1.Inthisworkjustareviewoftheconstructionofthesplineinterpolationsispresented,withoutrigurousproofs.ConsiderthefunctionsofDvariablesf:RD!R(1)andassumetheycanonlybecomputedonthegridG=DYj=1hjZ;(2)where0hj1arethegridconstantsandhZ=fhzz2Zg;agridcellisdenedas:C(z1;z2;:::;zD)DYj=1[hjzj;hj(zj+1)](3)Basically,thevaluesf(x)areonlyavailableifx2G|thiscanhappenforavarietyofreasons,thesimplestbeingthattheyaremeasuredfromanexperiment.Inthefollowingwewillsaythatwearecomputingapproximationsoff.Inamoreriguroussettingwewouldsaythatwearebuildinganarrayofsplinefunctionsthatconvergetoacertainlimitundercertainconditions;whenthefunctionfissucientlynice(atermthatstillneedsacleardenition),itwillbeequaltothatlimit.Forinstance,oneofthepropertiesofthislimitisthatanyofitsTaylorexpansionsconvergeeverywhere(forexampleanyfunctionwithanitediscreteFourierrepresentation).ApolynomialsplineisafunctiondenedonRD,thatisapolynomialoneachcellandithascontinuousderivativeseverywhereuptoacertainorder.Notethatthesetofgridconstantsfhjgisagiven,andthelimitspokenofbeforeisthelimitofverylargeordersofthepolynomialsenteringthespline.ThislimitwillalwaysexistaslongasfiswelldenedonthegridG;itistruethatforcertainfunctionsthelimitsfordierentgridswillnotbeequal,butthesefunctionsareirrelevanthere. 1C.C.Lalescu,B.Teaca,andD.Carati,\Implementationofhighordersplineinterpolationsfortrackingtestparticlesindiscretizedelds",submittedforpublication.2 AratherinterestingobservationisthatfortheHermitesplinesexactly2D(m+1)Dinputvaluesf(l1;:::;lD)(i1;:::;iD)arerequiredforagivenn=2m+1,whileforgridsplinesqDinputvaluesarerequiredforagivenq(andm2g,thusn2q�3).ThismeansthatgridsplinesgenerallyachieveagivendegreeofsmoothnessfromlessinformationthanaHermitespline|andthepriceisprobablyamuchlargererror.Also,thenumberofinputvaluesisthenumberoftermsinthesumtobecomputed,thusgridsplineswillbefastertocomputethanthecorrespondingHermitesplines(uptoamaximumqthatwilldependontheorder).Asanalnote,themixedderivativesoftheseeldsarealsosmooth:0@NYj=1@ @xjlj1As(n);0@NYj=1@ @xjlj1As(n;q)2Cm�maxfljg(21)(whenevermmaxfljg).4ImplementationofthegridsplinesInpracticetheHermitesplineswillprobablynotbeveryuseful,astheyrequiretoomuchmemory.Otherthanthat,theirimplementationshouldbesimilartothatofthegridsplines.Notethatwemention\parallelizedcodes"inthefollowing;thisrefersspecicallytocasesofcomputerprogramsthatrunonseveralprocessorsatonce,withthememorydividedbetweenthem,andtheseprogramsworkwithphysicalelds,eachprocessorkeepingasliceoftheseeldsinitsmemory.4.1Fullexpressionsofthe1DgridsplinesYouwillneedtouseacomputeralgebrasystem.Forq=4,denethefollowingfunctions:t(0;1)(h)=t0(0;1)+t1(0;1)h+t2(0;1)h2 2!(22)t(1;1)(h)=t0(1;1)+t1(1;1)h+t1(1;1)h2 2!(23)s(3;4)(x)=s(3;4)0+s(3;4)1x+s(3;4)2x2+s(3;4)3x3(24)s(5;4)(x)=s(5;4)0+s(5;4)1x+s(5;4)2x2+s(5;4)3x3+s(5;4)4x4+s(5;4)5x5(25)Afterwards,solvethefollowingsystemsofequations(symbolically):(t(0;1)(0)=f(0);t(0;1)(�1)=f(�1);t(0;1)(1)=f(1);t(1;1)(0)=f(1);t(1;1)(�1)=f(0);t(1;1)(1)=f(2);(26)8�&#x]TJ ;� -1;.93; Td;&#x [00;:s(3;4)(0)=t0(0;1);s(3;4)(1)=t0(1;1);d dxs(3;4)(x)x=0=t1(0;1);d dxs(3;4)(x)x=1=t1(1;1);(27)6 input:fractionoutput: i=(n;q)i()(arrayofdimensionq)2.gridcoordinates input:\normal"pointcoordinates(x1;:::;xD)algorithm:computeeach^xjbxj=hjcandeach~x1xj=hj�^xj.output:{thesetofintegers(^x1;:::;^xD){thesetoffractions(~x1;:::;~xD)3.splineformula input:{thesetoffractions(~x1;:::;~xD){thetypeofspline(n;q){apointer(orsimilarnotion)~ftoanarraycontainingtheinforma-tionaboutthelocaleld(thevaluesoftheeldonthenodesofthecellC(^x1;:::;^xD)andthenecessaryneighbouringcells),shiftedsuchthat~f(0;0;:::;0)=f(x1;:::;xD).algorithm:computethepolynomials(n;q)i(~xj)inthearray ij(bycallingthesplinepolynomialssubroutine),thencomputethesum^f=1+gXi1;:::;iD=0�g~f(i1;:::;iD)DYj=1 ij;(31)notethattestingshowsitismoreecienttointroduceaslittledowhileloopsaspossible|forour3Dimplementation,forq=4thesumiswritteninfullinthesourcecode,asintroducingdowhiletypeloopsslowsdownthecodeconsiderably.Forhighervaluesofqwejusthaveonedowhileloopforoneofthevariables.output:theapproximation^f.4.wrapper input:{\normal"pointcoordinates(x1;:::;xD){apointerftothearraycontainingtheinformationabouttheentireeld{thetypeofspline(n;q)algorithm:putthelocaleldinthearray~ffromfandthencomputetheapproximation^f(usingtheabovesubroutines)output:theapproximation^f8 Thewrapperisveryusefulinthecaseofaparallelizedcode.Alltheop-erationsrelatedtobringingtogetherinformationspreadonpossiblyseveralprocessorscanbeplacedinsidethewrapper,allowingforeasydebuggingandmaintenanceofthecode.Asanexampleofaparallelizedversion,inourimplementationa3Deldisdividedalongthezcoordinatebetweenprocessors.Theeldisperiodicinalldirections,andobtaining~fimpliesalittlecareinregardstothezcoordinate,butitbasicallyjustrequiresanormaluseoftheMODULOoperator.Weimposethateachprocessorhasatleastqnodesonthezdirectioninitsmemory,sothattheformulacontainsatmostinformationfromtwoprocessors(\low"and\up").Wethencompute(31)oneachprocessor,onlyforthenodesthatarein\itsdomain",andwethensumthetworesultingvalues^f=^flow+^fup;theamountofinformationpassedbetweenprocessorsisthuskepttoaminimum.AAppendix:UniquenessConstructionforgeneralcase:t(i;g)(u)=2gXk1;:::;kD=0tk(i;g)DYj=1ukjj(32)s(n;q)(x)=nXk1;:::;kD=0s(n;q)kDYj=1xkjj(33)withi1;:::;iD2f0;1g.TheTaylorexpansions(thecentereddierences)aregivenbythesystemofequationst(i;g)(v)=f(i+v)(34)wherev1;:::;vD= �g;g.Thissystemofequationshasauniquesolutiontk(i;g)f(k)(thereareasmanytk(i;g)astherearef(i+v),andthisisinfactasimpleLagrangeinterpolation).Thesplinesaregivenbythesystemofequations240@DYj=1d dxjlj1As(n;q)(x)35x=i=tl(i;g)(35)whichalsohasauniquesolution(andq2g+2),soformula(20)mustgivethisuniquesolution.Ageneralizationofthismethodwouldbetogotoothertypesofgrids.Forinstanceonecouldimaginethecasewhereinthreedimensionsthereisagridmadeupofequilateraltrianglesinthe(x;y)plane,andsquaresintheothertwoplanes.Whatchangesisthefactthatthesumsareabitmorecomplicated.ThecrucialpropertythathastobepreservedtohaveasmoothapproximationisthattheTaylorexpansionmustbethesamenomatterfromwhichcellweapproachagivennode.9