Thread Rating:
  • 0 Vote(s) - 0 Average
  • 1
  • 2
  • 3
  • 4
  • 5
Comparision: 4 different methods of interpolation
#1
Hi -

the question which is the best method to interpolate to fractional/complex tetration is still open, so I thought it would be nice to see this visually at "just-another-example".
I used base b=4, which requires a nasty (complex) fixpoint for the regular tetration and also it is boring to always gamble in the interval of convergence of the infinite tower.

So here is a short text about regular tetration, linear interpolation, polynomial interpolation (just the empirical diagonalization of the truncated matrix using your favourite software) and a log-polar interpolation.
What I'm missing are the computations of that example interval for the Cauchy-integral-method (could someone do this for me?) and the method derived from the slog-ansatz of Walker and Andy Robbins.
Here is the link


]Comparision of interpolations(pdf)

(I'll upload that file to our forum here when the text-version is stable)

Gottfried
Gottfried Helms, Kassel
Reply
#2
Quote:...
What I'm missing are the computations of that example interval for the Cauchy-integral-method (could someone do this for me?) and the method derived from the slog-ansatz of Walker and Andy Robbins.
Here is the link....
Gottfried

Hi Gottfried,

Your post is very interesting. The algorithm I have for generating the sexp(z) is also based on the idea of iterating logarithms, at different values for imag(z)=constant idelta, where the constant slowly moves closer to the real axis from one iteration to the next. This compares to your , ...

If you want a challenge, one thing you could do is generate the limiting behavior of iterating the logarithm, on the path starting with , and ending at . Count the number of times you iterate the logarithm. Then, subtract the fixed point "L". Then multiply by c^n. Then take the logarithm, base c. This will converge to a result, as the number of iterations of logarithms increases arbitrarily large. This resulting function can be analytic, if wrapped around the unit circle, or, equivalently, if represented as a 1-cyclic Fourier series. If it is analytic inside the unit circle, or the Fourier series is analytic for imag(z)>0, then the interpolation you were using will converges to the fixed point at +I*infinity. Also, you can use that 1-cyclic Fourier series correlating to your interpolation, and use that to extend your interpolation to other points that aren't on the initial interpolation line, whose imag(z) values are bigger than your interpolations , which with a little bit of trickery, can even be used to see how your interpolation behaves at values approaching the real axis, although the Fourier series has a singularity corresponding to log(z=0), and the Fourier series is not defined for imag(z)<0.

Here is the equation to use.




I'm also copying a 100 term Taylor series for , generated with my kneser.gp code, updated today in this post. Convergence for this series is best within a unit circle of radius of about "1", and it is accurate to about 32 decimal digits. This should be equivalent to the results for the Cauchy integral, in that the Cauchy integral also converges to the fixed point of L, L* at +/- infinity.
- Sheldon
Code:
sie sexp Taylor series approximation terms
0.9999999999999999999999999999999995070247187929191800883791929298
1.300874303182213070881446827688072813061714556754230656713299112
0.6134882416397192582109374555554626084484838898219689840815017709
0.4626136939158324624570419449320859356621508257631076067661781738
0.2580271665381400290406176995439707729428027624801088785511732678
0.1631972654995318557720687010383484813350480755688599657964371828
0.08962419078985280875791908367441636870291992358851208951306242480
0.05210943962785077456922545881576836319307194597998413171358219472
0.02787344355922481828308550101481237499554959951610037275785833015
0.01539917994225308398452361886392528717634511053121181516280489767
0.008038140872697811905389166736775937873687184851439611721620856142
0.004285814333364039540175648057057942915900950248681577696568660482
0.002190017732657062431267753342735334682355213685867116789645253986
0.001136799557106834835735373410210317366737445255609961738959439943
0.0005702933580844330487522396703324708581367111745349736234225003887
0.0002897924373405595552639524833821840895243704307010392362494910403
0.0001430633962592935882349730645760417255725999434193038877463683821
0.00007143535241514501365950250514888103430851484606325840796413359310
0.00003477052612699922375452114292328458686494995664991541145966486853
0.00001710763661062027438435026942770620065276202889441606559780515340
0.000008222668464823314272610019073460323839550881158819005935582905370
0.000003994852780504333014833969574749273443623562826957028502622836667
0.000001898398063874020939950524705115380873628494076288510587587700622
0.0000009122438146518910571971764421537032161885169426812387017113281077
0.0000004290385948155930958122182972265362537957940313092525320328638684
0.0000002041983306345472171355252468820480355653099911657003733656774073
0.00000009512318872976968042268016824201584840722705732172328549896514731
0.00000004489317526475872755381079755929962866192623484738624738947818941
0.00000002072720854761742214919996912392356563331344531270627679611102753
0.000000009709928580376567409288540894217012554522585199358692221947344054
0.000000004445457635977581093298134874172196352623391145821751830692478538
0.000000002069067831997830355680147212366306199858197057754070520301979538
0.0000000009396576972415515090385897649635544514935693430568984322768929411
0.0000000004348993806810866648462761422110304400679663399357988409106276443
1.959634169379625824533975324107307643182292461794426937648214767 E-10
9.026610708595697059589186571004657226483426982828780751255398719 E-11
4.035887940589495505962961813454430792050597579829137552712606253 E-11
1.851811158741082264579893345045055577900185993420077817637967014 E-11
8.215057781560342499234972804886894570184731791115949762884833338 E-12
3.758182099737354800100972442189710575470534180644706040162342618 E-12
1.653816116453905095315919495158619674267242509776579770628140753 E-12
7.551090292313603512182414144871894949960949212829628314979456250 E-13
3.294730854704294889645798366290620163251202126044696613532797910 E-13
1.503180202578246284601955177533010172432124301605887936254906187 E-13
6.498532414596416354625574128693528367693104525543242144376365141 E-14
2.966787694463398260101586711895208612540948702274735616760 E-14
1.269506066394563995980523451111737527416389957745425996754 E-14
5.809447211675512658406406010859664050117453892558603002430 E-15
2.456923887981265571022295372485400661348078556275414499022 E-15
1.129428858171867988405297967343673929768801862164662346633 E-15
4.711338066637788073516155078449022678560203083166264723447 E-16
2.181606517380487522129598303657716325016581211042840570040 E-16
8.951254371028589263572483822751652023271244878936893432056 E-17
4.190185379602400264475431039792092048187259950236125497032 E-17
1.684676325260767467952730786600274324353305534443332916510 E-17
8.00975310599499114034331211081169770466361095671141938590 E-18
3.139300680337191659721455825117368158842562177220748216374 E-18
1.525412448387190895963160857882397823035180144715226963088 E-18
5.787124927518397255730760899085273078103665520972826591387 E-19
2.897860085139928565881070500794655861108139200741166180865 E-19
1.0539017924966478755296435787405525807942497018701059912826 E-19
5.499741182310150685235453747877087354820770099024176682661 E-20
1.891815428233224354088005687028096685568736250862778647283 E-20
1.044661806434266648611246719755393417569252110062999131835 E-20
3.335460109848631173403946935187917579211177770612486976098 E-21
1.990391550180344368519802729546144052355230903443214604775 E-21
5.742544667618120422282442998171936954247274848707625156635 E-22
3.813938765209951555411152348110770741798315126591233113213 E-22
9.558229147437373582485505316131990712732375060961782289688 E-23
7.372258718528940476709714144299012707345703932341819279402 E-23
1.509540081900308541844775926021439334974435152344193470196 E-23
1.442355372479914264007309507475705715797347663296257830783 E-23
2.17261402936860554290691819576399424158155685466 E-24
2.86605801478900889251802391601744148442882364181 E-24
2.54262312111658143996213690853259003171756198089 E-25
5.80277732661624708084917829800408734781659878345 E-25
1.193043637496520775356974419572489627071440054140 E-26
1.200187232589791625741176000227785330653724342087 E-25
-6.14763319790035079365261034143156019979249851829 E-27
2.539858233552638308377817435553872085374218519626 E-26
-3.13057013130477445048407475360820250796497191757 E-27
5.50123125279568297609245444405153219242798923736 E-27
-1.056288294171983348089454164293866035364041344394 E-27
1.218417104876116850071141113015421401872978232018 E-27
-3.10037140662822825214492277765472433954612767950 E-28
2.75399762285912058468172127529495128841833289673 E-28
-8.49903533151774002213304461649241351169879536428 E-29
6.33607669365821705179359970435612411142282899585 E-29
-2.240662470687793097269648735303213434107446613212 E-29
1.479495589621324521501992703428978451151211123986 E-29
-5.76622626194605031140543603101348276541245427194 E-30
3.49645153513877337162203231356324697384970684598 E-30
-1.460807401340342751076050916111909985924317827766 E-30
8.34263457663592562247448858165024513769708270210 E-31
-3.66155034725519324324462302612346068105890371888 E-31
2.00634211056505813212940132182242432106754109740 E-31
-9.10436035774224627877942444074438826036859283602 E-32
4.86540080643117711358128043740583276363662229656 E-32
-2.243492597703563039982058352212732980162096265650 E-32
1.197927046648240273800047573989178344238550676938 E-32
Reply
#3
(10/14/2010, 11:14 PM)sheldonison Wrote: I'm also copying a 100 term Taylor series for , generated with my kneser.gp code, updated today in this post. Convergence for this series is best within a unit circle of radius of about "1", and it is accurate to about 32 decimal digits. This should be equivalent to the results for the Cauchy integral, in that the Cauchy integral also converges to the fixed point of L, L* at +/- infinity.
- Sheldon
Hi Sheldon -
I've tried your kneser.gp - many thanks: I'll have to study this in more detail later. I reproduced the coefficients for the powerseries for sexp successfully.


\r kneser.gp
init(4)
loop(14)
dumparray


Because I could not yet adapt the "initial value" from the (implicite) 1 to some other value (my z_k() values) (and didn't implement a binary search to find the relative height of the z_k to 1) I've to postpone that .

But just for a quick response ...

What I tried for the beginning was the curve for the complex heights of say h=exp(2*Pi*I *k) and then sexp(h ) for k=0..1 in a small stepwidth. I found that sexp could not handle the full complex circle of h
Here is the output:

pc = vectorv(65,r,sexp(2*Pi*I*(r-1)/64))
-----------
1.00000000000+0.E-105*I
0.994110914974+0.127276777925*I
0.976726543060+0.251971184839*I
0.948667822932+0.371671983801*I
0.911215082470+0.484288987138*I
0.865990762736+0.588164630170*I
0.814820935829+0.682137159840*I
0.759594358541+0.765553628742*I
0.702135585961+0.838238350728*I
0.644103926534+0.900427727697*I
0.586924205857+0.952684757593*I
0.531749806747+0.995806246947*I
0.479454236551+1.03073350627*I
0.430644959306+1.05847405350*I
0.385692356125+1.08003846045*I
0.344767068878+1.09639356436*I
0.307880151008+1.10843113947*I
0.274921939972+1.11694983442*I
0.245697026317+1.12264759756*I
0.219951093578+1.12611839242*I
0.196882860140+1.12728172854*I
0.104936081109+1.05034222516*I
-7.82699075484-7.02415700313*I
-709.298724923-690.557826748*I
-51957.5575671-48491.8027706*I
-3187869.14771-2853512.97028*I
-166166092.408-142895081.478*I
-7450029334.05-6164718982.34*I
-290496242538.-231637469246.*I
-9.94894155289E12-7.65497982943E12*I
-3.01936911598E14-2.24453542988E14*I
-8.18516413940E15-5.88555696918E15*I
-1.99639929780E17-1.39005088025E17*I
-4.40983394463E18-2.97626400168E18*I
-8.87445465596E19-5.81128395357E19*I
-1.63594655086E21-1.04032968135E21*I
-2.77630447925E22-1.71596276853E22*I
-4.35731956614E23-2.61966426185E23*I
-6.35111854280E24-3.71697245818E24*I
-8.63059373889E25-4.92042449645E25*I
-1.09734573161E27-6.09848329535E26*I
-1.30977172945E28-7.10017710659E27*I
-1.47207362265E29-7.78866883271E28*I
-1.56236696609E30-8.07288957626E29*I
-1.57003044882E31-7.92695731129E30*I
-1.49753982260E32-7.39194047837E31*I
-1.35892812751E33-6.56110592764E32*I
-1.17570899338E34-5.55507719216E33*I
-9.71776298482E34-4.49538064973E34*I
-7.68810138740E35-3.48354134756E35*I
-5.83218466008E36-2.58951669954E36*I
-4.24941962635E37-1.84960552474E37*I
-2.97850538389E38-1.27139182283E38*I
-2.01132601530E39-8.42282817727E38*I
-1.31035925744E40-5.38538829227E39*I
-8.24701803610E40-3.32755801839E40*I
-5.02048567638E41-1.98939873261E41*I
-2.95972451395E42-1.15216571225E42*I
-1.69161357840E43-6.47124093462E42*I
-9.38334654197E43-3.52855934257E43*I
-5.05661665865E44-1.86973058846E44*I
-2.64987303104E45-9.63707523115E44*I
-1.35160403892E46-4.83602667122E45*I
-6.71600116592E46-2.36473793281E46*I
-3.25363228864E47-1.12767556062E47*I


I'll give it more time later. Thanks so far!

Gottfried

Gottfried Helms, Kassel
Reply
#4
(10/15/2010, 02:07 PM)Gottfried Wrote: ....
Because I could not yet adapt the "initial value" from the (implicite) 1 to some other value (my z_k() values) (and didn't implement a binary search to find the relative height of the z_k to 1) I've to postpone that .

But just for a quick response ...

What I tried for the beginning was the curve for the complex heights of say h=exp(2*Pi*I *k) and then sexp(h ) for k=0..1 in a small stepwidth. I found that sexp could not handle the full complex circle of h
Here is the output:
...
Hi Gottfried,

I'm glad you could try the Kneser.gp program, I need to add some comments to it, and also explain the updated algorithm.

The problem you're having with imag(z)>2*I, is that sexp(z) implements the Taylor series, so it has a convergence radius of "2", in the imaginary direction, due to the singularity at sexp(-2). sexp(z) works for other values in the real direction, by using exp(TaylorSeries(z)), and log(TaylorSeries(z)). But the Taylor series itself has a radius of convergence of 2, since there is a singularity at sexp(-2). Also, the Taylor series convergence with 100 or so terms is best within a radius of about <1.

Try the following function, which should be reasonably accurate everywhere in the complex plane. It stitches together sexp(z), with the riemzprx(z), which converges best for values of imag(z)>=idelta, which is about 0.016i after 14 iterations. My algorithm works by iterating back and fourth, generating a more accurate riemaprx(z) from sexp(z), and then generating a more accurate sexp(z) from the updated riemaprx(z). Then it updates the array sizes and makes idelta a little smaller.
- Sheldon
Code:
stitchsexp(z) = {
  print(z);
  if (imag(z)>=1, return (riemaprx(z)));
  if (imag(z)<=-1, return (conj(riemaprx(conj(z)))));
  return(sexp(z));
}
pc = vectorv(65,r,stitchsexp(2*Pi*I*(r-1)/64))
Code:
for (s=1,65,print(pc[s]));
1.0000000000000000000000000000000 + 0.E-105*I
0.99411091497356251121660334054347 + 0.12727677792455698141545403455080*I
0.97672654305998578305917759789250 + 0.25197118483864734183151398762201*I
0.94866782293218727219984591701567 + 0.37167198380072132579505764662646*I
0.91121508246980523613033586612553 + 0.48428898713827736928291690621458*I
0.86599076273649904199094947918957 + 0.58816463017016779675563269986350*I
0.81482093582936119822241434165409 + 0.68213715983980310094848042368576*I
0.75959435854105266354457468711565 + 0.76555362874190903322261122971610*I
0.70213558596085847285750652219703 + 0.83823835072842774548779489678433*I
0.64410392653404100704979636884033 + 0.90042772769655387626345342340194*I
0.58692420585720595244283156324057 + 0.95268475759286478172829599485305*I
0.53174980674745764818795178029890 + 0.99580624694707393858489014630534*I
0.47945423655143147463739109030877 + 1.0307335062742955072331330682852*I
0.43064495930562204323634954278206 + 1.0584740535001761782448191001151*I
0.38569235612456845094921072388723 + 1.0800384604495544323525243042957*I
0.34476706887833815417334909473206 + 1.0963935643563012430808924865586*I
0.30788015100843650707563101146503 + 1.1084311394694043905329177454976*I
0.27492194001038214413926826693241 + 1.1169498344674768677723986077593*I
0.24569703876602253338971340665618 + 1.1226476131982804611669541038761*I
0.21995404404932586037347813521729 + 1.1261218964845432608621383877020*I
0.19740960493527573955655725416746 + 1.1278748910342390826684105768548*I
0.17776702936849884487293749742861 + 1.1283220397701643216380955188343*I
0.16073002863623299255628073431661 + 1.1278020154970139874441116391279*I
0.14601236024543245526018230481850 + 1.1265871324641023617450540637499*I
0.13334416296472653516342485574701 + 1.1248934309410459502330586563960*I
0.12247572704559135736508544106469 + 1.1228899870416972600103884620249*I
0.11317934826520773235503638826490 + 1.1207072178116520667083477139737*I
0.10524980412933908586090548559599 + 1.1184441017971748774371045992344*I
0.098503881541559116527366610557030 + 1.1161743321749868424037832109684*I
0.092779286895306502005151548157326 + 1.1139514767511686074259337461195*I
0.087933185984612679199410547995810 + 1.1118132484652623475840711442351*I
0.083840553159988453591219696326586 + 1.1097850007567059445554175501301*I
0.080392455728788036877818724515250 + 1.1078825612100714229445599015055*I
0.077494358790625486575004343963805 + 1.1061145092323597178189207359442*I
0.075064505315839933093265775124783 + 1.1044839924824197901153148522678*I
0.073032404200808376610851753152679 + 1.1029901645356753348872249262035*I
0.071337443393641346948052091906923 + 1.1016293141457599148697176131381*I
0.069927634407947893925195927960497 + 1.1003957451968517946114151542835*I
0.068758487367060749922974475255537 + 1.0992824563833147144830256779775*I
0.067792011158244554068295700263757 + 1.0982816609240691819657330364023*I
0.066995830571445197945752119851721 + 1.0973851791936958574531273262382*I
0.066342410885254327447809150203682 + 1.0965847309297133648948693130044*I
0.065808379828507206066064649667926 + 1.0958721485197513155662294392481*I
0.065373936888691168880042862182397 + 1.0952395286376436974032454744279*I
0.065022340344101677224898321355895 + 1.0946793360427426155349636051232*I
0.064739463015669870641047317468861 + 1.0941844705532570605920106190814*I
0.064513408463790668464414968424662 + 1.0937483059389388676853763144844*I
0.064334180126366560115336601834380 + 1.0933647076538990624484740937655*I
0.064193396661816161171945934128667 + 1.0930280348650106203942826351996*I
0.064084047497062011236948413414600 + 1.0927331310573251150080966931013*I
0.064000283269284013237609881379406 + 1.0924753065593200270793805495699*I
0.063937236483214214778898014606209 + 1.0922503155819538626738985873120*I
0.063890868279932516377395753476112 + 1.0920543297693504886801868088716*I
0.063857837728869286625789537680067 + 1.0918839097854470174685462083430*I
0.063835390514486963192534150305786 + 1.0917359760858963602760903651234*I
0.063821264296554837378338419121594 + 1.0916077797284012606331404037588*I
0.063813608382306521538321775082730 + 1.0914968738417688583645880251576*I
0.063810915664542902701493063497407 + 1.0914010861916872518284773701780*I
0.063811965056336580062258923458939 + 1.0913184931393879726606969700231*I
0.063815772894633215968042430315087 + 1.0912473951798114411841652328590*I
0.063821551995675476170516382327680 + 1.0911862941620766837225532123451*I
0.063828677228417059608758993813316 + 1.0911338722316975944866139433040*I
0.063836656631231378850492509629980 + 1.0910889724868448218285695781494*I
0.063845107235203121006324848070933 + 1.0910505813066151945143721194966*I
0.063853734876755975181544200369394 + 1.0910178122849990427062282998385*I


Reply
#5
(10/15/2010, 03:07 PM)sheldonison Wrote:
Code:
stitchsexp(z) = {
  print(z);
  if (imag(z)>=1, return (riemaprx(z)));
  if (imag(z)<=-1, return (conj(riemaprx(conj(z)))));
  return(sexp(z));
}
pc = vectorv(65,r,stitchsexp(2*Pi*I*(r-1)/64))
Code:
for (s=1,65,print(pc[s]));
1.0000000000000000000000000000000 + 0.E-105*I
0.99411091497356251121660334054347 + 0.12727677792455698141545403455080*I
0.97672654305998578305917759789250 +
(...)
0.063853734876755975181544200369394 + 1.0910178122849990427062282998385*I
Hi Sheldon -

upps - I think my Pari/GP-vector-call was meaningless.
Of course I wanted, that the "height" parameter circles 1 time in the complex plane, thus the last value of sexp should equal the first one in the vectorial output. Sorry. I should have written it in this way:

pc = vectorv(65,r,stitchsexp(exp(2*Pi*I*(r-1)/64)))

Then I get
Code:
4.00000000000+0.E-67*I
      3.87336468717+0.686320949087*I
       3.53286632591+1.25400277240*I
       3.07245040887+1.63627353165*I
       2.59063527737+1.83193139680*I
       2.15456857774+1.88240761546*I
       1.79318667065+1.84024780663*I
       1.50843156756+1.74933262494*I
       1.28927607041+1.63906429915*I
       1.12129007339+1.52628683915*I
      0.991249900203+1.41932158308*I
      0.888608334316+1.32147766690*I
      0.805497695424+1.23341451594*I
      0.736235931369+1.15453275928*I
      0.676770912155+1.08372321281*I
      0.624213185339+1.01973965792*I
     0.576484729285+0.961366385498*I
     0.532068845117+0.907479289431*I
     0.489837394148+0.857054039805*I
     0.448934299533+0.809148995558*I
     0.408699813497+0.762876526969*I
     0.368625526510+0.717369262745*I
     0.328334760078+0.671744585062*I
     0.287586855613+0.625070110831*I
     0.246306948051+0.576334576435*I
     0.204644547462+0.524432811667*I
     0.163062886522+0.468180840110*I
     0.122452775427+0.406386947333*I
    0.0842439806503+0.338012004684*I
    0.0504505580843+0.262443576481*I
    0.0235466506798+0.179859514724*I
  0.00607381145893+0.0915545538344*I
                  -9.01293432323E-55
  0.00607381145893-0.0915545538344*I
    0.0235466506798-0.179859514724*I
    0.0504505580843-0.262443576481*I
    0.0842439806503-0.338012004684*I
     0.122452775427-0.406386947333*I
     0.163062886522-0.468180840110*I
     0.204644547462-0.524432811667*I
     0.246306948051-0.576334576435*I
     0.287586855613-0.625070110831*I
     0.328334760078-0.671744585062*I
     0.368625526510-0.717369262745*I
     0.408699813497-0.762876526969*I
     0.448934299533-0.809148995558*I
     0.489837394148-0.857054039805*I
     0.532068845117-0.907479289431*I
     0.576484729285-0.961366385498*I
      0.624213185339-1.01973965792*I
      0.676770912155-1.08372321281*I
      0.736235931369-1.15453275928*I
      0.805497695424-1.23341451594*I
      0.888608334316-1.32147766690*I
      0.991249900203-1.41932158308*I
       1.12129007339-1.52628683915*I
       1.28927607041-1.63906429915*I
       1.50843156756-1.74933262494*I
       1.79318667065-1.84024780663*I
       2.15456857774-1.88240761546*I
       2.59063527737-1.83193139680*I
       3.07245040887-1.63627353165*I
       3.53286632591-1.25400277240*I
      3.87336468717-0.686320949087*I
              4.00000000000+0.E-67*I
which is at least by the general pattern what was expected. This looks much better... :-)

Gottfried

Gottfried Helms, Kassel
Reply
#6
@sheldonison

Notice:
You posted double same messages but your posts #4 and #5 are bit differents.

Smile
Reply
#7
Here is a plot of that computation:

   

Gottfried
Gottfried Helms, Kassel
Reply
#8
Two more pictures.
Here I compare the polynomial 32x32-interpolation with Sheldon's Kneser-sexp-implementation on the complex unit-circle of the height. The difference is not visible in the plot, it is about 1e-05 in numbers. (This difference is of the order of the error which the polynomial-method with 32x32 also produces, if tet(tet(x,0.5),0.5) - tet(x,1) is computed to check the consistency of two-times half-height-iteration vs one integer iteration: I get a difference of about 1e-5)

   

A plot of the differences (absolute difference between the evaluations of the two methods) shows two surprising *spikes*. Possibly there is a method- or an implementation-specific anomaly/bug?
(note: the h-parameter on the unit-circle 1..I..-1..-I..1 is replaced by the angle phi 0..2*Pi )

   

Gottfried
Gottfried Helms, Kassel
Reply
#9
Gottfried Wrote:A plot of the differences (absolute difference between the evaluations of the two methods) shows two surprising *spikes*. Possibly there is a method- or an implementation-specific anomaly/bug?

I see no evidence of any jumps in any of the first 25 derivatives of sexp(z), centered at z=I, which corresponds to the observed jump at exp(Pi*I/2). I generated these derivatives via Cauchy/Fourier integral, using a sample circle of radius=0.1 about sexp(z=I).
- Sheldon
Code:
0.576484729284574635281197 + 0.961366385497526303664192*I
0.467090029275013068362180 + 0.569916587535513312268876*I
-0.114761722576321231281567 + 0.465827218266975054414939*I
-0.161232170099668615555009 + 0.169045593187154316250322*I
-0.126613456383497199239067 + 0.0568918998688819763780698*I
-0.0762199970552966637382709 - 0.00833069158238119062194641*I
-0.0327098687746058093598951 - 0.0199959613710679405091179*I
-0.0114125765107396462800570 - 0.0167958355490362372839906*I
-0.00129636291030330078307152 - 0.0101328625577288779178557*I
0.00143930066431348953343801 - 0.00476584137840711770103100*I
0.00163769580449573434820048 - 0.00186185520615716786409671*I
0.00109070849328481275176054 - 0.000471938920586012785461553*I
0.000564655613106594209769527 + 0.00000181946046575202049960332*I
0.000245597595059366084519137 + 0.000111877991428776009534293*I
0.0000826274897755178426743768 + 0.0000942697925037007220291984*I
0.0000179529687173593648409477 + 0.0000554432825612102322217693*I
-0.00000285496852089677910122697 + 0.0000268624135777202968381695*I
-0.00000607171248036292730890827 + 0.0000106620295670742196028191*I
-0.00000445871359894296793470106 + 0.00000338692819941567600965410*I
-0.00000246334653238799809994874 + 0.000000592525326678551839351098*I
-0.00000111582497849715807791042 - 0.000000188196751582523640944798*I
-0.000000428067241111207786093519 - 0.000000273619144710082058056908*I
-0.000000126112874499856506266031 - 0.000000186858413233270239477061*I
-0.0000000196315391530546026715320 - 0.0000000976905679317030189981616*I
0.00000000815096620645684871922301 - 0.0000000425620288040458215752126*I
Reply
#10
(10/15/2010, 11:28 PM)sheldonison Wrote:
Gottfried Wrote:A plot of the differences (absolute difference between the evaluations of the two methods) shows two surprising *spikes*. Possibly there is a method- or an implementation-specific anomaly/bug?

I see no evidence of any jumps in any of the first 25 derivatives of sexp(z), centered at z=I, which corresponds to the observed jump at exp(Pi*I/2). I generated these derivatives via Cauchy/Fourier integral, using a sample circle of radius=0.1 about sexp(z=I).
- Sheldon

You're right.
It seems the implementation-artifact is in my routine for the polynomial interpolation.
The current version includes a separation of the height-parameter in the integer( hi = floor(real(h))) and the fractional & imaginary-part ( hf = h - hi) to have the fractional iteration only from a (comfortable) small interval of its parameter. At phi = pi/2 the real(h) changes sign and the floor-function steps. Because with the 32x32-matrices the k-fold-iterates of 1/k-height-stepwidth diverge from the 1-fold iteration of 1-height-stepwidt by about 1e-5 this introduces that error at this height-value. (Perhaps I should introduce some linear stretching for the polynomial method... at least for the comparision-graphs... <umpff>)

Gottfried

Gottfried Helms, Kassel
Reply




Users browsing this thread: 1 Guest(s)