C++ code calling GSL (GNU Scientific Library) functions with liboctave

Octave is amazing! I started writing C++ code using liboctave as a library. However, now and then, I need to call some GSL functions from my C++ code. In octave it’s pretty straight forward, you just need to install GSL package from octave-forge.  But, hey, we can extend to our C++ code as well. This is because the GSL from octave-forge is just a wrapper written in C++. That’s the beauty of open source.

So if you want to know more, download GSL from octave-forge. Then unzip them
$ tar xzvf gsl-1.0.7.tar.gz
$ cd gsl-1.0.7
$ ./configure
$ make

Wait, no need to do make install, you just need to go to the folder “src”. Make sure you read the “README” file. I put a paragraph below

There are several different templates to use depending
on the inputs to the special function.  For example, if
the function takes two doubles and returns a double, use
the double_double_to_double.cc template.

To see the available templates, use:

ls *template

As of this writing there are:

double fn(double)
double fn(double,double)
double fn(int)
double fn(int,double)
double fn(int,int,double)
double fn(double,mode)

so Octave provided some templates depending on your GSL function. If you have type make, then the build process will generate the .cc file from the template. Since I am interested in using the GSL special function, I need to look at “gsl_sf.cc”.

Copy this file to your folder and modify it. Here I give an example for computing Complete Integral Fermi Dirac function of order -1/2.

NOTE: ARgh, wordpress does not preserve the spacing and the tab. You can click on the links below to get the file:
testFermi.cc
mygsl_sf.h
fermi_dirac_mhalf.cc
Makefile
The main input file for testing is:

// testFermi.cc
#include
#include
#include "mygsl_sf.h"

int main(int argc,char* argv[])
{
double dbInput=2.0;
ColumnVector cvInput(10,2.0);
octave_value_list result=fermi_dirac_mhalf(octave_value(dbInput));
double resultScalar=result(0).scalar_value();
std::cout<<"Test scalar on fermi_dirac_mhalf(2.0): " << resultScalar <<"\n";
//modify 2nd and last element
cvInput(1)=3.0;
cvInput(9)=5.0;
std::cout<<
"Test fermi_dirac_mhalf with vector input\n Input vector : "<<
cvInput <<"\n";
result=fermi_dirac_mhalf(octave_value(cvInput));
ColumnVector resultVect(result(0).array_value());
std::cout<<"results : " << resultVect <<"\n";
return 0;
}

Then I declare my C++ to call GSL function in “mygsl_sf.h”

//mygsl_sf.h
#include
#include
#include
#include

octave_value_list fermi_dirac_mhalf(octave_value_list args);

After that, I write the function “fermi_dirac_mhalf.cc”

//fermi_dirac_mhalf.cc
#include
#include
#include
#include
#include "mygsl_sf.h"

octave_value_list fermi_dirac_mhalf(octave_value_list args)
{
int i;

NDArray x = args(0).array_value();
NDArray y(x.dims());
int lx = x.length();
// printf("length: %d\n", lx);
// printf("nargout: %d\n", nargout);
for(i = 0; i < lx; i++) {
y.xelem(i) = gsl_sf_fermi_dirac_mhalf(x.xelem(i));
}
return octave_value(y);

}

Now, let’s create the Makefile

OCTAVE_INCLUDE = /home/kurniawano/Programs/Octave3.0.1/include/octave-3.0.1/
GSL_INCLUD = /home/kurniawano/local/include
LIBS = -lgsl -lgslcblas -lm

all: testFermi

clean:
rm *.o testFermi

testFermi: testFermi.o fermi_dirac_mhalf.o
mkoctfile --link-stand-alone -o testFermi fermi_dirac_mhalf.o testFermi.o -I$(GSL_INCLUDE) $(LIBS)

testFermi.o: testFermi.cpp
g++ -c -I$(OCTAVE_INCLUDE) -I$(OCTAVE_INCLUDE)octave -o testFermi.o testFermi.cpp

fermi_dirac_mhalf.o: fermi_dirac_mhalf.cc
g++ -c -I$(OCTAVE_INCLUDE) -I$(OCTAVE_INCLUDE)octave -o fermi_dirac_mhalf.o fermi_dirac_mhalf.cc

Then, you can type make to try this codes. Enjoy!

Advertisements

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out / Change )

Twitter picture

You are commenting using your Twitter account. Log Out / Change )

Facebook photo

You are commenting using your Facebook account. Log Out / Change )

Google+ photo

You are commenting using your Google+ account. Log Out / Change )

Connecting to %s

%d bloggers like this: