Thread Subject: adding order matters for accuracy ?

Does the order in which we add numbers affect the accuracy of the result in Matlab?

For example, if I add

1/1 + 1/2 + ... + 1/100000 in this order

and, then, if I add

1/100000 + 1/2 + ... + 1/1 in this order

will the adding order generate errors from the point of view of accuracy?

Thank you.

"Ana-Maria Cuculescu" <amcuculescu@yahoo.com> wrote in message
news:h8mqh5$pgb$1@fred.mathworks.com...
> Does the order in which we add numbers affect the accuracy of the result
> in Matlab?
>
> For example, if I add
>
> 1/1 + 1/2 + ... + 1/100000 in this order
>
> and, then, if I add
>
> 1/100000 + 1/2 + ... + 1/1 in this order
>
> will the adding order generate errors from the point of view of accuracy?
>
> Thank you.

Yes.

Rules of thumb: avoid adding a small number to a very large number and avoid
subtracing 2 numbers close to each others in value.

In the first example you showed above, you'll end up adding 1/100000 to a
very large number by the time the sum gets to it.

Try to rearrange the sum order so that you are adding numbers close in
values to each others at each step, you'll avoid these problems. (may be by
reversing the sum order from that you showed in the first case above).

good luck

--Nasser

Thank you. It was a very good explanation.

Regards,
Ana.

Also, it is easy to check for yourself that the order matters in general:


N = 10000;
x = 1./(1:N);
s = sum(x); % Sum from largest to smallest.
t = sum(x(end:-1:1)); % Sum from smallest to largest.
u = sum(x(randperm(length(x)))); % Sum in a random order.
fprintf('%s%16.16f\n','s = ',s)
fprintf('%s%16.16f\n','t = ',t)
fprintf('%s%16.16f\n','u = ',u)

"Matt Fig" <spamanon@yahoo.com> wrote in message
news:h8piqr$gnu$1@fred.mathworks.com...
> Also, it is easy to check for yourself that the order matters in general:
>
>
> N = 10000;
> x = 1./(1:N);
> s = sum(x); % Sum from largest to smallest.
> t = sum(x(end:-1:1)); % Sum from smallest to largest.
> u = sum(x(randperm(length(x)))); % Sum in a random order.
> fprintf('%s%16.16f\n','s = ',s)
> fprintf('%s%16.16f\n','t = ',t)
> fprintf('%s%16.16f\n','u = ',u)

Just to add, to see the true sum, add the following:

x=sym(x);
true_sum=sum(x)

This is the output of the whole run

s = 9.7876060360443482
t = 9.7876060360443855
u = 9.7876060360443446

true_sum =

59731303408577589495741906390072043507599027161503501646678297383697896083301675990834893975828034989358435920082998720957994924884571851140401010658432623657075626196791205705114982585681763155357879419066228023462995179001815927755190029539622783747974164452254841118634890644221051242118078136299790687121414379471444514938982431376789461495970818866668619395322343934190599651481359018815262543960797884256526059820889405312462289377738461838086634382461221512147371413813678031551332103248022895431988668434221719412926685845621694138215890947593777607332120227872412112680993525752634465992057598344322743842636769273796539883437021681567009636559221413921809796946293545983221129067333217238694185641994697012288663188077422578333611895921339100126278583970871916569914605612060097954412538855422964220147309394611206634387804268132940387165100954361283005882917382102024566825241082114417785648506931516628284095236149907352445740757080270298970222298991424163619539379506123415305073900947308841264099732110132982632764632068198883523757157097282779105127947415908699046456803088783327339267892397218242149892531634819382214062724683836006845250776286884486715489190409737569114566129511883502649780640931649212129568035914234262923755559747412923111108520287130262209416913462005671669306269807934053781478249657278385565912518307525179105236506956901850157928214964371634050352901438577263415403688744876342581622281256030733222792180829791937510600332883455183961912573678574754990036711533802509672872155996863868388836727975729335797933388962851319194367493845850089824948702213450036620503830646107545091988313722406707058910017833653545414704249599650049775529073029215168015345293265606873008116247627933888624469681929209379851004445950547606319716967098816714100976792571808262515416127074639075705006903023146298225958490602857446499004364230098847223860901543346076376751962314336238774221150289209211751378432174019280769987477093794111640473517064798428096836658747328519765536656676298196590500678131337665197349138788681259861968351936981313054461987131681144483089340991856988043428025656850727976153708483896512667152931176462081303433355765022311774202512649985872100621715856418215625045464787964541932738773191357538604611479225542759007098027870358424560317443242553747238413952882517769076983690840766115761054346470457283501272775705272932015015574729960746787869871830967427872796160497277962209747924610215822414240041919981175183553527109962710692262381015919022945038699046125939158459536515698091554858912900025718007085057575800241172650965898950619696650228180742127488923350197416513983899495613680697847646864527374243212690287528183527042144807812892270503816902486104365285008915018816738494675706590503199629993335386307535817244728933946168295580830933976459037417011710209560438066324927310670987856025972703788030906730702690013353230858511164403003099629811642317147811098214186664826816112190853609874178781272023056290487159141104770242551077022922240000681046956385030670555583042397235810551992602706921203179878941865859375517849067548479812700691094322083291528056575524330694253212066757962691409448601584679990368892614397167682007332647730715826388640206065927409904606472258205780981624977146587435333894957725950134575030950005497539049808956325901221038754649703965036169657056546360742637792305376452012165429157765596363222949864277396100129549456545591305419783907745882377884078105054475788748702301081388038952336417940438634493102268143466861042836686693916770131202546582345064793779312444132692329868545519091171495526819030966607647613196032465924759505098159065109562700803475781336819122118833950986063496350089772134885342508815568074059916785775991250653654294450398794588358184383834078617089708262257087101187548998458958580044184666521405648523984415244952257156609263761417864479017935110784761636467201575450552556321009717575885944618018433778615054989100841212323410929694159577066483119474451653087883705259075470794023256828427805445735028837916654644021388362210940577930190391582636540850980470296616719034963423846654535422064177667222098347718656546479437040239468796487003250200473537404025025079546001003902525774329756640146471253407001427049473235336540768396996735533484560623411926098144095723014746718672904489789302271079738185744351573488375624241287/6102749046968969734217076666139513429511071629620266328234634310701547772196116277582252815944572761157246832065103066159966237375167796437093471691037743544384611658552391515520060941909772711257881920168722142811068438347932841397591751313659650734251625125355588271140981865532100203057929684374909750850294135717934229443018381215782610060696257186042504564900662727276932348143316811904118605684420663730185618156466173264054611402194624112434683350922762738194635045124739880919367348561980305461302784575159407780579594540237725618190720854571117117430409641454313268182386441724515028668319342120420276383396743219377030833563068623456264424366576076172793126144972015537648728844202912312926348642944309403351112490780691541179368116252636682475281595141012919611527536630624161973150852508905926127445413904848691275051791369203221673190687323458376552035455430124274740491512066468647499724417390270872744614706425172339537623801653979857053478213094618266036880598683915572713793895802586927933555089732433229729501480786211238129505770754815644359505863216675549078434075517322500148278177698501939081283332937002588168511793536597307115213293430360388908417395564156118123812419706746222746212519619912187337949353858977181685220241100864701579729776185987384175012736465341234330865258816791284461005988538717555694508542089506647811881934727406889353514106923440930659478427392491089609545898654359100814867781391463518222860447439866858821253988011770888869020790950787657356093141268624387996215126120657410862224334616258244046294083252008205940771730808325410452731072525298209664651667650971081505438998823327186832446448623344019194515907285346068598534738409522838846510777006412277468467685533963089646003224391621756302500028167816167491512333439104535113742164796271851483916412201713088118127605308970914828170372095093000270497899819271548137534369327374377685976999386818294996123274708217062691928244286999772760779210119211639987293670042546537974809372237292405002551943056143177190863305954914700394055397757522862183261870519110458743980507814952846559711287731182880641754507678600478635280385037672945457683481449989687582281362797528906617641032521523432058778719502580259536565711722226835349758782547581844967353715957763034679128377791394385763923576248256221301061237942707200720335738816554515108703456459227905819197946284295666497512051403098749778403559032544425494417530869808687900681529464389565519724651089950345171385897078379565176441437621284246985330237767312361787942957702866248805184065146444480606003227363241663995437164965158496312299539073963395485530129263865193240092940688049389692913434190027546908810881330469682713866817741818066431534264905737825582978642779307791592620860401336663875745544659038115844549650864590895315828292152699817264686815878164124524376060770338064363389602514676349022877972848977194464816173397269811488279861208702194624432312233525662059905917863068667902130460715669652937788779307051858430798065381840790015022294141326847837553724277532526314679651526216457436645252714558027527743004606975615467572142293313274618345604407999272252620467754145634782194603005226520916047839191067968605667419105324525236889457171366718693547352978309807155472914104321016276905182425406941142158471419167129388543776992620021438559699413531353860586941963934687604421868696227348094876292439630412228727150379368080136940567481066376842729563746554701416733924305098201904356234957654085636680528154361656716334983588260102973208260510882089920165577240155341050550475383807013601338435535783718057961847487996656107316204619731986967327198963787802206698922146471926047705303514591363406634361544152775214249697852015373126755545406545398085465570555645844691904262213915328037729341331940896530753657797438471315741036770153186758318883916115820892057395396795618893484129094048050596222795648559408659672397061567665463218095505770015238651624074909119968589623248826725106380303947966060659097097378869302555557178359416560485928816032985126405979389319280284777207457605421181461807138454400159907422045516533913313345710993828044114789191075875866114404717432236391996647027900922653043754944072697190599691082502222803704814458103130520035816265911562259346488897114507167127147224033663292066021873289732909336618094347985544288684637808103183718425600000

EDU>> vpa(true_sum)

ans =

9.7876060360443822641784779048516

You can see that the sum from the lowest to largest (called 't' above) was
the most accurate. To see this, align the digits

-------------------*
9.7876060360443822641784779048516 --> true
9.7876060360443855 -->t
9.7876060360443482 -->s

See that the digit marked '*', it went wrong first in the 's' sum, which is
the sum from large to small (it changed from 8 to 4). While the same digits
it still correct in the 't' sum, which is the sum from small to large.

--Nasser

"Nasser Abbasi" <nma@12000.org> wrote in message
> Just to add, to see the true sum, add the following:

.... Assuming one has the symbolic toolbox, that is. I do not, but thanks for posting the results for us to look at, Nasser. ;-)

"Matt Fig" <spamanon@yahoo.com> wrote in message
news:h8rrrk$2oo$1@fred.mathworks.com...
> "Nasser Abbasi" <nma@12000.org> wrote in message
>> Just to add, to see the true sum, add the following:
>

> .... Assuming one has the symbolic toolbox, that is. I do not, but thanks
> for posting the results for us to look at, Nasser. ;-)

opps, sorry, I did not see that the symbolic sum was so large, else I would
not have posted it. On the Matlab console, it was one line (but a very long
line, as I can see now).

--Nasser

"Nasser Abbasi" <nma@12000.org> wrote in message
> opps, sorry, I did not see that the symbolic sum was so large, else I would
> not have posted it. On the Matlab console, it was one line (but a very long
> line, as I can see now).
>
> --Nasser

No need to apologize, I wasn't complaining. I was simply pointing out that not everyone has the symbolic toolbox, and thanking you for sharing the light it sheds on this problem.

On 15 Sep, 03:27, "Ana-Maria Cuculescu" <amcucule...@yahoo.com> wrote:
> Does the order in which we add numbers affect the accuracy of the result in Matlab?
>
> For example, if I add
>
> 1/1 + 1/2 + ... + 1/100000 in this order
>
> and, then, if I add
>
> 1/100000 + 1/2 + ... + 1/1 in this order
>
> will the adding order generate errors from the point of view of accuracy?

Well, it is not the order the numbers are added
that generates the error. As you know from working
with pen and paper, the order the numbers are added
does not influence the result, assuming the numbers
are exact. The error is in the numbers themselves
because they are represented inside computer on
inexact floating point formats.

The order the numbers are added influence how
the already present error is magnified throughout
the computation.

Knowing how to handle such issues is an art,
and is among the details that separate brilliant
numeric programmers from mere good ones. You will
find these kinds of questions discussed in texts
on numerical methods.

Rune

"Nasser Abbasi" <nma@12000.org> wrote in message <03esm.40290$ec2.8916@newsfe13.iad>...
>
> "Matt Fig" <spamanon@yahoo.com> wrote in message
> news:h8piqr$gnu$1@fred.mathworks.com...
> > Also, it is easy to check for yourself that the order matters in general:
> >
> >
> > N = 10000;
> > x = 1./(1:N);
> > s = sum(x); % Sum from largest to smallest.
> > t = sum(x(end:-1:1)); % Sum from smallest to largest.
> > u = sum(x(randperm(length(x)))); % Sum in a random order.
> > fprintf('%s%16.16f\n','s = ',s)
> > fprintf('%s%16.16f\n','t = ',t)
> > fprintf('%s%16.16f\n','u = ',u)
>
> Just to add, to see the true sum, add the following:
>
> x=sym(x);
> true_sum=sum(x)
>
> This is the output of the whole run
>
> s = 9.7876060360443482
> t = 9.7876060360443855
> u = 9.7876060360443446
>
> true_sum =
>
>
597313034085775894957419063900720435075990271615035016466782973836978960833016759908348939758280349893584359200829987209579949248845718511404010106584326236570756261967912057051149825856817631553578794190662280234629951790018159277551900295396227837479741644522548411186348906442210512421180781362997906871214143794714445149389824313767894614959708188666686193953223439341905996514813590188152625439607978842565260598208894053124622893777384618380866343824612215121473714138136780315513321032480228954319886684342217194129266858456216941382158909475937776073321202278724121126809935257526344659920575983443227438426367692737965398834370216815670096365592214139218097969462935459832211290673332172386941856419946970122886631880774225783336118959213391001262785839708719165699146056120600979544125388554229642201473093946112066343878042681329403871651009543612830058829173821020245668252410821144177856
485069315166282840952361499073524457407570802702989702222989914241636195393795061234153050739009473088412640997321101329826327646320681988835237571570972827791051279474159086990464568030887833273392678923972182421498925316348193822140627246838360068452507762868844867154891904097375691145661295118835026497806409316492121295680359142342629237555597474129231111085202871302622094169134620056716693062698079340537814782496572783855659125183075251791052365069569018501579282149643716340503529014385772634154036887448763425816222812560307332227921808297919375106003328834551839619125736785747549900367115338025096728721559968638683888367279757293357979333889628513191943674938458500898249487022134500366205038306461075450919883137224067070589100178336535454147042495996500497755290730292151680153452932656068730081162476279338886244696819292093798510044459505476063197169670988167141009767925718082625154
161270746390757050069030231462982259584906028574464990043642300988472238609015433460763767519623143362387742211502892092117513784321740192807699874770937941116404735170647984280968366587473285197655366566762981965905006781313376651973491387886812598619683519369813130544619871316811444830893409918569880434280256568507279761537084838965126671529311764620813034333557650223117742025126499858721006217158564182156250454647879645419327387731913575386046114792255427590070980278703584245603174432425537472384139528825177690769836908407661157610543464704572835012727757052729320150155747299607467878698718309674278727961604972779622097479246102158224142400419199811751835535271099627106922623810159190229450386990461259391584595365156980915548589129000257180070850575758002411726509658989506196966502281807421274889233501974165139838994956136806978476468645273742432126902875281835270421448078128922705038
169024861043652850089150188167384946757065905031996299933353863075358172447289339461682955808309339764590374170117102095604380663249273106709878560259727037880309067307026900133532308585111644030030996298116423171478110982141866648268161121908536098741787812720230562904871591411047702425510770229222400006810469563850306705555830423972358105519926027069212031798789418658593755178490675484798127006910943220832915280565755243306942532120667579626914094486015846799903688926143971676820073326477307158263886402060659274099046064722582057809816249771465874353338949577259501345750309500054975390498089563259012210387546497039650361696570565463607426377923053764520121654291577655963632229498642773961001295494565455913054197839077458823778840781050544757887487023010813880389523364179404386344931022681434668610428366866939167701312025465823450647937793124441326923298685455190911714955268190309666076
47613196032465924759505098159065109562700803475781336819122118833950986063496350089772134885342508815568074059916785775991250653654294450398794588358184383834078617089708262257087101187548998458958580044184666521405648523984415244952257156609263761417864479017935110784761636467201575450552556321009717575885944618018433778615054989100841212323410929694159577066483119474451653087883705259075470794023256828427805445735028837916654644021388362210940577930190391582636540850980470296616719034963423846654535422064177667222098347718656546479437040239468796487003250200473537404025025079546001003902525774329756640146471253407001427049473235336540768396996735533484560623411926098144095723014746718672904489789302271079738185744351573488375624241287/610274904696896973421707666613951342951107162962026632823463431070154777219611627758225281594457276115724683206510306615996623737516779643709347169103774
354438461165855239151552006094190977271125788192016872214281106843834793284139759175131365965073425162512535558827114098186553210020305792968437490975085029413571793422944301838121578261006069625718604250456490066272727693234814331681190411860568442066373018561815646617326405461140219462411243468335092276273819463504512473988091936734856198030546130278457515940778057959454023772561819072085457111711743040964145431326818238644172451502866831934212042027638339674321937703083356306862345626442436657607617279312614497201553764872884420291231292634864294430940335111249078069154117936811625263668247528159514101291961152753663062416197315085250890592612744541390484869127505179136920322167319068732345837655203545543012427474049151206646864749972441739027087274461470642517233953762380165397985705347821309461826603688059868391557271379389580258692793355508973243322972950148078621123812950577075481
564435950586321667554907843407551732250014827817769850193908128333293700258816851179353659730711521329343036038890841739556415611812381241970674622274621251961991218733794935385897718168522024110086470157972977618598738417501273646534123433086525881679128446100598853871755569450854208950664781188193472740688935351410692344093065947842739249108960954589865435910081486778139146351822286044743986685882125398801177088886902079095078765735609314126862438799621512612065741086222433461625824404629408325200820594077173080832541045273107252529820966465166765097108150543899882332718683244644862334401919451590728534606859853473840952283884651077700641227746846768553396308964600322439162175630250002816781616749151233343910453511374216479627185148391641220171308811812760530897091482817037209509300027049789981927154813753436932737437768597699938681829499612327470821706269192824428699977276077921011921
163998729367004254653797480937223729240500255194305614317719086330595491470039405539775752286218326187051911045874398050781495284655971128773118288064175450767860047863528038503767294545768348144998968758228136279752890661764103252152343205877871950258025953656571172222683534975878254758184496735371595776303467912837779139438576392357624825622130106123794270720072033573881655451510870345645922790581919794628429566649751205140309874977840355903254442549441753086980868790068152946438956551972465108995034517138589707837956517644143762128424698533023776731236178794295770286624880518406514644448060600322736324166399543716496515849631229953907396339548553012926386519324009294068804938969291343419002754690881088133046968271386681774181806643153426490573782558297864277930779159262086040133666387574554465903811584454965086459089531582829215269981726468681587816412452437606077033806436338960251467
634902287797284897719446481617339726981148827986120870219462443231223352566205990591786306866790213046071566965293778877930705185843079806538184079001502229414132684783755372427753252631467965152621645743664525271455802752774300460697561546757214229331327461834560440799927225262046775414563478219460300522652091604783919106796860566741910532452523688945717136671869354735297830980715547291410432101627690518242540694114215847141916712938854377699262002143855969941353135386058694196393468760442186869622734809487629243963041222872715037936808013694056748106637684272956374655470141673392430509820190435623495765408563668052815436165671633498358826010297320826051088208992016557724015534105055047538380701360133843553578371805796184748799665610731620461973198696732719896378780220669892214647192604770530351459136340663436154415277521424969785201537312675554540654539808546557055564584469190426221391
5328037729341331940896530753657797438471315741036770153186758318883916115820892057395396795618893484129094048050596222795648559408659672397061567665463218095505770015238651624074909119968589623248826725106380303947966060659097097378869302555557178359416560485928816032985126405979389319280284777207457605421181461807138454400159907422045516533913313345710993828044114789191075875866114404717432236391996647027900922653043754944072697190599691082502222803704814458103130520035816265911562259346488897114507167127147224033663292066021873289732909336618094347985544288684637808103183718425600000
>
> EDU>> vpa(true_sum)
>
> ans =
>
> 9.7876060360443822641784779048516
>
> You can see that the sum from the lowest to largest (called 't' above) was
> the most accurate. To see this, align the digits
>
> -------------------*
> 9.7876060360443822641784779048516 --> true
> 9.7876060360443855 -->t
> 9.7876060360443482 -->s
>
> See that the digit marked '*', it went wrong first in the 's' sum, which is
> the sum from large to small (it changed from 8 to 4). While the same digits
> it still correct in the 't' sum, which is the sum from small to large.
>
> --Nasser
>

You are not really comparing apples to apples here. Start with these two lines:

N = 10000;
x = 1./(1:N);

This creates floating point approximations to the numbers 1/1, 1/2, 1/3, etc. For example, if you download my num2strexact function from the FEX (http://www.mathworks.com/matlabcentral/fileexchange/22239) and use it on x(3) you get the following:

num2strexact(x(3))
ans =
0.333333333333333314829616256247390992939472198486328125

So the actual value is of course only good to about 15-16 decimal digits, as we all know. But the sym command will "fix" that up:

x=sym(x);
x(3)
ans =
1/3

The same is true for the other numbers. So when you subsequently do the sum(x) command, you are summing up the numbers as if they were exactly equal to the original formula, and are not working with the actual floating point values that all of the floating point methods actually start with. So comparing the sym sum to the various floating point sums is not a valid comparison, particularly for an accuracy assessment. What one should do instead is add up the actual floating point numbers using some type of extended precision scheme (e.g. vpa) and compare *that* to the various floating point sums. i.e., make sure everybody starts with exactly the same numbers. E.g., use num2strexact on the floating point array, turn those into vpa numbers, and then add *those* up. e.g.,

digits 100
N = 10000;
x = 1./(1:N);
xx = num2strexact(x);
vx(N) = vpa(0);
for k=1:N
    vx(k) = vpa(xx{k});
end
sumvx = sum(vx)
 
9.787606036044382210194571627970283600461698370054364204406738281250

This number would be the best number to use for assessing accuracy of the various floating point methods, because it starts with exactly the same values as the other methods.

If one were really concerned about the most accurate strictly floating point method, I suppose one could sort the N array elements, replace the two smallest numbers with their sum, resort the N-1 elements, replace ... etc. But having said that, any program that would need to rely on such a scheme to get the *most* accurate answer is probably doomed from the beginning because of the accuracy of the original data.

James Tursa

"James Tursa" <aclassyguy_with_a_k_not_a_c@hotmail.com> wrote in message
news:h8tkq3$hu8$1@fred.mathworks.com...
> "Nasser Abbasi" <nma@12000.org> wrote in message
> <03esm.40290$ec2.8916@newsfe13.iad>...

*snip*

> The same is true for the other numbers. So when you subsequently do the
> sum(x) command, you are summing up the numbers as if they were exactly
> equal to the original formula, and are not working with the actual
> floating point values that all of the floating point methods actually
> start with. So comparing the sym sum to the various floating point sums is
> not a valid comparison, particularly for an accuracy assessment. What one
> should do instead is add up the actual floating point numbers using some
> type of extended precision scheme (e.g. vpa) and compare *that* to the
> various floating point sums. i.e., make sure everybody starts with exactly
> the same numbers. E.g., use num2strexact on the floating point array, turn
> those into vpa numbers, and then add *those* up. e.g.,
>
> digits 100
> N = 10000;
> x = 1./(1:N);
> xx = num2strexact(x);
> vx(N) = vpa(0);
> for k=1:N
> vx(k) = vpa(xx{k});
> end
> sumvx = sum(vx)
>
> 9.787606036044382210194571627970283600461698370054364204406738281250

If you have Symbolic Math Toolbox, you can use:

N = 10000;
x = sym(1./(1:N), 'f');
vpa(sum(x), 100)

sym(z, 'f') writes Z in the form +/- N*2^e with N, e integers -- which
should look a little familiar from the bottom of the first page of this
Cleve's Corner article:

http://www.mathworks.com/company/newsletters/news_notes/pdf/Fall96Cleve.pdf

--
Steve Lord
slord@mathworks.com
comp.soft-sys.matlab (CSSM) FAQ: http://matlabwiki.mathworks.com/MATLAB_FAQ

Rune Allnor wrote:
> On 15 Sep, 03:27, "Ana-Maria Cuculescu" <amcucule...@yahoo.com> wrote:
>> Does the order in which we add numbers affect the accuracy of the result in Matlab?
>>
>> For example, if I add
>>
>> 1/1 + 1/2 + ... + 1/100000 in this order
>>
>> and, then, if I add
>>
>> 1/100000 + 1/2 + ... + 1/1 in this order
>>
>> will the adding order generate errors from the point of view of accuracy?
>
> Well, it is not the order the numbers are added
> that generates the error. As you know from working
> with pen and paper, the order the numbers are added
> does not influence the result, assuming the numbers
> are exact. The error is in the numbers themselves
> because they are represented inside computer on
> inexact floating point formats.
>
> The order the numbers are added influence how
> the already present error is magnified throughout
> the computation.
...

One could construct it such that each of the floating point numbers are,
in fact, exactly representable and the order would still affect the
summation.

IOW, there are two sources of error in the case one is the actual
floating point accuracy and the other is algorithmic.

--

"Matt Fig" <spamanon@yahoo.com> wrote in message <h8piqr$gnu$1@fred.mathworks.com>...
> Also, it is easy to check for yourself that the order matters in general:
>
>
> N = 10000;
> x = 1./(1:N);
> s = sum(x); % Sum from largest to smallest.
> t = sum(x(end:-1:1)); % Sum from smallest to largest.
> u = sum(x(randperm(length(x)))); % Sum in a random order.
> fprintf('%s%16.16f\n','s = ',s)
> fprintf('%s%16.16f\n','t = ',t)
> fprintf('%s%16.16f\n','u = ',u)

There's something I still don't understand here. According to my numerical analysis books, the relative error does not accumulate when you add up numbers of the same sign. So, irrespective of the order of summation, if the relative error in x(i) is bounded by eps(x(i)), then the relative error in s,t, or u ought to be bounded by max(eps(x)).

The absolute error from truth should then be bounded by 10*max(eps(x)) where the 10 here is a bound on the true value of the sum. The differences amongst s,t, and u should therefore be bounded by 20*max(eps(x)).

However, in the modified code below, I am not seeing this bound obeyed. Where is my misconception...?


N = 10000;
x = 1./(1:N);
s = sum(x); % Sum from largest to smallest.
t = sum(x(end:-1:1)); % Sum from smallest to largest.
u = sum(x(randperm(length(x)))); % Sum in a random order.
fprintf('%s%16.16f\n','s = ',s)
fprintf('%s%16.16f\n','t = ',t)
fprintf('%s%16.16f\n','u = ',u)


errbound=2*max(eps(x))*10,
%4.4409e-015

abs(s-t)<=errbound, %false
abs(s-u)<=errbound, %false
abs(u-t)<=errbound, %false

"Matt " <xys@whatever.com> wrote in message <h8ttfl$sl6$1@fred.mathworks.com>...
> "Matt Fig" <spamanon@yahoo.com> wrote in message <h8piqr$gnu$1@fred.mathworks.com>...

>
> There's something I still don't understand here. According to my numerical analysis books, the relative error does not accumulate when you add up numbers of the same sign.

Not sure where this wrong statement you read from. Yes they sure can accumulate:

>> a=[2 repmat(eps,1,100)];
>> sum(a) % Exactly 2

ans =

     2

>> sum(fliplr(a))

ans =

    2.0000

>> sum(fliplr(a))-sum(a)

ans =

  2.2204e-014

>> 100*eps % <- error are accumulated 100 times

ans =

  2.2204e-014

As a final note, the SUM went multithread lately, and the sum order can no longer straight forward predicted (in one experiment on my PC, this happens with array larger than 88998 elements). The same code can give then different result when running on different computer just because of this reason.

% Bruno

"Bruno Luong" <b.luong@fogale.findmycountry> wrote in message <h8tuq7$qpr$1@fred.mathworks.com>...

> > There's something I still don't understand here. According to my numerical analysis books, the relative error does not accumulate when you add up numbers of the same sign.
>
> Not sure where this wrong statement you read from.

From "Introduction to Numerical Analysis" by Stoer and Bulirsch, but the idea is very simple. If the true values being summed are x(i) and have a relative error of u(i), then the computed sum S is given by,

S=sum( x.*(1+u) )

Since the true value of the sum is

S_true=sum(x)

this produces a relative error in S is

RelativeError= ( S- sum(x) )/sum(x)

= sum u.* x/sum(x)

But since all the x have the same sign the coefficients x(i)/sum(x) are a partition of unity, leading to

abs(RelativeError)<=max(abs(u))

So, as long as abs(u(i))<=eps(x(i)) my prediction follows.

"Bruno Luong" <b.luong@fogale.findmycountry> wrote in message
news:h8tuq7$qpr$1@fred.mathworks.com...

*snip*

> As a final note, the SUM went multithread lately, and the sum order can no
> longer straight forward predicted (in one experiment on my PC, this
> happens with array larger than 88998 elements). The same code can give
> then different result when running on different computer just because of
> this reason.

That was a bug that was fixed in the latest release:

http://www.mathworks.com/support/bugreports/532399

--
Steve Lord
slord@mathworks.com
comp.soft-sys.matlab (CSSM) FAQ: http://matlabwiki.mathworks.com/MATLAB_FAQ

"Bruno Luong" <b.luong@fogale.findmycountry> wrote in message
> Not sure where this wrong statement you read from. Yes they sure can accumulate:
>
> >> a=[2 repmat(eps,1,100)];
> >> sum(a) % Exactly 2
>
> ans =
>
> 2
>
> >> sum(fliplr(a))
>
> ans =
>
> 2.0000
>
> >> sum(fliplr(a))-sum(a)
>
> ans =
>
> 2.2204e-014
>
> >> 100*eps % <- error are accumulated 100 times
>
> ans =
>
> 2.2204e-014


I am not sure about that Bruno. sum(a) has no error because eps is smaller than eps(2) by definition. Isn't this only a special case? It seems to me that in the previous examples, for each element of x, x(ii), all other elements x(jj) were addable to x(ii) with eps(x(ii))<= abs(x(ii) - x(jj)) . I.e.,

>> a = [10 repmat(eps,1,4)];
>> sum(fliplr(a))-sum(a)
ans =
     0
>> a = [10 repmat(eps,1,5)];
>> sum(fliplr(a))-sum(a)
ans =
     1.77635683940025e-015


because

sum([eps eps eps eps eps]) > eps(10)/2 is true but not
sum([eps eps eps eps]) > eps(10)/2




Given this, wouldn't you want to use:

a=[2 repmat(eps(2),1,100)];

to make your point? Just wondering...

"Matt " <xys@whatever.com> wrote in message <h8u0tn$g7e$1@fred.mathworks.com>...
> "Bruno Luong" <b.luong@fogale.findmycountry> wrote in message <h8tuq7$qpr$1@fred.mathworks.com>...
>
> > > There's something I still don't understand here. According to my numerical analysis books, the relative error does not accumulate when you add up numbers of the same sign.
> >
> > Not sure where this wrong statement you read from.
>
> From "Introduction to Numerical Analysis" by Stoer and Bulirsch, but the idea is very simple. If the true values being summed are x(i) and have a relative error of u(i), then the computed sum S is given by,
>
> S=sum( x.*(1+u) )
>
> Since the true value of the sum is
>
> S_true=sum(x)
>
> this produces a relative error in S is
>
> RelativeError= ( S- sum(x) )/sum(x)
>
> = sum u.* x/sum(x)
>
> But since all the x have the same sign the coefficients x(i)/sum(x) are a partition of unity, leading to
>
> abs(RelativeError)<=max(abs(u))

The above is about finite precision error due to truncation (and their propagation), not error due when making the sum between two numbers that have widely different magnitudes.

The later can accumulate as various examples have showed.

Bruno

"Bruno Luong" <b.luong@fogale.findmycountry> wrote in message <h8u2ri$qtp$1@fred.mathworks.com>...

> The above is about finite precision error due to truncation (and their propagation), not error due when making the sum between two numbers that have widely different magnitudes.
======

Are these different? The reason summing numbers of different magnitudes is a problem is because one of the numbers needs to have its mantissa barrel-shifted and truncated in order to represent the two numbers with a common exponent.

"Matt Fig" <spamanon@yahoo.com> wrote in message <h8u1s7$kvb$1@fred.mathworks.com>...

>
>
> I am not sure about that Bruno. sum(a) has no error because eps is smaller than eps(2) by definition.

When it depends how each of us respectively define "error". But if you insist, here is another example:

a=[2 repmat(3*eps,1,100)]

s1 = sum(a)
s2 = sum(fliplr(a))
(s1-s2)/eps % <-100

% happy?

% Bruno

"Matt " <xys@whatever.com> wrote in message <h8u3gn$bc8$1@fred.mathworks.com>...
> "Bruno Luong" <b.luong@fogale.findmycountry> wrote in message <h8u2ri$qtp$1@fred.mathworks.com>...
>
> > The above is about finite precision error due to truncation (and their propagation), not error due when making the sum between two numbers that have widely different magnitudes.
> ======
>
> Are these different?

Yes, a subtle difference : in one case we truncate the sum result, in other we truncate the resultant.

Bruno

"Bruno Luong" <b.luong@fogale.findmycountry> wrote in message
> % happy?
>
> % Bruno

% Of Course!

"Steven Lord" <slord@mathworks.com> wrote in message <h8u14i$1oe$1@fred.mathworks.com>...
>
> "Bruno Luong" <b.luong@fogale.findmycountry> wrote in message
> news:h8tuq7$qpr$1@fred.mathworks.com...
>
> *snip*
>
> > As a final note, the SUM went multithread lately, and the sum order can no
> > longer straight forward predicted (in one experiment on my PC, this
> > happens with array larger than 88998 elements). The same code can give
> > then different result when running on different computer just because of
> > this reason.
>
> That was a bug that was fixed in the latest release:
>
> http://www.mathworks.com/support/bugreports/532399
>

Thank you for the info Steve. Indeed I played around witth SUM under 2009B by varying maxNumCompThreads (the soon-to-be-deleted command) and the result is consistent. When using on large array the behavior of SUM command has changed from previous Matlab versions. In short: it is no longer backward compatible. A small issue, nevertheless worth to notice.

Bruno

"Nasser Abbasi" <nma@12000.org> wrote in message <03esm.40290$ec2.8916@newsfe13.iad>...
>
> "Matt Fig" <spamanon@yahoo.com> wrote in message
> news:h8piqr$gnu$1@fred.mathworks.com...
> > Also, it is easy to check for yourself that the order matters in general:
> >

> > N = 10000;
> > x = 1./(1:N);
> > s = sum(x); % Sum from largest to smallest.
> > t = sum(x(end:-1:1)); % Sum from smallest to largest.
> > u = sum(x(randperm(length(x)))); % Sum in a random order.
> > fprintf('%s%16.16f\n','s = ',s)
> > fprintf('%s%16.16f\n','t = ',t)
> > fprintf('%s%16.16f\n','u = ',u)
http://www.wikio.com/article/bad-credit-payday-loans-176415445
> Just to add, to see the true sum, add the following:
>
> x=sym(x);
> true_sum=sum(x)
>
> This is the output of the whole run
>
> s = 9.7876060360443482
> t = 9.7876060360443855
> u = 9.7876060360443446
>
> true_sum =
>
>
597313034085775894957419063900720435075990271615035016466782973836978960833016759908348939758280349893584359200829987209579949248845718511404010106584326236570756261967912057051149825856817631553578794190662280234629951790018159277551900295396227837479741644522548411186348906442210512421180781362997906871214143794714445149389824313767894614959708188666686193953223439341905996514813590188152625439607978842565260598208894053124622893777384618380866343824612215121473714138136780315513321032480228954319886684342217194129266858456216941382158909475937776073321202278724121126809935257526344659920575983443227438426367692737965398834370216815670096365592214139218097969462935459832211290673332172386941856419946970122886631880774225783336118959213391001262785839708719165699146056120600979544125388554229642201473093946112066343878042681329403871651009543612830058829173821020245668252410821144177856
485069315166282840952361499073524457407570802702989702222989914241636195393795061234153050739009473088412640997321101329826327646320681988835237571570972827791051279474159086990464568030887833273392678923972182421498925316348193822140627246838360068452507762868844867154891904097375691145661295118835026497806409316492121295680359142342629237555597474129231111085202871302622094169134620056716693062698079340537814782496572783855659125183075251791052365069569018501579282149643716340503529014385772634154036887448763425816222812560307332227921808297919375106003328834551839619125736785747549900367115338025096728721559968638683888367279757293357979333889628513191943674938458500898249487022134500366205038306461075450919883137224067070589100178336535454147042495996500497755290730292151680153452932656068730081162476279338886244696819292093798510044459505476063197169670988167141009767925718082625154
161270746390757050069030231462982259584906028574464990043642300988472238609015433460763767519623143362387742211502892092117513784321740192807699874770937941116404735170647984280968366587473285197655366566762981965905006781313376651973491387886812598619683519369813130544619871316811444830893409918569880434280256568507279761537084838965126671529311764620813034333557650223117742025126499858721006217158564182156250454647879645419327387731913575386046114792255427590070980278703584245603174432425537472384139528825177690769836908407661157610543464704572835012727757052729320150155747299607467878698718309674278727961604972779622097479246102158224142400419199811751835535271099627106922623810159190229450386990461259391584595365156980915548589129000257180070850575758002411726509658989506196966502281807421274889233501974165139838994956136806978476468645273742432126902875281835270421448078128922705038
169024861043652850089150188167384946757065905031996299933353863075358172447289339461682955808309339764590374170117102095604380663249273106709878560259727037880309067307026900133532308585111644030030996298116423171478110982141866648268161121908536098741787812720230562904871591411047702425510770229222400006810469563850306705555830423972358105519926027069212031798789418658593755178490675484798127006910943220832915280565755243306942532120667579626914094486015846799903688926143971676820073326477307158263886402060659274099046064722582057809816249771465874353338949577259501345750309500054975390498089563259012210387546497039650361696570565463607426377923053764520121654291577655963632229498642773961001295494565455913054197839077458823778840781050544757887487023010813880389523364179404386344931022681434668610428366866939167701312025465823450647937793124441326923298685455190911714955268190309666076
47613196032465924759505098159065109562700803475781336819122118833950986063496350089772134885342508815568074059916785775991250653654294450398794588358184383834078617089708262257087101187548998458958580044184666521405648523984415244952257156609263761417864479017935110784761636467201575450552556321009717575885944618018433778615054989100841212323410929694159577066483119474451653087883705259075470794023256828427805445735028837916654644021388362210940577930190391582636540850980470296616719034963423846654535422064177667222098347718656546479437040239468796487003250200473537404025025079546001003902525774329756640146471253407001427049473235336540768396996735533484560623411926098144095723014746718672904489789302271079738185744351573488375624241287/610274904696896973421707666613951342951107162962026632823463431070154777219611627758225281594457276115724683206510306615996623737516779643709347169103774
354438461165855239151552006094190977271125788192016872214281106843834793284139759175131365965073425162512535558827114098186553210020305792968437490975085029413571793422944301838121578261006069625718604250456490066272727693234814331681190411860568442066373018561815646617326405461140219462411243468335092276273819463504512473988091936734856198030546130278457515940778057959454023772561819072085457111711743040964145431326818238644172451502866831934212042027638339674321937703083356306862345626442436657607617279312614497201553764872884420291231292634864294430940335111249078069154117936811625263668247528159514101291961152753663062416197315085250890592612744541390484869127505179136920322167319068732345837655203545543012427474049151206646864749972441739027087274461470642517233953762380165397985705347821309461826603688059868391557271379389580258692793355508973243322972950148078621123812950577075481
564435950586321667554907843407551732250014827817769850193908128333293700258816851179353659730711521329343036038890841739556415611812381241970674622274621251961991218733794935385897718168522024110086470157972977618598738417501273646534123433086525881679128446100598853871755569450854208950664781188193472740688935351410692344093065947842739249108960954589865435910081486778139146351822286044743986685882125398801177088886902079095078765735609314126862438799621512612065741086222433461625824404629408325200820594077173080832541045273107252529820966465166765097108150543899882332718683244644862334401919451590728534606859853473840952283884651077700641227746846768553396308964600322439162175630250002816781616749151233343910453511374216479627185148391641220171308811812760530897091482817037209509300027049789981927154813753436932737437768597699938681829499612327470821706269192824428699977276077921011921
163998729367004254653797480937223729240500255194305614317719086330595491470039405539775752286218326187051911045874398050781495284655971128773118288064175450767860047863528038503767294545768348144998968758228136279752890661764103252152343205877871950258025953656571172222683534975878254758184496735371595776303467912837779139438576392357624825622130106123794270720072033573881655451510870345645922790581919794628429566649751205140309874977840355903254442549441753086980868790068152946438956551972465108995034517138589707837956517644143762128424698533023776731236178794295770286624880518406514644448060600322736324166399543716496515849631229953907396339548553012926386519324009294068804938969291343419002754690881088133046968271386681774181806643153426490573782558297864277930779159262086040133666387574554465903811584454965086459089531582829215269981726468681587816412452437606077033806436338960251467
634902287797284897719446481617339726981148827986120870219462443231223352566205990591786306866790213046071566965293778877930705185843079806538184079001502229414132684783755372427753252631467965152621645743664525271455802752774300460697561546757214229331327461834560440799927225262046775414563478219460300522652091604783919106796860566741910532452523688945717136671869354735297830980715547291410432101627690518242540694114215847141916712938854377699262002143855969941353135386058694196393468760442186869622734809487629243963041222872715037936808013694056748106637684272956374655470141673392430509820190435623495765408563668052815436165671633498358826010297320826051088208992016557724015534105055047538380701360133843553578371805796184748799665610731620461973198696732719896378780220669892214647192604770530351459136340663436154415277521424969785201537312675554540654539808546557055564584469190426221391
5328037729341331940896530753657797438471315741036770153186758318883916115820892057395396795618893484129094048050596222795648559408659672397061567665463218095505770015238651624074909119968589623248826725106380303947966060659097097378869302555557178359416560485928816032985126405979389319280284777207457605421181461807138454400159907422045516533913313345710993828044114789191075875866114404717432236391996647027900922653043754944072697190599691082502222803704814458103130520035816265911562259346488897114507167127147224033663292066021873289732909336618094347985544288684637808103183718425600000
>
> EDU>> vpa(true_sum)
>
> ans =
>
> 9.7876060360443822641784779048516
>
> You can see that the sum from the lowest to largest (called 't' above) was
> the most accurate. To see this, align the digits
>
> -------------------*
> 9.7876060360443822641784779048516 --> true
> 9.7876060360443855 -->t
> 9.7876060360443482 -->s
>
> See that the digit marked '*', it went wrong first in the 's' sum, which is
> the sum from large to small (it changed from 8 to 4). While the same digits
> it still correct in the 't' sum, which is the sum from small to large.
>
> --Nasser
>

Hi!

> N = 10000;
> x = 1./(1:N);
> s = sum(x); % Sum from largest to smallest.
> t = sum(x(end:-1:1)); % Sum from smallest to largest.
> u = sum(x(randperm(length(x)))); % Sum in a random order.
> fprintf('%s%16.16f\n','s = ',s)
> fprintf('%s%16.16f\n','t = ',t)
> fprintf('%s%16.16f\n','u = ',u)

I've experimeted with different error compensations for the sum: Kahan, Knuth, 2-level Knuth and using long doubles. E.g. Knuth accumulates the numbers with the same quality as with quadrupel precision.
  x = 1 ./ (1:10000);
  s = XSum(x, 'Knuth');
  t = XSum(x(end:-1:1), 'Knuth');
  u = XSum(x(randperm(N)), 'Knuth');
==> s == t == u

If the values to be added are normally distributed with mean 0 and std 1, an error compensation is important. Sorting increases the error and much slower:
  x = randn(1, 1e6);
  exact = XSum(x, 2, 'Knuth');
  sum(x) - exact => 2.09e-011
  sum(sort(x, 2, 'ascend')) - exact => -2.2e-007
  sum(sort(x, 2, 'descend')) - exact => -2.2e-007
  [dum, index] = sort(abs(x));
  sum(x(index)) - exact => 3.29e-12
  XSum(x, 2, 'Knuth2') - exact => 0.0 % as sixtupel precision

Conclusions:
- Sorting (ascending, descending, absolute values) can be better or worse than random order depending on the values.
- Sorting of large arrays needs a lot of time and temporary memory.
- Knuth's error compensation has a higher accuracy, is fast and does not need temporary memory.
- And as C-mex the speed is ok:
  http://www.mathworks.com/matlabcentral/fileexchange/26800
- On modern computers vectors can use Gigs of RAM and error compensating SUM methods get more important.

Thanks to D. E. Knuth, Ogita, Rump and Oishi !
Kind regards, Jan