I need to use a Fermi Dirac Integral function of order -3/2 which is not available in Octave. Octave only has order -1/2, 1/2, and 3/2 from its GSL package. The source code is available in Fortran from

A. Trellakis, A. J. Galick, U. Ravaioli, Rational Chebyshev approximation for the Fermi-Dirac integral F_{−3/2}(x). Solid-State Electronics, Volume 41, Issue 5, May 1997, Pages 771-773.

To create an Octave loadable function. Follow the instruction in this Octave Wiki. Basically, we need to create a C++ function, say *fm3half.cc*

#include

#include

#include "f77-fcn.h"

extern "C"

{

int F77_FUNC (fm3hf77, FM3HF77) (double& vret, const double& x);

}

DEFUN_DLD (fm3half, args, ,

"- Loadable Function: [FM3HALF] = fm3half ()\n\

\n\

Returns the fermi dirac integral of order -3/2.")

{

octave_value_list retval;

const double x=args(0).double_value(); //=args(0).double_value();

double vret;

F77_XFCN (fm3hf77, FM3HF77, (vret,x));

if (f77_exception_encountered)

{

error ("unrecoverable error in fm3half");

return retval;

}

retval(0) = octave_value (vret);

return retval;

}

which load a Fortran function called *fm3hf77.f*. You can take a look at the fortran function from this sourcecodesfm3half PDF file.

Note that the Fortran function is written as a subroutine and the arguments are passed by reference. The C++ function takes the first argument *x* and convert it from octave_value to double type. Then it passes to the Fortran function. The Fortran function saves the result in the first parameter *vret*, which is then converted to octave_value type.

Filed under: Uncategorized | Tagged: fortran, octave | 1 Comment »