Setup. Suppose we are given a function $$f:\mathbb N\to\{\text{False},\text{True}\}$$ such that $$f(n)=\text{True}\implies f(n+1)=\text{True}$$ and such that $f(n)=\text{True}$ for some $n$ large enough.

In natural language. The function $f$ imposes a condition on the natural numbers which is fulfilled once $n$ is large enough.

My question 1. How can I, for any given $f$, find the smallest $n$ such that $f(n)=\text{True}$?

A first idea would be to start with $n=1$ and to increment $n$ by one until $f(n)$ is True. However, this is fairly slow. The next step is a "binary search" algorithm (see below).

My question 2. How can I implement a good algorithm in Mathematica?

Here is an example of the first algorithm that I thought of, implemented in Python:

def bin_search(cond):
    n = 1
    while not cond(n):
        n *= 2
    lower_bound = n//2
    upper_bound = n
    middle = (lower_bound + upper_bound)//2
    while upper_bound - lower_bound > 1:
        if cond(middle):
            upper_bound = middle
            lower_bound = middle
        middle = (lower_bound + upper_bound)//2
    return upper_bound

For example, one such condition would be $$f(n)=[H_n\geq 10],$$

where $$H_n=\sum_{i=1}^n \frac 1i$$ is the $n$th harmonic number.

EDIT: A remark about the harmonic numbers: In fact I wanted to compare the runtime of a good brute force method against an intelligent method for this problem. For solving $H_n\geq 10$, it suffices to find the smallest integer $n$ such that $$\ln(n)+\gamma\geq 10,$$ where $\gamma$ denotes the Euler-Mascheroni constant. So we can find $$n=[ \exp(10-\gamma)],$$ where $[\cdot]$ denotes the nearest integer function.

  • 2
    $\begingroup$ 1. You are asking two different questions: (Q1) an algorithm derivation one, and (Q2) a code translation one. Which one do you want? 2. If Q1, then there is no need for ugly Python code... :) $\endgroup$ – Anton Antonov Dec 3 '19 at 13:26
  • $\begingroup$ @AntonAntonov Since I am looking at general $f$, my question is: What is a good implementation in Mathematica of a good algorithm which is able to find the smallest $n$ such that ...? I don't want to use any specific properties of $f$ except the two mentioned above... My Python code was just an example of the simplest algorithm I could come up with, namely, a sort of binary search in order to find the minimal $n$ $\endgroup$ – Maximilian Janisch Dec 3 '19 at 13:32
  • $\begingroup$ @AntonAntonov By the way, what is wrong with the Python implementation? 😯 $\endgroup$ – Maximilian Janisch Dec 3 '19 at 18:59
  • $\begingroup$ You are aware, of course, that the problem has no solution in general, right? One can always construct a sequence growing slower than the $2^n$ built in the algorithm. $\endgroup$ – yarchik Dec 3 '19 at 19:15
  • 2
    $\begingroup$ @yarchik Actually, based on your Input, I will consider asking a part of this question on the StackExchange for theoretical Computer Science. That might reduce the confusion $\endgroup$ – Maximilian Janisch Dec 3 '19 at 19:34

Here is one implementation of your algorithm:

firstTrue[f_] := next[f, 1, Infinity]

next[f_, last_, cur_] := If[last + 1 == cur,
    If[f[last], last, cur],
    If[f[last], next[f, dec[last], last], next[f, inc[last], cur]]

inc[n_] := With[{e = IntegerExponent[n, 2]},
    If[n == 2^e, 2^(e+1), n + 2^(e-1)]

dec[n_] := With[{e = IntegerExponent[n, 2]},
    If[n == 2^e, 2^(e-1) + 2^(e-2), n - 2^(e-1)]

Your example:

firstTrue[HarmonicNumber[#] >= 10&]


Of course, for analytic functions like HarmonicNumber there are much better algorithms to do this.

  • $\begingroup$ Thank you! About HarmonicNumber: I will add some motivation; of course you are right $\endgroup$ – Maximilian Janisch Dec 3 '19 at 18:20
  • $\begingroup$ Ok I have added some motivation to the Harmonic Number problem. In my opinion it is a very nice solution using Euler-Mascheroni $\endgroup$ – Maximilian Janisch Dec 3 '19 at 18:23

As Carl Woll says:

Of course, for analytic functions like HarmonicNumber there are much better algorithms to do this.

Here is some code related to that statement -- it uses one of the universal WL functions:

H[n_] := Sum[1/i, {i, n}];

FindArgMin[{H[n], H[n] >= 10, n > 0}, n]
(* {12366.5} *)

(* 10. *)

(* 9.99996 *)
  • $\begingroup$ This is interesting because it simply uses numerical methods. A caveat is the innacuracy: For example, FindArgMin[{H[n], H[n] >= 30, n > 0}, n, WorkingPrecision -> 1000, AccuracyGoal -> 1000, PrecisionGoal -> 1000] returns $$6.000022445096458984375\cdot 10^{12}$$ and the correct result is $$6000022499693$$ But it is indeed a very quick way to get close to the correct solution $\endgroup$ – Maximilian Janisch Dec 3 '19 at 18:55

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.