Ask Your Question
1

plot of ratio of polynomials with coefficients varying from extremely big to extremely small values

asked 2018-03-19 13:52:13 +0200

joaoff gravatar image

updated 2018-03-19 14:13:32 +0200

Does anyone know how to plot a ratio of polynomials where the coefficients of the polynomials vary from extremely big to extremely small values? e.g.

var("omega")

numerator = 1/434425210923785365644729072814805003095647594571477081143688112179*omega^82 - 1/282168800363107583713154019070571550253203986968116550714*omega^80 + 1/376705284799486258761328410091603792725536189789*omega^78 - 1/775984649289050778795273868963329671283*omega^76 + 1/2194552498783346688198647992959*omega^74 - 1/7997055354595006175132*omega^72 + 37/1335393102730856*omega^70 - 455327/89425453412*omega^68 + 6451507166/8156449*omega^66 - 41921205646451/398*omega^64 + 12148488559833861120*omega^62 - 1222797178046083079207911425*omega^60 + 107973546113555797789307958819028992*omega^58 - 8389589284510832786372039497137019252899840*omega^56 + 574047573454150250336814935212839009181133950681089*omega^54 - 34513862440990797385476079081710948473169710438951140982784*omega^52 + 1811420748484519385045717526620321156600093813148380753880114266113*omega^50 - 81787566131583729215998508705538459221433180523307219172710272077796474880*omega^48 + 3074808395164144704086261320345404352244943366045903809072779420014923038874140672*omega^46 - 88158128188124577566928559996352418914343887707417923271162886096982268067866251309875201*omega^44 + 1276001990602982330745542293770646665590092983251428019005953641290986839752143036376864694206464*omega^42 + 50369501995070658072835654251040011409949189464195973029581067274496453807222446814723919181288107409408*omega^40 - 5422646273721796408605794927447900588248745194267196182599017580889204367989136820721458419743182469561086640129*omega^38 + 294713197016069575830983289812182224150427565530453919249784823279683263369786981460098799522297673811884617147560755200*omega^36 - 12215223571528656750794949703498245908596636883090964706566516953450365851729599110977560326879247539894799131557796203328438273*omega^34 + 421566491928721680419007829094184826541933844529415542073312144440199053433647609288595755917957129357084944852204800077564039202340864*omega^32 - 12516212267147404429229647731387721025123388989961149553475986047447463987490660724352358293390873442017645025300898526431390960031953078714369*omega^30 + 324384318879686142437160140587334223782579714788774469706503221987046551849501681825766358071172441792448781390801119733389149415560400234273231077377*omega^28 - 7387292185075462418862140718878310563299525368703143556507183065676300148539804750553952220382528962008497259504339924871164528693646839499632967763533758464*omega^26 + 148147576549670480176424742570505113382496592810105724265276479687408049874080107426078378616519348224098107322066761943275443363035890238181101932475433564088827904*omega^24 - 2614352983825249603402079949177727066634757002976206933577141540374477238146660413804057344751020839726456625121445634542631825250862617604413287932677361615994705332404224*omega^22 + 40472930950614252427322391199345022952821847319764623404423487895464586861905495460361611937362853552737369167685282072188371812469915846107436301498874685121753362553033067593728*omega^20 - 546796479239276734969817343152988884071617215095039987059486720343626501580823236262129279129106406271648455070886215057036604477324486364130961502454958245222249278481585821602157690881*omega^18 + 6398345482963010405651688543654226237482870669667751595384527999380305228718755002303831885277585496704965845755258790817010267468004294456270161640983423358987394850223629580900096084482719744*omega^16 - 64180410520272272525286908376473230431734138441676083665323637745196749056460187115412218416398327877779227771800677121061088325929088938579331150090894803951427168305021791509246069685553904083795968*omega^14 + 544208653291391185309536180524159975422536018745994149714204580484188377200450473903397939463023201727188280315535659623140071365231759008260639937606202290659210954860546302434710221600085853152285286203393*omega^12 - 3827148869236131645057753748864256011706470094795063609207751917613788002030507913841732628668035491051373157906167394821657884872636999755793875871427355836181313542579221345353770527185046940717229548007726252032*omega^10 + 21731619560104735623412248166557020474247938972624063082336680903151973764168130050518970368891893871902596734207676031052134545024607407417007348957518735674156914150659282214614029610623838986362932572240851152089907201*omega^8 - 95767367248890998859080840076211432090249689089635000943797128834097631391862494890486062534305250474327151750740126931966007985728703285377423534457675359002872120478001359450002524654691233467116529177706330479483962206978048*omega^6 + 307398380835349522251122973685605867590098291648854481500785164487397373217556109794667977004956100249926660979665856367573287860824766046271883357679682890511797573992243503579120101379656959151097962039219471247899783414197876948992*omega^4 - 639479912680222769699574758361470606739916016079184564965894672578774414095139913044892342897844958391466586738992899235402376325166165958853889423782954279875630911656947991004581824214687332076761443029188086281653700846286211570814943233*omega^2 + 647222323629001353412618991139040623627834152189367973808816854381809885796462656684768497116303794227646539449954073178816871620703525349577186183599165579353490438287333256004323407793268413622639886629815730954228627165293634675091625356558336

denominator = 1/137761281603135848379924667497137087691869145590694743389457616256325*omega^84 - 1/83086834766162983807691746766150178643427527433316275522827*omega^82 + 1/102667767044170877623975342730336352676592625615621*omega^80 - 1/195052973429231983176184037372982355474872*omega^78 + 1/506764970405918852128607114944932*omega^76 - 1/1689091026941888449775972*omega^74 + 1/6938475729210507*omega^72 - 5933/202766817954*omega^70 + 210652963/41684674*omega^68 - 3785514892410/5023*omega^66 + 98176715539323913*omega^64 - 11274496021929013098315776*omega^62 + 1149760802785547342617432310153216*omega^60 - 104740316526399710014129038296835497656321*omega^58 + 8564675098852218158330202009453158486665628483585*omega^56 - 631107412233357996092481381992488681322377443105002487809*omega^54 + 42040826550209669791480148216500925719160340847286518196157284352*omega^52 - 2538156657414023370036118838246125316621096397799706645672301006268923904*omega^50 + 139157701207028251044471300427394350671149356489494088577874141239880561868144640*omega^48 - 6938806439383168134484769650927631938020954059414471869084090238768620526122595210231809*omega^46 + 314993119151231424855218501630778642353808903970689645637119382608673180637105584990275091038209*omega^44 - 13026311952157375634499439761001799003899804334363328258128277741022032909434938363115295924632479072256*omega^42 + 490832913386158964971261411043934330175709099657034165101518631974970770766908644352744356524979896994217590784*omega^40 - 16848063598754041432195740457548489087523338761447189453678864215130767690468257118676488451717823178276208643556769792*omega^38 + 526508123680113455403918388001376421311837933201943096641237980895995131957184498765311703601193482744568171270729639711997953*omega^36 - 14964040137529112789336002405814742574969978116283981052553508948071227427385569689635300763344329290377090567018871249668589630783489*omega^34 + 386219588613305226787305963970336593407234250843368709064728339440768265801980913517208370322485059415982349821971022988796755219950914764801*omega^32 - 9034392239290768500236248643616567747019665018495660272110199314250181111897744491977335480100194146015876816462940326122858266556363974484342865921*omega^30 + 191046534415310850754936618207854871774036324685618690113988679137516259702053348085287489422888909156122197081115795038125721730271097600962795307182063616*omega^28 - 3640615079002674868268803209292339284671635529598000125401708718259709095257942353116876576750774280418512074696423958832428090501175351634961736513120496483565568*omega^26 + 62273172374254986415765227778440797173097273754339193891373421481657115467189617254708606757292829726836650836046370485734292092085898801412029454712950542649407572017153*omega^24 - 951531652629767939573304547715201546281156310628865493286337971797183608774930095570943761071381252056585979625290081231116086193892640991438658207371842803101206913382977896449*omega^22 + 12911217909000424455313473021628455155346806719690986449907489732516140526838114148296344376049550653353703159891539620541282581173653621727840150941148476570112554868892393249471201280*omega^20 - 154437621072641585160192364210764638535150742299453909620719685765830187029475768884884341432582322941793749721471314800162159259791070934517719047186939224330486171808759549778279712277135361*omega^18 + 1613671375800221558063007173458743660146871072783243100013792246050123478339713407322883872334136082537546559536378103103080560054853188334529829307580805332330325585126522606938881345808027181121537*omega^16 - 14559082441673569008977537112513409409617326807467488077449500018188832706180435122319535717730339986042614839268533169311805347564162239812492671767653840792592722234547429968326169980798063761096993406976*omega^14 + 111744222840860258586904336704874254832742992254602728119554334607230626003726011474866600712329730965417565026295984240724114851734359742329870679353469856198907573647755676148914916780030398787721120780818317313*omega^12 - 715267855855396257045094528442779490633060111715529941383146691210931540467322126692565592976853593736823767984638431457105145557082338533900534849060590000205707733135840156253287305930305666272715830450207103255576577*omega^10 + 3714905633019119124314145982488666373006597207416225659385733815586350823207206537470310618757702145455955559354122274345149874879416613750368567586179362486479162207672443597112579110941483842098201059445524219729129391521792*omega^8 - 15039534985098642489127495070878079816478699521964932372348795943530923555616039497846344032659763466508345351341047837858908515818598488134274908478846308392419344474084655780759561686415523039911831647254805083163507636676075716609*omega^6 + 44523108397807453960903723278683589626375425344092590311786839673418629587601620660013463353874988948445949661864697766839069544771099155966119222987528998667273718783412135043865598080899844899439778708053188505321146926956138446907768832*omega^4 - 85727466500775256683870209515231431180659141489998596286522465643370098508698223577125979654943088313884931697767434884070649253899301725049084279835311596901096918650207437145768941287635219222105620194647347139209573654378699720526506978443265*omega^2 + 80567082948458938117841918763525355285225446281873107502306332907427746631501737918800317997386306466900211456674973441765747419882727073257243150087944844059450542039822520151117935401484171322088826759710569301983550807117385762337391493721136037889

#numerator = 1/434425210923785365644729072814805003095647594571477081143688112179*omega^82 + 647222323629001353412618991139040623627834152189367973808816854381809885796462656684768497116303794227646539449954073178816871620703525349577186183599165579353490438287333256004323407793268413622639886629815730954228627165293634675091625356558336

#denominator = 1/137761281603135848379924667497137087691869145590694743389457616256325*omega^84 + 80567082948458938117841918763525355285225446281873107502306332907427746631501737918800317997386306466900211456674973441765747419882727073257243150087944844059450542039822520151117935401484171322088826759710569301983550807117385762337391493721136037889

gd = numerator/denominator
plt3 = plot(gd, (2*pi*0.98*10^3, 2*pi*1.02*10^3))
plt3.show(gridlines=True)

As you can perceive, I tried to rationalize the coefficients, even though, I get the following message, with an empty plot

verbose 0 (3757: plot.py, generate_plot_points) WARNING: When plotting, failed to evaluate function at 200 points.
verbose 0 (3757: plot.py, generate_plot_points) Last error message: '(34, 'Numerical result out of range')'

The plot should resemble an "M" letter centered at 2pi10^3.

edit retag flag offensive close merge delete

1 Answer

Sort by ยป oldest newest most voted
1

answered 2018-03-26 15:58:30 +0200

slelievre gravatar image

Not sure exactly how to explain the problem you get when plotting, but here is a way to get an actual plot.

After defining gd as in the question, the following gets you a plot:

tau = 2 * RDF.pi()
plt3 = plot(lambda x: gd.subs({omega: QQ(x)}), (tau*0.98e3, tau*1.02e3))
plt3.show(gridlines=True)

With this, I get a plot which looks gently curved (not an M shape).

Notes

The trick with RDF.pi() avoids keeping a symbolic pi lying around.

The other trick is to substitute rational values into gd, so that gd(x) is computed as a rational, returning a rational number (computed in arbitrary precision); this rational can then be transformed into a floating point number for plotting.

The naive plotting uses a floating-point omega, so that the various terms in the numerator and denominator are computed separately as floating point numbers, and approximation errors get out of control with the high powers of omega in numerator and denominator.

edit flag offensive delete link more

Comments

Well, it is embarrassing, this was not the result that I was expecting, but thank you! I will revise my computations and try to simplify the polynomials and improve the calculations in some way.

Interesting to see the lambda function used in here. It is not the first time it solves my problems.

By the way, regarding the notes, I had tried to rationalize the bounding values of "omega" in the plot function, but nothing changed (sorry for not reporting this). In this regard, I haven't understood the RDF.pi() trick. I am not sure why you have transformed the symbolic pi to floating point and what is the problem of having symbolic pi's laying around.

joaoff gravatar imagejoaoff ( 2018-03-26 22:10:44 +0200 )edit

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

1 follower

Stats

Asked: 2018-03-19 13:52:13 +0200

Seen: 167 times

Last updated: Mar 26 '18