# When I use Plot to plot an interpolating function which is the solution of a big ODEs solved by NDSolve,I found it takes a long time

I'm a Chinese student.I can't speak English very well.If there are any mistakes,please forgive me. I use NDSolve to solve a big ODEs.

tmax = 8640000*4;
planetcounts=33;
sol = NDSolve[equs,Flatten[{positionvectors, velocityvectors}], {t, 0, tmax} ];
evectors = {Table[Subscript[x, n, m][t], {n, 1, planetcounts, 1}, {m, 1, 3, 1}], Table[Subscript[vx, n, m][t], {n, 1, planetcounts, 1}, {m, 1, 3, 1}]} /. sol;
epositionvectors = evectors[[1, 1]];
evelocityvectors = evectors[[1, 2]];
DumpSave["planetorbit.mx", {epositionvectors, evelocityvectors}];


Then I save the solution,there is 33*3*2 interpolating functions in total.The size of the file(planetorbit.mx) is about 430MB. Then I plot one of the interpolating functions.

Get["planetorbit.mx"];
tmax = 8640000*4;
Plot[epositionvectors[[4, 2]], {t, 0, tmax}]


I can plot the graph of the function,but it takes a little long time.It taks about tens of seconds. Then I lookup the document of InterpolatingFunction.I use N to make a new InterpolatingFunction with numerical values of all the data,then I plot the new InterpolatingFunction.

Get["planetorbit.mx"];
tmax = 8640000*4;
if = N[epositionvectors[[4, 2]]];
Plot[if, {t, 0, tmax}]


Now it's much faster than before,it takes less than a second. But when I want to plot lots of interpolating functions faster,I try to turn all the InterpolatingFunctions in a Table into new InterpolatingFunctions with numerical values of all the data. I try to use N directly on a Table.

Get["planetorbit.mx"];
tmax = 8640000*4;
eepositionvectors=N[epositionvectors];
Plot[eepositionvectors[[4, 2]], {t, 0, tmax}]


Or use N on every element in the Table.

Get["planetorbit.mx"];
tmax = 8640000*4;
planetcounts=33;
eepositionvectors = Table[N[epositionvectors[[n, m]]], {n, 1,planetcounts, 1}, {m, 1, 3, 1}];
Plot[eepositionvectors[[4, 2]], {t, 0, tmax}]


But both of them doesn't works,it still takes tens of seconds. So,is there a way to turn all the InterpolatingFunctions in a Table into new InterpolatingFunctions with numerical values of all the data?So that I can plot them much faster.

There is the program I used to generate planetorbit.mx

iinput = {{{2570440309792.98,-136860685075.315,-1482235526242.25},{-1294.1476850743,-1789.16293429772,3025.07567697685}},{{-339554588574.549,67683927452.2663,-660150670223.064},{10698.6145003882,-1368.97734524111,-447.310042542305}},{{-301241541803.123,-27923669989.7742,-224389650127.63},{5036.67534881656,-579.109182034442,-15757.3273997457}},{{5408324630930.4,194446178532.619,884475126419.986},{-2164.60069863904,855.72835769034,4478.89130642752}},{{270411152308.477,-4179468045.56236,-265088975592.963},{18946.6311071931,-2178.26039613655,14527.1972179053}},{{973140539318.019,-120131674589.468,1048204522461.09},{-6331.71647706687,-3172.06930083491,-3082.32951360313}},{{-519434624.835267,89477127.874234,1269609497.13248},{-16889.0045086371,3507.51531741377,-23677.7277738523}},{{973542357944.268,-120204767197.587,1048518373662.5},{-18550.9339217503,6020.36854395739,7835.33840532535}},{{-340810560587.677,67814435915.4726,-660356579138.102},{15217.8565127673,-2254.93453229351,4300.44629877502}},{{-341839584134.79,67768478371.4825,-659144033618.068},{629.275969072511,782.576557522181,-12449.101033493}},{{975715403326.939,-121237035948.568,1045727875589.62},{-4819.1637362358,533.559607723971,8292.81650778346}},{{-341538123827.744,67891697717.6915,-660546386105.968},{27484.0501782225,-1499.72373665824,-12477.4163199288}},{{-341444789525.154,67839331912.7634,-660136681384.9},{10698.7829553519,-513.696747250918,-8545.78072185835}},{{270391294130.213,-4187280427.6221,-265098725843.189},{19579.8026803021,-2306.64980315873,13340.8046052286}},{{162212886066.433,-3898411371.18696,-112078410586.023},{34218.5085798306,-2018.07989836899,47380.8510588956}},{{973554253730.111,-120230070717.025,1048458921357.1},{-19685.1274600069,6970.09516895006,10479.6051378158}},{{2570363160927.42,-136999954181.344,-1482271711006.45},{396.218755612714,-6777.16061476125,5785.49164348138}},{{-751001169.372561,147200989.996129,966576840.349661},{-16120.8370736851,3564.99893446096,-24339.0537300266}},{{1700960076115.86,-354393160837.101,4147792755844.68},{-5010.47273859918,437.781400473648,1920.6463361279}},{{2570172714623.67,-137540773374.887,-1482317757505.45},{5145.47599241715,-2309.63396077755,7939.19849209034}},{{270382932217.381,-4190198376.84808,-265101540144.164},{20325.4759829074,-2596.51737135682,11338.2330398433}},{{5408316398840.54,194433857384.827,884487950663.342},{-2105.81082501867,990.826415194224,4646.44417128792}},{{973349896457.649,-120039324168.892,1048752827960.96},{-14859.0665210035,2757.18216076937,2943.41072048841}},{{973485110524.187,-120256084320.206,1048292136162.41},{-7372.58107363601,466.792171877495,6225.52683140201}},{{123017008329.084,-6016868826.22613,-87204388878.7599},{16536.9037211291,-3414.47391977056,23414.5902100918}},{{973584679708.929,-120213811906.906,1048566216437.7},{-17064.417928125,5719.39476591529,8936.62631149546}},{{972751990087.101,-120249988765.447,1047340821926.26},{-3538.01727721721,-2284.70133409074,3086.24008611723}},{{2570643838902.96,-136658694207.763,-1482119568576.2},{647.610552650719,901.51776029384,4002.23178792455}},{{1700683780251.11,-354615317091.889,4147780019670.37},{-4429.73246153819,-527.720276533523,6163.47908897866}},{{2570263767310.07,-136954179444.37,-1482351320720.03},{1460.38115756068,-5045.50668599653,5659.11274687601}},{{2570455255553.85,-137045867032.889,-1482192269433.8},{3123.35991290122,-697.503045523929,6139.36870330027}},{{102455434171.624,-8993398535.93708,19108136199.6033},{-34337.2525123439,2141.95757647127,-6341.25391312778}},{{444130984675.576,-29528782328.4646,-292355888246.344},{9016.32912034143,148.865757329006,15431.1496711581}}};
rawallaxes = {{0,-2.88487317448016*^-05,2.12253045001054*^-12},{0,-4.35332458437188*^-06,-1.2973918272173*^-13},{0,-0.000192343941307627,-2.5078800147682*^-11},{0,-1.13675769171095*^-05,0},{0,-5.76004531467333*^-05,5.14988138938799*^-12},{0,-2.65377038886072*^-05,7.9088509829725*^-13},{0,-7.29211460566148*^-05,-4.34643944099977*^-12},{0,-5.2905576012563*^-05,1.57670903243312*^-12},{0,-2.04663683689432*^-05,9.14917881039701*^-13},{0,-1.01579316833522*^-05,-1.43796758354797*^-12},{0,-9.16841258913337*^-07,-1.36619993882378*^-14},{0,-4.10561769967899*^-05,-3.05892377400196*^-12},{0,-0.000175851819221862,7.20608989146498*^-12},{0,-7.08821607986465*^-05,-2.11245300765772*^-12},{0,-1.24001303447585*^-06,0},{0,-7.6770840678364*^-05,1.1439746719355*^-12},{0,-5.14248822582886*^-05,0},{0,-2.66627057499136*^-06,0},{0,-0.000108330175862648,-7.13799202037291*^-12},{0,-5.40129440196324*^-06,-4.82913405444596*^-13},{0,-0.000227858836296946,-2.8011733557709*^-11},{0,-1.13675096145016*^-05,6.77556318343098*^-13},{0,-1.60886793310056*^-05,1.43843986201431*^-12},{0,-0.000165121018653736,-1.45082713204348*^-11},{0,-2.86532963400532*^-06,0},{0,-3.8443598896265*^-05,0},{0,-4.55995768788853*^-06,-1.35897329088704*^-13},{0,-8.35256923892302*^-06,2.48925961223931*^-13},{0,1.23743166113854*^-05,0},{0,-1.75462064362364*^-05,5.22917700893771*^-13},{0,-0.000101237666967791,0},{0,2.99244959478528*^-07,-2.2295484347081*^-14},{0,-0.000326718291034922,-1.94739294134827*^-11}};
soibody = "Sun";
bios = {0, 0, 0, 864000};
G = 6.67384*10^(-11);
tmax = 8640000*4;
planetname = {"Ariel",
"Callisto",
"Ceres",
"Charon",
"Deimos",
"Dione",
"Earth",
"Europa",
"Ganymede",
"Iapetus",
"Io",
"Jupiter",
"Mars",
"Mercury",
"Mimas",
"Miranda",
"Moon",
"Neptune",
"Oberon",
"Phobos",
"Pluto",
"Rhea",
"Saturn",
"Sun",
"Tethys",
"Titan",
"Titania",
"Triton",
"Umbriel",
"Uranus",
"Venus",
"Vesta"};
planetdata =
{{8.346344431770477*^10, 0, 0},
{7.179289361397270*^12, 0, 0},
{6.26325000000000*^10, 0, 0},
{1.058799888601881*^11, 0, 0},
{9.615569648120313*^4, 0, 0},
{7.311636648732*^10, 0, 0},
{3.9860043543609598*^14, 1.08262544999999997*^-3, 6378136.3},
{7.211454165826*^9, 0, 0},
{3.202738774922892*^12, 0, 0},
{9.887834453334144*^12, 0, 0},
{1.205120887033*^11, 0, 0},
{5.959916033410404*^12, 0, 0},
{1.266865349218008*^17, 1.46956199999999995*^-2, 71492000},
{4.282837362069909*^13, 1.95660915940866617*^-3, 0.3396000000000000*^7},
{2.2031780000000021*^13, 4.40443532482049831*^-6, 2.4400000000000000*^6},
{2.503524000000*^9, 0, 0},
{4.319516899232100*^9, 0, 0},
{4.9028000661637961*^12, 2.03215684649525711*^-4, 1738000},
{6.835099502439672*^15, 3.40843000000000005*^-3, 25225000},
{2.053234302535623*^11, 0, 0},
{7.087546066894452*^5, 0, 0},
{8.696138177608748*^11, 0, 0},
{1.539424643535*^11, 0, 0},
{3.793120749865224*^16, 1.62907099999999999*^-2, 60330000},
{1.3271244004193938*^20, 2.11060885327268404*^-7, 696000000},
{4.121107782641*^10, 0, 0},
{8.978138376543*^12, 0, 0},
{2.269437003741248*^11, 0, 0},
{1.427598140725034*^12, 0, 0},
{8.509338094489388*^10, 0, 0},
{5.793951322279009*^15, 3.51067999999999996*^-3, 25559000},
{3.2485859200000006*^14, 4.40443532482049831*^-6, 6.05100000000000*^6},
{0.1728824496930000*^11, 7.10608919544419154*^-2,0.2650000000000000*^6}};
constaxe = {"Sun", "Earth"};
referencepoint = "Earth";
planetcounts = Length[planetname];
allaxes = Table[Normalize[rawallaxes[[n]]], {n, 1, planetcounts, 1}];
input = iinput + Table[{{0, 0, 0}, If[planetname[[n]] == soibody, {0, 0, 0}, iinput[[Position[planetname, soibody][[1, 1]], 2]]]}, {n, 1, planetcounts, 1}];
set = Sum[planetdata[[n, 1]]*input[[n]], {n, 1, planetcounts, 1}]/Sum[planetdata[[n, 1]], {n, 1, planetcounts, 1}];
inti = input - Table[set, {n, 1, planetcounts, 1}];
positionvectors = Table[{Subscript[x, n, 1][t], Subscript[x, n, 2][t], Subscript[x, n, 3][t]}, {n, 1, planetcounts, 1}];
velocityvectors = Table[{Subscript[vx, n, 1][t], Subscript[vx, n, 2][t], Subscript[vx, n, 3][t]}, {n, 1, planetcounts, 1}];
norm[r_] := Sqrt[Sum[r[[n]]^2, {n, 1, Length[r], 1}]];
gravitationalparameterfunction[r1_, r0_, data_, axes_] := data[[1]]*(1 + data[[2]]*(data[[3]]/norm[r0 - r1])^2*(3 ((r1 - r0).axes/norm[r0 - r1])^2 - 1)/2);
V = -Sum[gravitationalparameterfunction[positionvectors[[n]], positionvectors[[m]], planetdata[[n]], allaxes[[n]]]*gravitationalparameterfunction[positionvectors[[n]], positionvectors[[m]], planetdata[[m]], allaxes[[m]]]/(G*norm[positionvectors[[n]] - positionvectors[[m]]]), {n, 1, planetcounts, 1}, {m, n + 1, planetcounts, 1}];
odes = Table[{Subscript[vx, n, m]'[t]*planetdata[[n, 1]]/G == -\!$$\*SubscriptBox[\(\[PartialD]$$, $$\(\*SubscriptBox[\(x$$, $$n, m$$]\)[t]\)]V\), Subscript[vx, n, m][t] == Subscript[x, n, m]'[t]}, {n, 1, planetcounts, 1}, {m, 1, 3, 1}];
bds = Table[{Subscript[vx, n, m][0] == inti[[n, 2, m]], Subscript[x, n, m][0] == inti[[n, 1, m]]}, {n, 1, planetcounts, 1}, {m, 1, 3, 1}];
equs = Flatten[{bds, odes}];
sol = NDSolve[equs, Flatten[{positionvectors, velocityvectors}], {t, 0, tmax} ];
evectors = {Table[Subscript[x, n, m][t], {n, 1, planetcounts, 1}, {m, 1, 3, 1}], Table[Subscript[vx, n, m][t], {n, 1, planetcounts, 1}, {m, 1, 3, 1}]} /. sol;
epositionvectors = evectors[[1, 1]];
evelocityvectors = evectors[[1, 2]];
DumpSave["planetorbit.mx", {epositionvectors, evelocityvectors}];


Maybe there are some bugs in version 11.0.I can run it faster in version 11.3.

• Please give complete data. equs, positionvectors, and velocityvectors are undefined. – Henrik Schumacher Aug 21 '18 at 12:21
• @HenrikSchumacher Ok,I added the program I used to generate 'planetorbit.mx' .It takes about several minutes to run. – YiYun Wang Aug 21 '18 at 12:50
• Can't reproduce the issue with tmax=tmax/100 in v9.0.1. Did you test the code only once? – xzczd Aug 21 '18 at 13:04
• Hm. I cannot reproduce that behavior with version 11.3 on macOS. ParametricPlot3D[Evaluate@epositionvectors, {t, 0, tmax}] plot all the curves within 1.2 seconds. – Henrik Schumacher Aug 21 '18 at 13:05
• Well, also Plot[Evaluate@epositionvectors[[4,2]], {t, 0, tmax}] is plotts alsmost immediately. Notice the Evaluate. – Henrik Schumacher Aug 21 '18 at 13:14