SERVING THE QUANTITATIVE FINANCE COMMUNITY

• 1
• 2
• 3
• 4
• 5
• 17

outrun
Posts: 4573
Joined: April 29th, 2016, 1:40 pm

### Re: Looking for hardware recommendations

Does this look good?

1) Bessel I function with complex args line 458 and text below, from the C source code used in the R package Bessel (pdf)
2) the confluent hypergeometric is more difficult to find! I found a power series sum routine

I know too little about these functions..

Alan
Topic Author
Posts: 10015
Joined: December 19th, 2001, 4:01 am
Location: California
Contact:

### Re: Looking for hardware recommendations

Re a C/C++ version, I thought a little about that several days ago. The main bottleneck there is I need two special functions with all complex arguments: the Bessel I(nu,z) and the confluent hypergeometric M(a,b,z). Googling, I could only locate one library supporting those and it was not really set up for Windows, my current environment. If anybody knows such a library that can easily be hooked in to Visual Studio on Windows, I'd appreciate a link. (Some usual suspects, Boost and GSL, don't have the complex argument support, as far as I could tell).
Yeah, typically a Fortran era function that has not (yet) been ported to C++.
So, given these functions should be available in Fortran, perhaps I should try this:
https://en.wikipedia.org/wiki/Intel_Fortran_Compiler

outrun
Posts: 4573
Joined: April 29th, 2016, 1:40 pm

### Re: Looking for hardware recommendations

That's a great idea.

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

### Re: Looking for hardware recommendations

Maybe TOMS (NSWC)
http://jblevins.org/mirror/amiller/

Now, Code Block for FORTAN is free

I actually wrote the BS formula in it just to see I could remember my Fortran 66/77 from some 35 years earlier.
It worked!
ACHTUNG: don't go beyond line 72 because it annoys the compiler
N.B. And as far as I remember install was kinda easy.
http://www.datasimfinancial.com
http://www.datasim.nl

Every Time We Teach a Child Something, We Keep Him from Inventing It Himself
Jean Piaget

Alan
Topic Author
Posts: 10015
Joined: December 19th, 2001, 4:01 am
Location: California
Contact:

### Re: Looking for hardware recommendations

Yes, I have downloaded TOMS 707 and will see if I can get it compiled after I get done with my Mathematica run.

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

### Re: Looking for hardware recommendations

linkable to C/C++ code compiled under Visual Studio on Windows?  (It seems to have many dependencies on unix-type stuff).

It's unfortunate that many otherwise useful software libraries fall at this last (configuration) hurdle.
Boost is nice because it is mostly header-only.
http://www.datasimfinancial.com
http://www.datasim.nl

Every Time We Teach a Child Something, We Keep Him from Inventing It Himself
Jean Piaget

outrun
Posts: 4573
Joined: April 29th, 2016, 1:40 pm

### Re: Looking for hardware recommendations

For Linux users it's fine.

Boost is also fine for windows users because it's cross platform and serves windows too ...but a lot of scientific tools don't bother with windows, being cross platform gives a lot of overhead.
The solution is simple: run a Linux VM, it'll save you a lot of trouble.

outrun
Posts: 4573
Joined: April 29th, 2016, 1:40 pm

### Re: Looking for hardware recommendations

Not sure if "setting up a Debian Linux VM" looks complicated or not?

https://williamjturkel.net/2013/05/31/i ... -a-vm/amp/

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

### Re: Looking for hardware recommendations

This article, eq (4) shows CHF as an absolutely convergent series, so now that C++11 and Boost support the needed building blocks I thought it could be programmed on a lazy Sunday afternoon. I tested against A&S, C++ std::exp(z) and the famous $e^5$. It's about 6 lines of code. I just tried it to kick the tyres. Either it works or it doesn't. And it seems feasible.
Larger values if z ==> more terms in (4) and you might also need multiprecision at some stage.
https://dspace.mit.edu/bitstream/handle ... sequence=1

Have I produced the first mainstream/accessible results in C++ for CHF?

I'll post on $e^5$ after I review + take some more tests.
http://www.datasimfinancial.com
http://www.datasim.nl

Every Time We Teach a Child Something, We Keep Him from Inventing It Himself
Jean Piaget

Alan
Topic Author
Posts: 10015
Joined: December 19th, 2001, 4:01 am
Location: California
Contact:

### Re: Looking for hardware recommendations

This article, eq (4) shows CHF as an absolutely convergent series, so now that C++11 and Boost support the needed building blocks I thought it could be programmed on a lazy Sunday afternoon. I tested against A&S, C++ std::exp(z) and the famous $e^5$. It's about 6 lines of code. I just tried it to kick the tyres. Either it works or it doesn't. And it seems feasible.
Larger values if z ==> more terms in (4) and you might also need multiprecision at some stage.
https://dspace.mit.edu/bitstream/handle ... sequence=1

Have I produced the first mainstream/accessible results in C++ for CHF?

I'll post on $e^5$ after I review + take some more tests.
Thanks for experimenting. Here's a paper with a few tests you might want to try (Table 1, 3):
http://www.cs.upc.edu/~argimiro/mypaper ... MS2016.pdf

If I am reading it correctly, Table 1 gives the error in the real and imaginary part, but not the actual function value. If you post your corresponding function values for M(a,b,z) for Table 1 examples and run-times, when my Mathematica run is done, I can post the same from Mathematica and we can compute your errors. To put the pressure on, I'll wait on your numbers before posting mine.
Last edited by Alan on October 16th, 2016, 4:41 pm, edited 1 time in total.

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

### Re: Looking for hardware recommendations

OK.
I am now first doing sanity checks for benign values of z. I am taking a fixed N terms in the series but I might modify the stopping criterion. I get correct sanity results for

exp(z) and C++11
A&S (real)
Kummer's first formula
3-term recurrence relation

For kick-off I get stuff like

M(2,3, Complex(1/2, 1/2)) == (1.31289906741, 0.4746562395)
M(2,2, Complex(1/2, 1/2)) == (1.4468890364, 0.7904390832)
M(1,1, Complex(10, 10)) == std::exp(Complex(10, 10)) to 9 decimal places.

===
N.B. Above using double precision.
Last edited by Cuchulainn on October 16th, 2016, 8:20 pm, edited 5 times in total.
http://www.datasimfinancial.com
http://www.datasim.nl

Every Time We Teach a Child Something, We Keep Him from Inventing It Himself
Jean Piaget

Alan
Topic Author
Posts: 10015
Joined: December 19th, 2001, 4:01 am
Location: California
Contact:

### Re: Looking for hardware recommendations

Not sure if "setting up a Debian Linux VM"  looks complicated or not?

https://williamjturkel.net/2013/05/31/i ... -a-vm/amp/
It's a good idea among many good ones in this thread. So far we have, in what seems a rough order of difficulty/complexity for me personally (easiest to hardest):

-  Run my current Mathematica on newer/better hardware with more cores.
-  Port my algorithm to C/C++ (Windows/VisualStudio) with "roll your own" special functions (Daniel's simple series sum).
-  Port the algorithm to FORTRAN (Windows/VisualStudio), using Intel's compiler and perhaps parallelization tools.
- Port the algorithm to C (unix-type system, perhaps virtual) in order to use some special function libs for that.

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

### Re: Looking for hardware recommendations

This article, eq (4) shows CHF as an absolutely convergent series, so now that C++11 and Boost support the needed building blocks I thought it could be programmed on a lazy Sunday afternoon. I tested against A&S, C++ std::exp(z) and the famous $e^5$. It's about 6 lines of code. I just tried it to kick the tyres. Either it works or it doesn't. And it seems feasible.
Larger values if z ==> more terms in (4) and you might also need multiprecision at some stage.
https://dspace.mit.edu/bitstream/handle ... sequence=1

Have I produced the first mainstream/accessible results in C++ for CHF?

I'll post on $e^5$ after I review + take some more tests.
Thanks for experimenting. Here's a paper with a few tests you might want to try (Table 1, 3):
http://www.cs.upc.edu/~argimiro/mypaper ... MS2016.pdf

If I am reading it correctly, Table 1 gives the error in the real and imaginary part, but not the actual function value. If you post your corresponding function values for M(a,b,z) for Table 1 examples and run-times, when my Mathematica run is done, I can post the same from Mathematica and we can compute your errors. To put the pressure on, I'll wait on your numbers before posting mine.
Table 1 is not all that clear. Anyways, I take the first three examples using C++ 100-decimal digit data type (double precision gives overflow for largish z):
M(1,4,50i) -> (0.00241259etc,,, 0.059983..)
M(3,10,30 + 100i) -> (-13811.6.., -953.439..)
M(15,20,200i)-> (-2.82866E-6, -3.27839E-6)
==
M(4,1,50i) -> (-9044.44, -18975.1)
M(10,3,30+100i) -> (-1.8468E+20, 9.77876E+21)
M(4000, 4200,50i) -> (-0868055, -0.46886).

As a sanity check I took M(a,a,z) for various z and checked against std::exp(z). Accuracy is O(E-17).

For some of the other values I get overflow even with 100 digits; this demands more research to complement the series solution, such as asymptotics, term rewriting using Euler's formula and checking for overflow. I am also not 100% how Multiprecision handles maths functions.

I haven't yet paid attention to speed but it's of the order of [10,20] seconds in Release(!) mode. It is trivial to parallelise and almost linear speedup should be possible. For later maybe.

I'll wait on feedback regarding above. I only printed the first few digits.

// This exercise is a good test of new features in C++11 and Boost. A lot of stuff being used.
// page 509 of A&S is very impressive.
http://www.datasimfinancial.com
http://www.datasim.nl

Every Time We Teach a Child Something, We Keep Him from Inventing It Himself
Jean Piaget

Alan
Topic Author
Posts: 10015
Joined: December 19th, 2001, 4:01 am
Location: California
Contact:

### Re: Looking for hardware recommendations

Mtest[a_,b_,z_,trials_]:=
Timing[Do[M[a,b,z],{trials-1}];
M[a,b,z]]
Thanks again for your experiments. I just tested with the above Mathematica code. The result gives the timing (in secs) to do the computation 'trials' times, and then the M value.

Taking 10,000 trials:

Mtest[1.,4.,50 I,10000]
{1.185608,0.00241259 +0.0599983 I}

Mtest[3.,10.,30+100. I,10000]
{1.747211,-13811.6-953.439 I}

Mtest[15.,20.,200. I,10000]
{1.326008,-2.82866*10^-6-3.27839*10^-6 I}

Mtest[4.,1.,50. I,10000]
{0.046800,-9044.44-18975.1 I}

Mtest[10.,3.,30+100. I,10000]
{0.078001,-1.8468*10^20+9.77876*10^21 I}

Mtest[4000.,4200.,50. I,10000]
{0.249602,-0.868055-0.468886 I}

So, roughly, the results are [0.05 - 1.75] secs to do 10,000 evals. To make a port to C/C++ (or FORTRAN or whatever) worth the effort, I would need a factor of 10x improvement on these times or better, as well as a similar speedup with the Bessel $I_{\nu}(z)$ function. (Also, $M(a,b,z)$ needs to admit complex 'a' and 'z', and $I_{\nu}(z)$ needs complex z). I know you haven't paid attention to speed, but that's where I am on the ballpark requirements.

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

### Re: Looking for hardware recommendations

There is no silver bullet; I think the requirements need to be scoped and even then you could spend a lot of time on the different regimes. My series approach works for abs(z) < 50 (see Pearson et al section 3.8) and it is very fast. After that other approaches are needed as discussed here

- range of a,b,z
- real, complex a,b,z
- accuracy desired
- etc.

https://arxiv.org/pdf/1407.7786.pdf
http://www.datasimfinancial.com
http://www.datasim.nl

Every Time We Teach a Child Something, We Keep Him from Inventing It Himself
Jean Piaget