I'm trying to generate combinations of M integers from a set of M^2 consecutive integers, then exclude some combinations based on some (currently irrelevant) criteria, and trying to store them dynamically in memory so that I may be able to access them later.
The combinations are being generated using the 'gsl_combination' structure in the GNU Scientific Library (docs here). As you can see from the code snippet below, the combinations are being generated okay (and exclusions are as desired).
Then I need to append each combination structure object to memory for later retrieval. I'm trying to use glib pointer arrays to do this.
Code: gptrwierd.c
#include <glib.h>
#include <gsl/gsl_combination.h>
//Hardcoded the size of the rigol lattice
//MUST BE ODD !!!
#define M 3
//Returns GSL_SUCCESS if the combination contains a forbidden site, otherwise returns GSL_FAILURE
int
is_forbidden (const gsl_combination * test)
{
int i, forbidden = GSL_FAILURE; //Not forbidden by default
size_t k = test->k;
size_t *data = test->data;
int halfway = ((M * M - 1) / 2); //The halfway point in the 2D grid
...
... //Some code for updating 'forbidden'
...
return forbidden;
}
#undef __FUNCT__
#define __FUNCT__ "main"
int
main (int argc, char **argv)
{
long count, istride;
//Generate the basis states
gsl_combination *c = gsl_combination_calloc (M * M, M);
//An array that will contain selected basis vectors.
GPtrArray *basis = g_ptr_array_new ();
count = 0;
do
{
//c is the ith combination of M integers from M^2 integers
//Get all the non forbidden elements
if (is_forbidden (c) == GSL_FAILURE)
{ //If the site is NOT forbidden, append c to basis
g_ptr_array_add (basis, c); //Appends c to basis
{
printf ("count %ld {", count);
gsl_combination_fprintf (stdout,
g_ptr_array_index (basis, count), " %u");
printf (" }\n");
}
count++;
}
}
while (gsl_combination_next (c) == GSL_SUCCESS);
printf("\n\n");
//Now, access each basis element from GArray
for (istride = 0; istride < basis->len; istride++)
{
printf ("istride %ld {", istride);
gsl_combination_fprintf (stdout, g_ptr_array_index (basis, istride),
" %u");
printf (" }\n");
}
gsl_combination_free (c);
g_ptr_array_free (basis, TRUE);
return 0;
}
Makefile:
GSL_FLAGS = gsl-config --cflags
GLIB_LOC = pkg-config --cflags --libs glib-2.0
CFLAGS = -Wall -O3 `$(GSL_FLAGS)` `$(GLIB_LOC)`
GSL_LIBS = gsl-config --libs-without-cblas
GLIB_LIBS = pkg-config --libs glib-2.0
gptrwierd:
gcc $(CFLAGS) -c gptrwierd.c
gcc -o $@ gptrwierd.o -lgsl -lgslcblas -lm `$(GLIB_LIBS)`
allclean: clean
${RM} gptrwierd *.dat
Compilation and running is done thus (M is hardcoded to 3):
$make gptrwierd
gcc -Wall -O3 `gsl-config --cflags` `pkg-config --cflags --libs glib-2.0` -c gptrwierd.c
gcc -o gptrwierd gptrwierd.o -lgsl -lgslcblas -lm `pkg-config --libs glib-2.0`
$./gptrwierd > out.dat
contents of output file are pasted below:
count 0 { 1 2 3 }
count 1 { 1 2 4 }
count 2 { 1 2 6 }
count 3 { 1 2 7 }
count 4 { 1 2 8 }
count 5 { 1 3 4 }
count 6 { 1 3 6 }
count 7 { 1 3 7 }
count 8 { 1 3 8 }
count 9 { 1 4 6 }
count 10 { 1 4 7 }
count 11 { 1 4 8 }
count 12 { 1 6 7 }
count 13 { 1 6 8 }
...
...
count 28 { 3 6 7 }
count 29 { 3 6 8 }
count 30 { 3 7 8 }
count 31 { 4 6 7 }
count 32 { 4 6 8 }
count 33 { 4 7 8 }
count 34 { 6 7 8 }
istride 0 { 6 7 8 }
istride 1 { 6 7 8 }
istride 2 { 6 7 8 }
istride 3 { 6 7 8 }
istride 4 { 6 7 8 }
istride 5 { 6 7 8 }
istride 6 { 6 7 8 }
istride 7 { 6 7 8 }
istride 8 { 6 7 8 }
istride 9 { 6 7 8 }
istride 10 { 6 7 8 }
istride 11 { 6 7 8 }
istride 12 { 6 7 8 }
...
...
istride 30 { 6 7 8 }
istride 31 { 6 7 8 }
istride 32 { 6 7 8 }
istride 33 { 6 7 8 }
istride 34 { 6 7 8 }
As you can see, it prints the combinations correctly in the do...while loop where the combinations are iterated, but only prints the lexicographically final combination when the gptrarray elements are iterated in the subsequent loop.
What am I doing wrong?
Update: One possible workaround seems to be to do away with glib altogether and manually copy the combination data to a dynamically allocated array:
gsl_combination *c = gsl_combination_calloc (M * M, M);
long **basis;
long dim = 0;
//Evaluate the dimensionality of the Hilbert space (# of allowed combinations)
do
{
//c is the ith combination of M integers from M^2 integers
//Get all the non forbidden elements
if (is_forbidden (c) == GSL_FAILURE)
{ //If the site is NOT forbidden, append dim
dim++;
}
}
while (gsl_combination_next (c) == GSL_SUCCESS);
//Now, generate the actual basis
//basis is an array of arrays that will contain the
//selected basis vectors. Each basis vector is an array of M integers
basis = (long **) malloc (dim * sizeof (long *));
for (count = 0; count < dim; count++)
basis[count] = (long *) malloc (M * sizeof (long));
//Reset the combination to the one that is lexicographically first
gsl_combination_init_first (c);
count = 0;
//Now, append all allowed combinations to the basis
do
{
//c is the ith combination of M integers from M^2 integers
//Get all the non forbidden elements
if (is_forbidden (c) == GSL_FAILURE)
{ //If the site is NOT forbidden, append data in c to basis
for (istride = 0; istride < M; istride++)
basis[count][istride] = c->data[istride];
count++;
}
}
while (gsl_combination_next (c) == GSL_SUCCESS);
//Now, access each basis element from GArray
//basis[i] is the ith combination. jth combination element is basis[i][j]
printf (
"Printing all allowed combinations of sites in basis\n"
"---------------------------------------------------\n");
for (istride = 0; istride < dim; istride++)
{
printf ("stride # %3ld:\t{", istride);
for (count = 0; count < M; count++)
printf (" %ld", basis[istride][count]);
printf (" }\n");
}
/*Do whatever*/
...
...
free(basis);