Sparse Poisson problem in eigen

OwlgenBack to scientific computing. Lately, I have been using the Eigen libraries for linear algebra. They were recommended by the CGAL project, and they indeed share a common philosophy. Thanks to the rather high C++ level they can accomplish this sort of magic:

  int n = 100;
  VectorXd x(n), b(n);
  SpMat A(n,n);



  // solve Ax = b
  Eigen::BiCGSTAB<SpMat> solver;
  //ConjugateGradient<SpMat> solver;


  x = solver.solveWithGuess(b,x0);

Notice that A is a sparse matrix! I am next describing how to use this in order to solve the 1D Poisson equation.

Continue reading

CGAL from scratch

So, you need computational geometry. It could be you need some nice Voronoi division of some domain. Or you need to mesh a surface. Whatever your needs, you may know there is this CGAL project that has these things implemented in a very professional manner.

Here is a list of steps in order to get these routines working for you.

Continue reading

Parameters and include files

Actually, scrap what I said yesterday. The common situation still is: you have a set of parameters that are basically global to the program. But, of course, it is not good to go totally global. OK, I think the best idea is to define a class in a header file. Then, declare one of these classes in the header file as external. Then declare (for real) just once somewhere. Ideally a singleton class should be used, but this is a bit over my head, just remember to declare one of these classes (and just once). (Unless you do want several sets of parameters, who knows…)

The members of the funtion are all public, and contain the parameters. Plain static members don’t make much sense (since there will only be one of these classes), but constant data should be declared static const, and declared somewhere (globally).

So, for example:

file “defs.h”


using std::endl;
using std::cout;

class dta {

double f;

static const int i;


extern dta data;

Then, define stuff somewhere:

file “defs.cpp”


dta data;

const int dta::i=10;

void defs(void) {  data.f=0; }


file “main.cpp”


int main(void) {

cout << data.i << endl;
cout << data.f << endl;
cout << “Done” << endl;
return 1;

Just compile with g++ main.cpp defs.cpp.

Namespaces and include files in c++

This is a common situation (at least for me): you have a set of parameters that are basically global to the program. But, of course, it is not good to go totally global, so one may try a namespace, thus:

namespace fun {
double f(double);
int i=5;

This way, one can keep simple names (“f” and “i”), but still call them fun::f, and fun::i. This is ok as long as one has one source cpp file (or one include file, included only once). What happens when one splits the code into several files, as one should eventually? f() is fine, but i get defined twice… error!

The way to go: simple, just use extern!. Examples follow.

Continue reading

The CGAL project

CGAL logo

I have been following and using this project for some years now, on and off. In my opinion, a very remarkable cooperative effort of several research centers, coordinated at INRIA Sophia Antipolis. They provide very robust implementations of algorithms on computational geometry: Delaunay triangulations, Voronoi diagrams ski sums, meshes, alpha shapes, convex hulls, all that. Warning: very abstract C++ is needed (I actually learned C++ because of these libraries).

Continue reading