[Tutorial] Non-unimodal ternary search
Difference between en3 and en4, changed 190 character(s)
Hello everyone. Today I would like to introduce you to a new technique — ternary search on non-unimodal functions. I came up with it during ROI 2024 and it helped me to become a prize winner.↵

The motivation↵
------------------↵
There are many examples in CP where we need to maximise or minimise the value of some function $f$. If the function is unimodal, there is no problem using ternary search. But if it's not, there's nothing to do. Of course, if $f$ is random, we should check all its values. Fortunately, intuition tells us that if the function is often close to unimodal, we can also use ternary search on a large part of its area of definition.↵

### The intuition↵
Consider the unimodal function $0.03x^2$ after adding some noise ($sin (x)$) to it:↵

<spoiler summary="Graph">↵
![ ](/predownloaded/43/c2/43c2bf81135659543d66aba9c8802cd0d876161c.png)↵
</spoiler>↵

Intuition suggests that if the noise is not strong, we can still use a ternary search away from the global minimum.↵

In this example, we can put it more formally &mdash; the derivative of the function is $0.06x + cos(x)$, and when $|x| > \frac{1}{0.06} \approx 16.67$, adding cosine does not change the sign of the derivative.↵
#### The first optimisation↵
So, the first idea is to run the ternary search using not `while (r - l > eps)` but `while (r - l > C)` and then bruteforce all values between $l$ and $r$ with some precision. In many cases when $f$ takes an integer argument there will be no precision issue at all.↵
#### The second optimisation↵
I
 shouldt is mentioned in [this blog](https://codeforces.me/blog/entry/60702). It tells us a similar idea of splitting all argument values into blocks and applying ternary search on them.↵

This is the only thing related to the blog that I found. I tried googling and asking people involved in CP, but none of them knew about it before.↵

### Testing↵
The function from the example is boring, so consider a more interesting function: [Weierstrass function](https://www.desmos.com/calculator/lae5b4wdhq)↵

We zoom in and find that the maximum is about $1.162791$↵

<spoiler summary="Spoiler">↵
![ ](/predownloaded/ee/bb/eebb6dc44e4209d597b4d81fb23318cc51a2c16d.png)↵
</spoiler>↵

We will search for the maximum on the interval (-1, 1).↵

<spoiler summary="Code">↵

~~~~~↵
#include <bits/stdc++.h>↵
#define M_PI 3.14159265358979323846↵
using namespace std;↵
#define ld long double↵
ld Weierstrass(ld x) {↵
    ld y = 0;↵
    for (int i = 0; i < 30; i++) y += pow(0.14, i) * cos(pow(51, i) * M_PI * x);↵
    return y;↵
}↵
ld Classical_ternary(ld f(ld), ld l = -1, ld r = 1, ld eps = 1e-7) {↵
    while (r - l > eps) {↵
        ld midl = l + (r - l) / 3;↵
        ld midr = r - (r - l) / 3;↵
        if (f(midl) > f(midr)) r = midr;↵
        else l = midl;↵
    }↵
    return f(l);↵
}↵
int main() {↵
    cout << Classical_ternary(Weierstrass) << endl;↵
}↵
~~~~~↵

</spoiler>↵

This gives us $1.12881$. Changing $eps$ will change this value slightly.↵

Let's split the arguments into blocks. Since the arguments are real, we will not actually split them explicitly into blocks, we will take the maximum in some range near the argument.↵

<spoiler summary="Blocks">↵

~~~~~↵
ld Blocks(ld f(ld), ld l = -1, ld r = 1, ld eps = 1e-7) {↵
    int iter = 10000;↵
    while (r - l > eps) {↵
        ld midl = l + (r - l) / 3;↵
        ld midr = r - (r - l) / 3;↵
        ld vall = -2e9, valr = -2e9;↵
        for (ld x = midl, i = 0; i < iter; x += eps, i++) vall = max(vall, f(x));↵
        for (ld x = midr, i = 0; i < iter; x += eps, i++) valr = max(valr, f(x));↵
        if (vall > valr) r = midr;↵
        else l = midl;↵
    }↵
    return f(l);↵
}↵
~~~~~↵
</spoiler>↵

It gives $1.15616$, which is quite good. We can optimise it by taking the maximum of all the values of $f$ we have ever calculated:↵

<spoiler summary="Take maximum">↵
~~~~~↵
ld Blocks(ld f(ld), ld l = -1, ld r = 1, ld eps = 1e-7) {↵
    int iter = 10000;↵
    ld ans = -2e9;↵
    while (r &mdash; l > eps) {↵
        ld midl = l + (r &mdash; l) / 3;↵
        ld midr = r &mdash; (r &mdash; l) / 3;↵
        ld vall = -2e9, valr = -2e9;↵
        for (ld x = midl, i = 0; i < iter; x += eps, i++) vall = max(vall, f(x));↵
        for (ld x = midr, i = 0; i < iter; x += eps, i++) valr = max(valr, f(x));↵
        ans = max({ans, vall, valr});↵
        if (vall > valr) r = midr;↵
        else l = midl;↵
    }↵
    return max(ans, f(l));↵
}↵
~~~~~↵
</spoiler>↵

This gives us $1.16278$, which is very close to the expected $1.162791$. It seems that we have found a good approach. But let's go further.↵

##The third optimisation↵
Let's change the constant $3$ in the code. We will call it 
C$C$. It is not new to experienced people, often it is good to choose C$C$ equal to 2$2$ (binary search by the derivative) or $\frac{\sqrt5+1}{2}$ (golden ratio).↵
As we are cutting out $\frac{1}{C}$ part of our interval on each iteration, the probability of the maximum being indise the cut-off part decreases as C increases.↵

Let's choose 
$C = 100$ and run ternary search. As we already know, taking maximum of all function calls is very good optimisation, so I have already added it.↵

<spoiler summary="Code">↵
~~~~~↵
ld find_maximum(ld f(ld), int C = 3, ld l = -1, ld r = 1, ld eps = 1e-7) {↵
    ld ans = -2e9;↵
    while (r &mdash; l > eps) {↵
        ld midl = l + (r &mdash; l) / C;↵
        ld midr = r &mdash; (r &mdash; l) / C;↵
        ld vall = f(midl), valr = f(midr);↵
        ans = max({ans, vall, valr});↵
        if (vall > valr) r = midr;↵
        else l = midl;↵
    }↵
    return max(ans, f(l));↵
}↵
~~~~~↵
</spoiler>↵

If we run it with $C = 3$, we get 
$1.13140$ (the algo is the same as the classic ternary search, but we take the maximum, so the answer is much better).↵
Now let's now increase $C$ and watch the answer increases:↵

We run it with $C = 30$ and get $1.16275$.↵
We run it with $C = 100$ and get ... $1.15444$↵

In fact, increasing $C$ does not guarantee a better answer.↵
### The fourth optimisation↵
Let's bruteforce all values of $C$ from 3 to 100 and take the best answer:↵

`for (int C = 3; C < 100; C++) res = max(res, find_maximum(Weierstrass, C));`↵

It gives us $1.16279$ and works faster than block splitting. To compare them further, we need to change the function, because both methods return values that are very close to the answer.↵

Let's use $a = 0.88$ and $b = 51$ for the Weierstrass function. Note that it is impossible to see the actual maximum of the function in Desmos.↵

I will compare the above 2 approaches on Codeforces Custom test.↵

<spoiler summary="C < 100">↵

~~~~~↵
Blocks res: 6.6809↵
Blocks time: 3.645↵
100-ternary res: 6.27581↵
100-ternary time: 0.753↵
~~~~~↵
</spoiler>↵


<spoiler summary="C < 200">↵

~~~~~↵
Blocks res: 6.6809↵
Blocks time: 3.365↵
200-ternary res: 6.48924↵
200-ternary time: 2.725↵
~~~~~↵


</spoiler>↵


<spoiler summary="C < 400">↵

~~~~~↵
Blocks res: 6.6809↵
Blocks time: 3.442↵
400-ternary res: 7.05544↵
400-ternary time: 10.751↵
~~~~~↵
</spoiler>↵

I tried combining these methods together, but is performs worse than both techniques (because I used $1000$ iterations instead of $10000$ and bruteforced $C < 10$ not to run out of time).↵

Conclusion↵
------------------↵
On the function I considered, both approaches are good enough. On real contests I only used my approach.↵

From recent problems, [problem:2035E], I have used this technique successfully
: [submission:288346163].. Here is the implementation for integer arguments: [submission:288346163].↵

If you know how to find the maximum better (not by simulated annealing), please write.↵

History

 
 
 
 
Revisions
 
 
  Rev. Lang. By When Δ Comment
en8 English polosatic 2024-11-04 15:28:13 26 Tiny change: 'y approach.\n\nFrom ' -> 'y approach (the fourth optimisation).\n\nFrom '
ru1 Russian polosatic 2024-11-04 15:27:30 7676 Первая редакция перевода на Русский
en7 English polosatic 2024-11-04 15:24:31 40
en6 English polosatic 2024-11-04 15:23:43 24
en5 English polosatic 2024-11-04 14:59:46 6 Tiny change: 'nction is often close to ' -> 'nction is close to '
en4 English polosatic 2024-11-04 14:54:37 190 Tiny change: 'isation\nI should mention [this bl' -> 'isation\nIt is mentioned in [this bl' (published)
en3 English polosatic 2024-11-04 14:22:09 2719 Tiny change: 's, C));`\nIt gives' -> 's, C));`\n\nIt gives'
en2 English polosatic 2024-11-04 12:50:52 1362 Tiny change: 'the worst.' -> 'the worst.\n\nThe third optimisation ------------------'
en1 English polosatic 2024-11-04 12:09:48 3911 Initial revision (saved to drafts)