[saumavadey]

I have written two functions for Modified Bessel Functions of First and Second Kind of Order 0 as besselI0 and besselK0 respectively.

The function besselK0 calls the function besselI0 within itself and has been included as a header file.

The codes are as follows:

namespace Foam
{
scalar besselI0(scalar X)
{
const scalar A1 = 1.0;
const scalar A2 = 3.5156229;
const scalar A3 = 3.0899424;
const scalar A4 = 1.2067492;
const scalar A5 = 0.2659732;
const scalar A6 = 0.0360768;
const scalar A7 = 0.0045813;

const scalar B1 = 0.39894228;
const scalar B2 = 0.01328592;
const scalar B3 = 0.00225319;
const scalar B4 = -0.00157565;
const scalar B5 = 0.00916281;
const scalar B6 = -0.02057706;
const scalar B7 = 0.02635537;
const scalar B8 = -0.01647633;
const scalar B9 = 0.00392377;

scalar T;
scalar I;
if (mag(X) < 3.75)
{
T = sqr(X/3.75);
I = A1 + T*(A2 + T*(A3 + T*(A4 + T*(A5 + T*(A6 + T*A7)))));
}
else
{
T = 3.75/mag(X);
I = B1 + T*(B2 + T*(B3 + T*(B4 + T*(B5 + T*(B6 + T*(B7 + T*(B8 + T*B9)))))));
I = (I*exp(mag(X)))/sqrt(mag(X));
}
return I;
}
}

===============================================
#include "besselI0.H"

namespace Foam
{
scalar besselK0(scalar X)
{
const scalar A1 = -0.57721566;
const scalar A2 = 0.42278420;
const scalar A3 = 0.23069756;
const scalar A4 = 0.03488590;
const scalar A5 = 0.00262698;
const scalar A6 = 0.00010750;
const scalar A7 = 0.00000740;

const scalar B1 = 1.25331414;
const scalar B2 = -0.07832358;
const scalar B3 = 0.02189568;
const scalar B4 = -0.01062446;
const scalar B5 = 0.00587872;
const scalar B6 = -0.00251540;
const scalar B7 = 0.00053208;

scalar T;
scalar K;
if (X <= 2.0)
{
T = sqr(X/2.0);
K = -(log(X/2.0)*besselI0(X)) + A1 + T*(A2 + T*(A3 + T*(A4 + T*(A5 + T*(A6 + T*A7)))));
}
else
{
T = 2.0/X;
K = B1 + T*(B2 + T*(B3 + T*(B4 + T*(B5 + T*(B6 + T*B7)))));
K = K/(exp(X)*sqrt(X));
}
return K;
}
}
===============================================
While compiling I get the following error:

In file included from besselK0.H:11:0, besselI0.H: In function ‘Foam::scalar Foam::besselI0(Foam::scalar)’: besselI0.H:13:9: error: redefinition of ‘Foam::scalar Foam::besselI0(Foam::scalar)’ scalar besselI0(scalar X) ^~~~~~~~

In file included from besselK0.H:11:0, besselI0.H:13:9: note: ‘Foam::scalar Foam::besselI0(Foam::scalar)’ previously defined here scalar besselI0(scalar X)
===============================================

I am not being able to figure out where am I getting it wrong. I seek help in recognizing the error and how to get out of it.

Regards,

Saumava

[raumpolizei]

Hey, could you share the whole header file BesselI0.H? Is the header file included elsewhere? Is BesselIO.H including the preprocessor directive ->

#ifndef BESSELI0_H
#define BESSELI0_H
//your awesome function here
#endif

[saumavadey]

You got this right. Adding the header guards worked.

Thanks a lot.

[olesen]

Originally Posted by saumavadey

I have written two functions for Modified Bessel Functions of First and Second Kind of Order 0 as besselI0 and besselK0 respectively.

The function besselK0 calls the function besselI0 within itself and has been included as a header file.

Although you have things working, if the function body is actually defined in the header it should likely also be declared "inline".

To help the compiler out, the various A, B constants should be constexpr. Finally, reduce the use of additional (non-const) variables to help convince the compiler some more.

For example,

namespace Foam
{
	inline scalar besselI0(const scalar X)
	{
		constexpr scalar A1 = 1.0;
		constexpr scalar A2 = 3.5156229;
		...
		
		if (mag(X) < 3.75)
		{
			const scalar T = sqr(X/3.75);

			return A1 + T*(A2 + T*(A3 + T*(A4 + T*(A5 + T*(A6 + T*A7)))));
		}
		else
		{
			const scalar T = 3.75/mag(X);

			return (B1 + T*(B2 + T*(B3 + T*(B4 + T*(B5 + T*(B6 + T*(B7 + T*(B8 + T*B9)))))))) 
		           * exp(mag(X)/sqrt(mag(X);
		}
	}			
}