Ask Your Question
1

Finding the maximum on an interval and find root

asked 2016-09-03 21:26:33 +0100

vukov gravatar image

I have a function f that is the truncation of the fourier transform of a sequence (i.e. f = sum( cos(a*x) for a in sequence), where sequence has 10,000 terms). I want to find the absolute maximum on an interval.

The documentation suggest getting the maximum of the plot, which works. I would like to know where this maximum occurs, though. As the maximum is 3272.9689885743473, f - 3200 should have a root, but find_root doesn't find a root. I guess this is technically a bug, but not very surprising given the definition of f.

My main question is, how do I find out where this maximum occurs? I can squint at the plot, but this isn't accurate enough.

sequence = [3, 4, 6, 9, 10, 17, 18, 25, 30, 32, 37, 44, 45, 52, 59, 60, 71, 72, 73, 86, 87, 94, 99, 106, 107, 114, 121, 122, 134, 135, 142, 149, 156, 161, 163, 168, 175, 176, 183, 190, 197, 204, 211, 218, 225, 226, 237, 238, 245, 252, 253, 265, 266, 273, 280, 287, 288, 300, 301, 314, 315, 322, 327, 334, 335, 342, 349, 350, 362, 363, 370, 377, 384, 391, 396, 403, 404, 411, 418, 419, 431, 432, 439, 446, 453, 454, 465, 466, 473, 480, 481, 493, 494, 501, 508, 515, 516, 528, 529, 542, 543, 550, 555, 562, 563, 570, 577, 584, 591, 598, 605, 612, 619, 620, 631, 632, 639, 646, 647, 659, 660, 667, 674, 681, 682, 694, 695, 708, 709, 721, 722, 729, 736, 743, 744, 756, 757, 764, 771, 772, 784, 785, 798, 805, 812, 813, 825, 826, 833, 840, 841, 848, 853, 860, 867, 874, 875, 887, 888, 895, 902, 903, 910, 915, 922, 923, 930, 937, 944, 949, 956, 957, 964, 971, 972, 984, 985, 992, 999, 1006, 1013, 1018, 1026, 1033, 1040, 1041, 1053, 1054, 1061, 1068, 1069, 1080, 1081, 1082, 1095, 1102, 1103, 1115, 1116, 1123, 1130, 1137, 1138, 1150, 1157, 1158, 1165, 1172, 1177, 1184, 1185, 1192, 1199, 1200, 1212, 1213, 1220, 1227, 1234, 1235, 1246, 1247, 1254, 1261, 1262, 1269, 1274, 1281, 1282, 1289, 1296, 1297, 1309, 1310, 1323, 1324, 1331, 1336, 1343, 1344, 1351, 1358, 1365, 1372, 1379, 1386, 1393, 1400, 1401, 1412, 1413, 1420, 1427, 1428, 1440, 1441, 1448, 1455, 1462, 1463, 1475, 1482, 1489, 1490, 1502, 1503, 1510, 1517, 1524, 1525, 1537, 1538, 1551, 1552, 1559, 1566, 1571, 1578, 1579, 1586, 1593, 1594, 1606, 1607, 1614, 1621, 1628, 1629, 1640, 1641, 1648, 1655, 1656, 1668, 1669, 1676, 1683, 1690, 1691, 1703, 1704, 1717, 1718, 1725, 1730, 1737, 1738, 1745, 1752, 1759, 1766, 1773, 1780, 1787, 1794, 1795, 1806, 1807, 1814, 1821, 1822, 1834, 1835, 1842, 1849, 1856, 1857, 1869, 1870, 1883, 1884, 1896, 1897, 1904, 1911, 1918, 1919, 1931, 1932, 1939, 1946, 1953, 1960, 1965, 1973, 1980, 1987, 1988, 2000, 2001, 2008, 2015, 2022, 2023, 2035, 2042, 2049, 2050, 2062, 2063, 2070, 2077, 2078, 2085, 2090, 2097, 2104, 2111, 2112, 2119, 2124, 2131, 2132, 2139, 2146 ...
(more)
edit retag flag offensive close merge delete

1 Answer

Sort by ยป oldest newest most voted
1

answered 2016-09-04 19:50:13 +0100

nbruin gravatar image

I'm pretty sure your function attains its maximum at x=0, where each of the terms will have its maximal (positive) value 1. I assume your series arises as a truncation of an infinite series? It doesn't look like this series will be pointwise convergent, so I'm not so sure that the extremal values of the truncated series will tell you something useful.

In general, I'd expect the algorithms in scipy for optimizing functions to be pretty good. Use for instance:

sage: import scipy.optimize as opt sage: f=fast_callable(sum( cos(a*x) for a in sequence),vars=(x,),domain=RDF) sage: opt.minimize_scalar(f)

although you might want to double check if the numerics are trustworthy. Your expression has an awful lot of terms, so there's plenty of room for precision loss.

edit flag offensive delete link more

Comments

scipy is giving me the wrong answer. Also, I am looking for the x-value of the giant spike that occurs at about 5.56, not the maximum on [0, 2pi]. (There are many local maximum around there, so find local maximum doesn't work.)

vukov gravatar imagevukov ( 2016-09-04 23:25:14 +0100 )edit

Scipy's optimize routines also allow for bracketing. If you give enough information, I would expect that one of the algorithms implemented in scipy will do the job. If the precision issues are so bad that the methods don't work, you may have to look into using higher precision (i.e., RealField(100) or so), in whcih case you probably can't use the scipy inplementations anymore.

nbruin gravatar imagenbruin ( 2016-09-05 03:06:15 +0100 )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: 2016-09-03 21:26:33 +0100

Seen: 1,192 times

Last updated: Sep 04 '16