Mercurial > hg > camir-aes2014

annotate toolboxes/FullBNT-1.0.7/graph/triangulate.c @ 0:e9a9cd732c1e tip

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
first hg version after svn
author wolffd
date Tue, 10 Feb 2015 15:05:51 +0000
parents
children
rev   line source
wolffd@0 1 /* triangulate.c written by Ilya Shpitser */
wolffd@0 2
wolffd@0 3 #include <stdlib.h>
wolffd@0 4
wolffd@0 5 #ifdef UNIX
wolffd@0 6 #include "matlab.h"
wolffd@0 7 #endif
wolffd@0 8
wolffd@0 9 #include "matrix.h"
wolffd@0 10 #include "mex.h"
wolffd@0 11
wolffd@0 12 #include "elim.h"
wolffd@0 13 #include "map.h"
wolffd@0 14 #include "misc.h"
wolffd@0 15
wolffd@0 16 void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]){
wolffd@0 17
wolffd@0 18 int dims [2];
wolffd@0 19 int i, j, k, m, n;
wolffd@0 20 long index;
wolffd@0 21 double * G_pr;
wolffd@0 22 double * stage_pr;
wolffd@0 23 double * answer_G_pr, * fill_ins_pr;
wolffd@0 24 double * matlab_clique_pr;
wolffd@0 25 mxArray * matlab_clique;
wolffd@0 26 Elimination e;
wolffd@0 27 float ** adj_mat;
wolffd@0 28 int ** order = (int **) NULL;
wolffd@0 29 Iterator iter, iter2;
wolffd@0 30 word w, w2;
wolffd@0 31 int ** fill_ins;
wolffd@0 32 Map cliques;
wolffd@0 33 Map clique;
wolffd@0 34 mxArray * fill_ins_mat;
wolffd@0 35 int * nodes;
wolffd@0 36 mxArray * full;
wolffd@0 37
wolffd@0 38 // (original) full = mlfFull((mxArray *) prhs[0]);
wolffd@0 39 full = (mxArray *) mlfFull((mxArray *) prhs[0]); // added typecasting
wolffd@0 40 /* Obtain graph matrix information. */
wolffd@0 41 m = mxGetM(full);
wolffd@0 42 n = mxGetN(full);
wolffd@0 43 G_pr = mxGetPr(full);
wolffd@0 44
wolffd@0 45 if(n < 1 || m < 1){
wolffd@0 46 return;
wolffd@0 47 }
wolffd@0 48
wolffd@0 49 /* Allocate and populate the log weight adjacency matrix corresponding
wolffd@0 50 to the input graph. */
wolffd@0 51 adj_mat = (float **) malloc(sizeof(float *) * m);
wolffd@0 52 adj_mat[0] = (float *) malloc(sizeof(float) * m * n);
wolffd@0 53 for(i = 1; i < m; i++){
wolffd@0 54 adj_mat[i] = adj_mat[i - 1] + n;
wolffd@0 55 }
wolffd@0 56 /* We no longer have log weight info, but we have a (total) ordering on
wolffd@0 57 the nodes already, so we do not need this information. */
wolffd@0 58 for(i = 0; i < m; i++){
wolffd@0 59 for(j = 0; j < n; j++){
wolffd@0 60 index = j * m + i;
wolffd@0 61 if(G_pr[index] > 0){
wolffd@0 62 adj_mat[i][j] = 1;
wolffd@0 63 } else {
wolffd@0 64 adj_mat[i][j] = 0;
wolffd@0 65 }
wolffd@0 66 }
wolffd@0 67 }
wolffd@0 68
wolffd@0 69 /* Convert the total elimination ordering into a partial order argument
wolffd@0 70 for the elimination routine. The elimination routine's purpose in this
wolffd@0 71 mode of operation is to return cliques and fill-in edges. */
wolffd@0 72 if(nrhs > 1){
wolffd@0 73 order = (int **) malloc(sizeof(int *) * m);
wolffd@0 74 order[0] = (int *) malloc(sizeof(int) * m * n);
wolffd@0 75 for(i = 1; i < m; i++){
wolffd@0 76 order[i] = order[i - 1] + n;
wolffd@0 77 }
wolffd@0 78 for(i = 0; i < m; i++){
wolffd@0 79 for(j = 0; j < n; j++){
wolffd@0 80 order[i][j] = 0;
wolffd@0 81 }
wolffd@0 82 }
wolffd@0 83 stage_pr = mxGetPr(prhs[1]);
wolffd@0 84 for(i = 0; i < mxGetN(prhs[1]) - 1; i++){
wolffd@0 85 order[(int) stage_pr[i] - 1][(int) stage_pr[i + 1] - 1] = 1;
wolffd@0 86 }
wolffd@0 87 }
wolffd@0 88
wolffd@0 89 /* Find the elimination ordering. */
wolffd@0 90 e = find_elim(n, adj_mat, order, -1);
wolffd@0 91
wolffd@0 92 /* Allocate memory for the answer, and set the answer. */
wolffd@0 93 plhs[0] = mxCreateDoubleMatrix(m, n, mxREAL);
wolffd@0 94 answer_G_pr = mxGetPr(plhs[0]);
wolffd@0 95 cliques = get_cliques(e);
wolffd@0 96 /*
wolffd@0 97 dims[0] = 1;
wolffd@0 98 dims[1] = get_size_Map(cliques);
wolffd@0 99 plhs[1] = mxCreateCellArray(2, (const int *) dims);*/
wolffd@0 100 plhs[1] = mxCreateCellMatrix(get_size_Map(cliques), 1);
wolffd@0 101 fill_ins = get_fill_ins(e);
wolffd@0 102 fill_ins_mat = mxCreateDoubleMatrix(m, n, mxREAL);
wolffd@0 103 fill_ins_pr = mxGetPr(fill_ins_mat);
wolffd@0 104
wolffd@0 105 for(i = 0; i < n; i++){
wolffd@0 106 for(j = 0; j < m; j++){
wolffd@0 107 index = j * m + i;
wolffd@0 108 answer_G_pr[index] = G_pr[index];
wolffd@0 109 if(fill_ins[i][j] > 0){
wolffd@0 110 answer_G_pr[index] = 1;
wolffd@0 111 fill_ins_pr[index] = 1;
wolffd@0 112 }
wolffd@0 113 }
wolffd@0 114 }
wolffd@0 115 mxDestroyArray(full);
wolffd@0 116 // (original) plhs[2] = mlfSparse(fill_ins_mat, NULL, NULL, NULL, NULL, NULL);
wolffd@0 117 plhs[2] = (mxArray *) mlfSparse(fill_ins_mat, NULL, NULL, NULL, NULL, NULL); // added typecasting
wolffd@0 118 mxDestroyArray(fill_ins_mat);
wolffd@0 119 nodes = (int *) malloc(sizeof(int) * n);
wolffd@0 120 k = 0;
wolffd@0 121 iter = get_Iterator(cliques);
wolffd@0 122 while(!is_empty(iter)){
wolffd@0 123 w = next_key(iter);
wolffd@0 124 clique = (Map) w.v;
wolffd@0 125 matlab_clique = mxCreateDoubleMatrix(1, get_size_Map(clique), mxREAL);
wolffd@0 126 matlab_clique_pr = mxGetPr(matlab_clique);
wolffd@0 127 for(i = 0; i < n; i++){
wolffd@0 128 nodes[i] = 0;
wolffd@0 129 }
wolffd@0 130 iter2 = get_Iterator(clique);
wolffd@0 131 while(!is_empty(iter2)){
wolffd@0 132 w2 = next_key(iter2);
wolffd@0 133 nodes[w2.i] = w2.i + 1;
wolffd@0 134 }
wolffd@0 135 j = 0;
wolffd@0 136 for(i = 0; i < n; i++){
wolffd@0 137 if(nodes[i] > 0){
wolffd@0 138 matlab_clique_pr[j++] = nodes[i];
wolffd@0 139 }
wolffd@0 140 }
wolffd@0 141 mxSetCell(plhs[1], k++, matlab_clique);
wolffd@0 142 }
wolffd@0 143 free(nodes);
wolffd@0 144
wolffd@0 145 /* Finally, free the allocated memory. */
wolffd@0 146 destroy_Elimination(e);
wolffd@0 147 if(adj_mat){
wolffd@0 148 if(adj_mat[0]){
wolffd@0 149 free(adj_mat[0]);
wolffd@0 150 }
wolffd@0 151 free(adj_mat);
wolffd@0 152 }
wolffd@0 153 if(order){
wolffd@0 154 if(order[0]){
wolffd@0 155 free(order[0]);
wolffd@0 156 }
wolffd@0 157 free(order);
wolffd@0 158 }
wolffd@0 159 }