SERVING THE QUANTITATIVE FINANCE COMMUNITY

Cuchulainn
Posts: 54228
Joined: July 16th, 2004, 7:38 am
Location: Amsterdam
Contact:

### Re: exp(5) = $e^5$

Instead of Newton, bisection, analytical approximations (too many assumptions that can break during stormy weather) we do least-squares minimisation of the unimodal(!) function $(log(x) - 5.0)^2$  using the (most?) efficient three-point interval search on [32,243] with the least number of function evaluations. Does it do better than other algos for iv?

Here is the output (tol = 1.0e-5). Nice thing is the answer always lies in an interval.

Minimum: x, f(x) [84.75,190.25], 0.187801
Minimum: x, f(x) [111.125,163.875], 0.0467709
Minimum: x, f(x) [137.5,163.875], 0.00782745
Minimum: x, f(x) [144.094,157.281], 0.00212025
Minimum: x, f(x) [144.094,150.688], 0.000551828
Minimum: x, f(x) [147.391,150.688], 0.000139543
Minimum: x, f(x) [147.391,149.039], 3.27546e-05
Minimum: x, f(x) [147.803,148.627], 9.52942e-06
Minimum: x, f(x) [148.215,148.627], 1.93004e-06
Minimum: x, f(x) [148.318,148.524], 4.84548e-07
Minimum: x, f(x) [148.369,148.472], 1.23176e-07
Minimum: x, f(x) [148.395,148.447], 3.28335e-08
Minimum: x, f(x) [148.395,148.421], 8.72959e-09
Minimum: x, f(x) [148.408,148.421], 1.95913e-09
Minimum: x, f(x) [148.411,148.418], 5.47337e-10
Minimum: x, f(x) [148.411,148.414], 1.22007e-10
Minimum: x, f(x) [148.412,148.414], 3.37671e-11
Minimum: x, f(x) [148.413,148.414], 7.74401e-12
Minimum: x, f(x) [148.413,148.413], 2.22907e-12
Minimum: x, f(x) [148.413,148.413], 4.62369e-13
Minimum: x, f(x) [148.413,148.413], 1.17684e-13
Minimum: x, f(x) [148.413,148.413], 3.15131e-14
Minimum: x, f(x) [148.413,148.413], 8.20017e-15
Minimum: x, f(x) [148.413,148.413], 1.90432e-15
Bracket  [148.413,148.413]
Bracket [2.78939e-15,1.01925e-15]

For OP, tol = 1.0e-1 gives

Minimum: x, f(x) [84.75,190.25], 0.187801
Minimum: x, f(x) [111.125,163.875], 0.0467709
Minimum: x, f(x) [137.5,163.875], 0.00782745
Minimum: x, f(x) [144.094,157.281], 0.00212025
Minimum: x, f(x) [144.094,150.688], 0.000551828
Minimum: x, f(x) [147.391,150.688], 0.000139543
Minimum: x, f(x) [147.391,149.039], 3.27546e-05
Minimum: x, f(x) [147.803,148.627], 9.52942e-06
Minimum: x, f(x) [148.215,148.627], 1.93004e-06
Minimum: x, f(x) [148.318,148.524], 4.84548e-07
Minimum: x, f(x) [148.369,148.472], 1.23176e-07
Bracket  [148.369,148.472]
Bracket [8.70207e-08,1.59332e-07]
http://www.datasimfinancial.com

What is 2^77232917 -1?

Cuchulainn
Posts: 54228
Joined: July 16th, 2004, 7:38 am
Location: Amsterdam
Contact:

### Re: exp(5) = $e^5$

We all know by now that the answer is 148.413 and that the function to be minimised is simple. For more complicated functions things can become more tricky and then we want to make C++ look more like mathematics (functional programming languages do it before breakfast). C++11 helps a bit by lambda functions so that you can build a vector space of functions to assemble functions from simpler ones. You would kind of expect that a library builder out there has done this (knowledge of maths is needed).
// TestHOF101.cpp
//
// Testing lambda functions; build complex functions using algebra.
// A vector space of functions.
// In  this case compute x - e^5 = 0 using least squares and Brent's method
// x == 148.413
//
// DD

#include <functional>
#include <cmath>
#include <iostream>
#include <boost/math/tools/minima.hpp>

template <typename T>
std::function<T (T)> operator * (std::function<T (T)>& f, std::function<T (T)>& g)

return [=](T x)
{
return  f(x) * g(x);
};
}

template <typename T>
std::function<T (T)> operator - (std::function<T (T)>& f, T a)

return [=](T x)
{
return  f(x) - a;
};
}

template <typename T>
std::function<T (T)> log(std::function<T (T)>& f)
{ // log

return [=](T x)
{
return  std::log(f(x));
};
}
template <typename T>
std::function<T(T)> exp(std::function<T(T)>& f)
{ // log

return [=](T x)
{
return  std::exp(f(x));
};
}

double funcGlobal(double x)
{
return (std::log(x) - 5.0)* (std::log(x) - 5.0);
}

int main()
{

// Define interval [min, max] where search will take place
double min = 5.0;
double max = 300.0;

int bits = 50;
boost::uintmax_t maxIter = 10000;

// Standard function call
std::pair<double, double> result = boost::math::tools::brent_find_minima(funcGlobal, min, max, bits, maxIter);
std::cout << "Abscissa, value f(x) : " << result.first << ", " << result.second << std::endl;

// Now use lambda and algebra; assemble functions incrementally
std::function<double(double)> f1 = [](double x) { return std::log(x); };
decltype(f1) f2 = f1 - 5.0;
decltype(f2) funcHOF = f2*f2;

// Stress testing
std::function<double(double)> f0 = [](double x) { return x; };
decltype(f0) funcExtreme = exp(log(exp(log(f0))));
std::cout << "? " << funcExtreme(5) << '\n';
decltype(f0) f3 = log(funcExtreme);
decltype(f0) f4 = exp(funcExtreme);
std::cout << "? " << f3(5) << '\n';
std::cout << "? " << f4(5) << '\n';
//

std::pair<double, double> result2 = boost::math::tools::brent_find_minima(funcHOF, min, max, bits, maxIter);
std::cout << "Abscissa, value f(x) HOF: " << result2.first << ", " << result2.second << std::endl;

}

//
http://www.datasimfinancial.com

What is 2^77232917 -1?

Cuchulainn
Posts: 54228
Joined: July 16th, 2004, 7:38 am
Location: Amsterdam
Contact:

### Re: exp(5) = $e^5$

Until now, we used double precision data types for the minimisation using interval searching. Will the algorithm work with multiprecision types? It would seem yes at the moment. We use 100-digit type and we get 148.413..... to 100 digits that we checked against Boost and that big online calculator.

And the run-time performance doesn't seem to be all that bad. No test done.

The output is

3-point interval Search method
148.4131 5910 2576 6034 2111 5580 0405 5227 9623 4876 6759 3878 9890 4675 2845 1109 1206 4820 9585 7607 9688 4094 5989 90211

Exact (yuge precision)
148.4131 5910 2576 6034 2111 5580 0405 5227 9623 4876 6759 3878 9890 4675 2845 1109 1206 4820 9585 7607 9688 4094 5989 9021 1

\-> etc.412 9280 8270 6663 2605 2909 9262 3769 4548 1827 8044 9112 4354 2477 7223 6978 7450 8409 7439 3799 2123 5516 6686 3241 8855 9697 6672 2232 1326 5623 4323 8300 9773 3581 8269

http://www.datasimfinancial.com

What is 2^77232917 -1?

Cuchulainn
Posts: 54228
Joined: July 16th, 2004, 7:38 am
Location: Amsterdam
Contact:

### Re: exp(5) = $e^5$

Paul wrote:
What if the entire notion that a transition probability even exists is false?

If you'd read PWOQF2 then you'd know that people might quote transition matrices for a period of a year for which there is no infinitesimal-period matrix, because you can't invert e^M (and by "invert" I don't mean one "one over"!). The implication being either the transition matrix being quoted is rubbish or that whoever generates the matrix is making some sophisticated assumption about time dependence. (No, I don't believe the latter either!) Yet another reason why the proper defintion of exp(M) is important and the other one is pretty pointless!
P

One way to persuade infidels is to define $e^A$ as the solution of a matrix ODE.

$dB/dt = AB, B(0) = I$

The solution is $B(t) = e^{At}$. So plug in  $t = 1$ and Bob's your uncle.

Now, the next step is to approximate the matrix ODE by a numeric ODE solver, easy in C++.
http://www.datasimfinancial.com

What is 2^77232917 -1?

Cuchulainn
Posts: 54228
Joined: July 16th, 2004, 7:38 am
Location: Amsterdam
Contact:

### Re: exp(5) = $e^5$

Here's the C++ code in Boost odeint. Few lines of code.

using value_type = double;
typedef boost::numeric::ublas::matrix<value_type> state_type;

class MatrixOde
{
private:
// dB/dt = A*B, B(0) = C;
boost::numeric::ublas::matrix<value_type> A_;
boost::numeric::ublas::matrix<value_type> C_;
public:

MatrixOde(const boost::numeric::ublas::matrix<value_type>& A, const boost::numeric::ublas::matrix<value_type>& IC)
: A_(A), C_(IC) {}

void operator()( const state_type &x , state_type &dxdt , double t ) const
{

for( std::size_t i=0 ; i < x.size1();++i )
{
for( std::size_t j=0 ; j < x.size2(); ++j )
{
dxdt(i, j) = 0.0;
for (std::size_t k = 0; k < x.size2(); ++k)
{
dxdt(i, j) += A_(i,k)*x(k,j);
}
}
}
}

};
http://www.datasimfinancial.com

What is 2^77232917 -1?

JGJPR
Posts: 12
Joined: September 7th, 2017, 3:33 pm

### Re: exp(5) = $e^5$

Is a quant a calculation engine?
Very strange question to me...

Cuchulainn
Posts: 54228
Joined: July 16th, 2004, 7:38 am
Location: Amsterdam
Contact:

### Re: exp(5) = $e^5$

JGJPR wrote:
Is a quant a calculation engine?
Very strange question to me...

Good point. Some want closed-form solution.
http://www.datasimfinancial.com

What is 2^77232917 -1?

Cuchulainn
Posts: 54228
Joined: July 16th, 2004, 7:38 am
Location: Amsterdam
Contact:

### Re: exp(5) = $e^5$

Here is another solution; all you need is to use your office wall as a dart board.

$\int_{0}^{5} \exp(x) dx = \exp(5) -1$

Then use Monte Carlo Hit or Miss Integration. Draw the curve and throw darts at the wall. How many hits? Maybe get more people to join in? Just be careful when throwing those darts!
http://www.datasimfinancial.com

What is 2^77232917 -1?

Cuchulainn
Posts: 54228
Joined: July 16th, 2004, 7:38 am
Location: Amsterdam
Contact:

### Re: exp(5) = $e^5$

The only solution left is tie solve this problem using (a hand-crafted?) ML algorithm?

Can I just look at the diagram and give the answer?

http://www.datasimfinancial.com

What is 2^77232917 -1?

Posts: 23371
Joined: September 20th, 2002, 8:30 pm

### Re: exp(5) = $e^5$

"Paper and pencil" would let us construct our own slide rule either for multiplication or for direct exponentiation. But it seems like a labor-intensive solution in that: 1) it requires assigning accurate multi-digit decimal values to specific marks on the slide rule; 2) it seems to imply that we must solve the more general problem of $e^x$. It also might take a lot of paper if that if the marks on the paper are accurate to only O(1 mm), then we might need a piece of paper O(14841 mm) long. We can use a roll of adding machine paper to do this.

Cuchulainn
Posts: 54228
Joined: July 16th, 2004, 7:38 am
Location: Amsterdam
Contact:

### Re: exp(5) = $e^5$

Cuchulainn wrote:
Here's the C++ code in Boost odeint. Few lines of code.

using value_type = double;
typedef boost::numeric::ublas::matrix<value_type> state_type;

class MatrixOde
{
private:
// dB/dt = A*B, B(0) = C;
boost::numeric::ublas::matrix<value_type> A_;
boost::numeric::ublas::matrix<value_type> C_;
public:

MatrixOde(const boost::numeric::ublas::matrix<value_type>& A, const boost::numeric::ublas::matrix<value_type>& IC)
: A_(A), C_(IC) {}

void operator()( const state_type &x , state_type &dxdt , double t ) const
{

for( std::size_t i=0 ; i < x.size1();++i )
{
for( std::size_t j=0 ; j < x.size2(); ++j )
{
dxdt(i, j) = 0.0;
for (std::size_t k = 0; k < x.size2(); ++k)
{
dxdt(i, j) += A_(i,k)*x(k,j);
}
}
}
}

};

A sanity check is to take the scalar case du/dt = u. I think we all agree that a scalar can be seen as a 1x1 matrix. So, nothing shocking really but it is OP.
148.413.
// Sanity check; 1X1 matrix case dB/t = B on (0,5), B(0) = I
{ // dB/dt = A*B, B(0) = C

std::size_t N = 1;
ublas::identity_matrix<value_type> C(N);
ublas::matrix<value_type> A(N, N);
A(0, 0) = 1.0;
ublas::matrix<value_type> B(N, N);

// Initialise the solution matrix B
B = C;

MatrixOde ode(A, C);
double L = 0.0;
double T = 5.0;
double dt = 0.001;

std::size_t steps = Bode::integrate(ode, B, L, T, dt, write_out);

std::cout << "Value at time: " << T << '\n';
std::cout << B << '\n';
}

http://www.datasimfinancial.com

What is 2^77232917 -1?

Posts: 23371
Joined: September 20th, 2002, 8:30 pm

Put $1 in a bank account yielding exactly 2% continuously compounded and then check the account balance at noon in 252 years and 179 days. Cuchulainn Posts: 54228 Joined: July 16th, 2004, 7:38 am Location: Amsterdam Contact: ### Re: exp(5) = $e^5$ One glaring omission is that no one has mentioned how $e$ came to be. Everyone blames Euler (some blame APL) but its roots are probably due to a 17th century banker or trader who stumbled on compound interest: $\displaystyle\lim_{n\to\infty} (1 + r/n)^{nt} \rightarrow e^{rt}$ Then the special case gives us $\displaystyle\lim_{n\to\infty} (1 + 1/n)^{n} == e$ There you have it. Even Einstein agrees on $e^5$. Last edited by Cuchulainn on November 5th, 2017, 11:37 am http://www.datasimfinancial.com What is 2^77232917 -1? Cuchulainn Posts: 54228 Joined: July 16th, 2004, 7:38 am Location: Amsterdam Contact: ### Re: exp(5) = $e^5$ Traden4Alpha wrote: Put$1 in a bank account yielding exactly 2% continuously compounded and then check the account balance at noon in 252 years and 179 days.

Very good answer. I had not read this until this morning. But funnily the same idea crossed my mind yesterday.

So, the question was: How long does it take to earn $148.413 by investing$1 at a fixed rate of 2%?
http://www.datasimfinancial.com

What is 2^77232917 -1?

Posts: 23371
Joined: September 20th, 2002, 8:30 pm

### Re: exp(5) = $e^5$

Cuchulainn wrote:
Put $1 in a bank account yielding exactly 2% continuously compounded and then check the account balance at noon in 252 years and 179 days. Very good answer. I had not read this until this morning. But funnily the same idea crossed my mind yesterday. So, the question was: How long does it take to earn$148.413 by investing \$1 at a fixed rate of 2%?
Yep!

I also pondered if other exponential processes might be adapted to answer this question. But base-2 doubling requires an inconvenient 5/ln(2) = 7.213475 doublings to get the answer. And looking for the right seed for a Fibonacci process seemed to yield an answer of starting with 7457 and 12067 which after 10 cycles reaches 10000*e^5 with requisite accuracy but it feel like a hack!

Wilmott.com has been "Serving the Quantitative Finance Community" since 2001. Continued...

 JOBS BOARD

Looking for a quant job, risk, algo trading,...? Browse jobs here...