Add MEX'ed edition of 'partitionUI'. The MEX function only
implements the first edition of the 'partitionUI' algorithm (i.e., before afg added sub-grid support in r2349), but it is nevertheless a useful experiment and facilitates testing a C/C++ implementation of the MsMFE method.
This commit is contained in:
Bård Skaflestad
parent
0b8b26e6c2
commit
2a30852257
155
mex_partition_ui.c
Normal file
155
mex_partition_ui.c
Normal file
@ -0,0 +1,155 @@
|
|||||||
|
#include <stddef.h>
|
||||||
|
#include <string.h>
|
||||||
|
|
||||||
|
#include <mex.h>
|
||||||
|
|
||||||
|
#include "partition.h"
|
||||||
|
|
||||||
|
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
static int
|
||||||
|
args_ok(int nlhs, int nrhs, const mxArray *prhs[])
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
{
|
||||||
|
int i, ok;
|
||||||
|
|
||||||
|
ok = nlhs == 1;
|
||||||
|
|
||||||
|
for (i = 0; ok && (i < nrhs); i++) {
|
||||||
|
ok = mxIsDouble(prhs[i]) || mxIsInt32(prhs[i]);
|
||||||
|
}
|
||||||
|
|
||||||
|
return ok && (nrhs == 3);
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
static void
|
||||||
|
release_rhs(int *idx, int *fine_d, int *coarse_d)
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
{
|
||||||
|
if (coarse_d != NULL) { mxFree(coarse_d); }
|
||||||
|
if (fine_d != NULL) { mxFree(fine_d); }
|
||||||
|
if (idx != NULL) { mxFree(idx); }
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
static void
|
||||||
|
copy_intvec(const mxArray *M_i, int *i)
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
{
|
||||||
|
size_t e, ne;
|
||||||
|
|
||||||
|
int *pi;
|
||||||
|
double *pd;
|
||||||
|
|
||||||
|
ne = mxGetNumberOfElements(M_i);
|
||||||
|
|
||||||
|
if (mxIsDouble(M_i)) {
|
||||||
|
pd = mxGetPr(M_i);
|
||||||
|
|
||||||
|
for (e = 0; e < ne; e++) { i[e] = pd[e]; }
|
||||||
|
} else {
|
||||||
|
pi = mxGetData(M_i);
|
||||||
|
memcpy(i, pi, ne * sizeof *i);
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
static void
|
||||||
|
extract_rhs(const mxArray *prhs[],
|
||||||
|
int *nidx, int *ndims,
|
||||||
|
int **idx, int **fine_d, int **coarse_d)
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
{
|
||||||
|
int ni, nd, p, *i, *f, *c;
|
||||||
|
|
||||||
|
ni = mxGetNumberOfElements(prhs[0]);
|
||||||
|
nd = mxGetNumberOfElements(prhs[1]);
|
||||||
|
|
||||||
|
if (mxGetNumberOfElements(prhs[2]) != nd) {
|
||||||
|
mexErrMsgTxt("'coarseDim' must have same number of "
|
||||||
|
"elements as 'fineDim'.");
|
||||||
|
}
|
||||||
|
|
||||||
|
i = mxMalloc(ni * sizeof *i);
|
||||||
|
f = mxMalloc(nd * sizeof *f);
|
||||||
|
c = mxMalloc(nd * sizeof *c);
|
||||||
|
|
||||||
|
if ((i == NULL) || (f == NULL) || (c == NULL)) {
|
||||||
|
release_rhs(i, f, c);
|
||||||
|
|
||||||
|
*nidx = -1; *ndims = -1;
|
||||||
|
*idx = NULL; *fine_d = NULL; *coarse_d = NULL;
|
||||||
|
} else { /* Adjust for 1-based indexing */
|
||||||
|
copy_intvec(prhs[0], i); for (p = 0; p < ni; p++) { i[p]--; }
|
||||||
|
copy_intvec(prhs[1], f);
|
||||||
|
copy_intvec(prhs[2], c);
|
||||||
|
|
||||||
|
*nidx = ni; *ndims = nd;
|
||||||
|
*idx = i; *fine_d = f; *coarse_d = c;
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
static void
|
||||||
|
assign_lhs(const int *p, mxArray *lhs)
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
{
|
||||||
|
size_t e, ne;
|
||||||
|
|
||||||
|
int *pi;
|
||||||
|
double *pd;
|
||||||
|
|
||||||
|
ne = mxGetNumberOfElements(lhs);
|
||||||
|
|
||||||
|
if (mxIsDouble(lhs)) {
|
||||||
|
pd = mxGetPr(lhs);
|
||||||
|
|
||||||
|
for (e = 0; e < ne; e++) { pd[e] = p[e] + 1; }
|
||||||
|
} else {
|
||||||
|
pi = mxGetData(lhs);
|
||||||
|
|
||||||
|
for (e = 0; e < ne; e++) { pi[e] = p[e] + 1; }
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
/*
|
||||||
|
* p = mex_partition_ui(i, fineDim, coarseDim)
|
||||||
|
*/
|
||||||
|
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
void
|
||||||
|
mexFunction(int nlhs, mxArray *plhs[],
|
||||||
|
int nrhs, const mxArray *prhs[])
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
{
|
||||||
|
int nidx, ndims, *idx, *p, *fine_d, *coarse_d;
|
||||||
|
char errmsg[1023 + 1];
|
||||||
|
|
||||||
|
if (args_ok(nlhs, nrhs, prhs)) {
|
||||||
|
extract_rhs(prhs, &nidx, &ndims, &idx, &fine_d, &coarse_d);
|
||||||
|
|
||||||
|
if (nidx > 0) {
|
||||||
|
plhs[0] = mxDuplicateArray(prhs[0]);
|
||||||
|
p = mxMalloc(nidx * sizeof *p);
|
||||||
|
|
||||||
|
partition_unif_idx(ndims, nidx, fine_d, coarse_d, idx, p);
|
||||||
|
|
||||||
|
assign_lhs(p, plhs[0]);
|
||||||
|
|
||||||
|
mxFree(p); release_rhs(idx, fine_d, coarse_d);
|
||||||
|
}
|
||||||
|
} else {
|
||||||
|
sprintf(errmsg,
|
||||||
|
"Calling sequence is\n"
|
||||||
|
"\tp = %s(ix, fineDim, coarseDim)\n",
|
||||||
|
mexFunctionName());
|
||||||
|
|
||||||
|
mexErrMsgTxt(errmsg);
|
||||||
|
}
|
||||||
|
}
|
||||||
42
mex_partition_ui.m
Normal file
42
mex_partition_ui.m
Normal file
@ -0,0 +1,42 @@
|
|||||||
|
function varargout = mex_partition_ui(varargin)
|
||||||
|
%Partition grid uniformly in logical space using compiled C code.
|
||||||
|
%
|
||||||
|
% SYNOPSIS:
|
||||||
|
% p = mex_partition_ui(ix, fineDim, coarseDim)
|
||||||
|
%
|
||||||
|
% PARAMETERS:
|
||||||
|
% ix - Indices, e.g., G.cells.indexMap.
|
||||||
|
%
|
||||||
|
% fineDim - Cartesian dimensions of underlying index space (e.g.
|
||||||
|
% G.cartDims)
|
||||||
|
%
|
||||||
|
% coarseDim - Cartesian dimensions of requested index space. Typically,
|
||||||
|
% all(coarseDim <= fineDim).
|
||||||
|
%
|
||||||
|
% RETURNS:
|
||||||
|
% p - Partition vector of SIZE(ix).
|
||||||
|
%
|
||||||
|
% SEE ALSO:
|
||||||
|
% partitionUI.
|
||||||
|
|
||||||
|
%{
|
||||||
|
#COPYRIGHT#
|
||||||
|
%}
|
||||||
|
|
||||||
|
% $Date$
|
||||||
|
% $Revision$
|
||||||
|
|
||||||
|
buildmex CFLAGS="\$CFLAGS -Wall -Wextra -ansi -pedantic ...
|
||||||
|
-Wformat-nonliteral -Wcast-align -Wpointer-arith ...
|
||||||
|
-Wbad-function-cast -Wmissing-prototypes -Wstrict-prototypes ...
|
||||||
|
-Wmissing-declarations -Winline -Wundef -Wnested-externs ...
|
||||||
|
-Wcast-qual -Wshadow -Wconversion -Wwrite-strings ...
|
||||||
|
-Wno-conversion -Wchar-subscripts -Wredundant-decls" ...
|
||||||
|
...
|
||||||
|
-g -largeArrayDims -DCOMPILING_FOR_MATLAB=1 ...
|
||||||
|
...
|
||||||
|
mex_partition_ui.c partition.c
|
||||||
|
|
||||||
|
% Call MEX'ed edition.
|
||||||
|
[varargout{1:nargout}] = mex_partition_ui(varargin{:});
|
||||||
|
end
|
||||||
131
partition.c
Normal file
131
partition.c
Normal file
@ -0,0 +1,131 @@
|
|||||||
|
#include <assert.h>
|
||||||
|
#include <stddef.h>
|
||||||
|
#include <stdlib.h>
|
||||||
|
|
||||||
|
#include "partition.h"
|
||||||
|
|
||||||
|
|
||||||
|
#define MAX(a,b) (((a) > (b)) ? (a) : (b))
|
||||||
|
|
||||||
|
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
static void
|
||||||
|
partition_coord_idx(int ndims, int idx, const int *size, int *cidx)
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
{
|
||||||
|
int i;
|
||||||
|
|
||||||
|
for (i = 0; i < ndims; i++) {
|
||||||
|
cidx[i] = idx % size[i];
|
||||||
|
idx /= size[i];
|
||||||
|
}
|
||||||
|
|
||||||
|
assert (idx == 0);
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
static int
|
||||||
|
partition_lin_idx(int ndims, const int *size, const int *cidx)
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
{
|
||||||
|
int i, idx;
|
||||||
|
|
||||||
|
idx = cidx[ndims - 1];
|
||||||
|
for (i = ndims - 2; i >= 0; i--) {
|
||||||
|
idx = cidx[i] + size[i]*idx;
|
||||||
|
}
|
||||||
|
|
||||||
|
return idx;
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
/* Load-balanced linear distribution.
|
||||||
|
*
|
||||||
|
* See Eric F. Van de Velde, Concurrent Scientific Computing,
|
||||||
|
* 1994, Springer Verlag, p. 54 (Sect. 2.3) for details. */
|
||||||
|
static void
|
||||||
|
partition_loadbal_lin_dist(int ndims, const int *size, const int *nbins,
|
||||||
|
int *idx)
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
{
|
||||||
|
int i, L, R, b1, b2;
|
||||||
|
|
||||||
|
for (i = 0; i < ndims; i++) {
|
||||||
|
L = size[i] / nbins[i]; /* # entities per bin */
|
||||||
|
R = size[i] % nbins[i]; /* # bins containing one extra entity */
|
||||||
|
|
||||||
|
b1 = idx[i] / (L + 1);
|
||||||
|
b2 = (idx[i] - R) / L ;
|
||||||
|
|
||||||
|
idx[i] = MAX(b1, b2);
|
||||||
|
}
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
int
|
||||||
|
partition_unif_idx(int ndims, int nc,
|
||||||
|
const int *fine_d, const int *coarse_d, const int *idx,
|
||||||
|
int *p)
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
{
|
||||||
|
int c, ret, *ix;
|
||||||
|
|
||||||
|
ix = malloc(ndims * sizeof *ix);
|
||||||
|
|
||||||
|
if (ix != NULL) {
|
||||||
|
for (c = 0; c < nc; c++) {
|
||||||
|
partition_coord_idx(ndims, idx[c], fine_d, ix);
|
||||||
|
|
||||||
|
partition_loadbal_lin_dist(ndims, fine_d, coarse_d, ix);
|
||||||
|
|
||||||
|
p[c] = partition_lin_idx(ndims, coarse_d, ix);
|
||||||
|
}
|
||||||
|
|
||||||
|
ret = nc;
|
||||||
|
} else {
|
||||||
|
ret = -1;
|
||||||
|
}
|
||||||
|
|
||||||
|
free(ix);
|
||||||
|
|
||||||
|
return ret;
|
||||||
|
}
|
||||||
|
|
||||||
|
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
int
|
||||||
|
partition_compress(int n, int *p)
|
||||||
|
/* ---------------------------------------------------------------------- */
|
||||||
|
{
|
||||||
|
int ret, i, max, *compr;
|
||||||
|
|
||||||
|
max = -1;
|
||||||
|
for (i = 0; i < n; i++) {
|
||||||
|
assert (0 <= p[i]); /* Only non-neg partitions (for now?). */
|
||||||
|
max = MAX(max, p[i]);
|
||||||
|
}
|
||||||
|
|
||||||
|
compr = calloc(max + 1, sizeof *compr);
|
||||||
|
|
||||||
|
if (compr != NULL) {
|
||||||
|
for (i = 0; i < n; i++) { compr[p[i]]++; }
|
||||||
|
|
||||||
|
compr[0] = -1 + (compr[0] > 0);
|
||||||
|
for (i = 1; i <= max; i++) {
|
||||||
|
compr[i] = compr[i - 1] + (compr[i] > 0);
|
||||||
|
}
|
||||||
|
|
||||||
|
for (i = 0; i < n; i++) { p[i] = compr[p[i]]; }
|
||||||
|
|
||||||
|
ret = compr[max];
|
||||||
|
} else {
|
||||||
|
ret = -1;
|
||||||
|
}
|
||||||
|
|
||||||
|
free(compr);
|
||||||
|
|
||||||
|
return ret;
|
||||||
|
}
|
||||||
14
partition.h
Normal file
14
partition.h
Normal file
@ -0,0 +1,14 @@
|
|||||||
|
#ifndef PARTITION_H_INCLUDED
|
||||||
|
#define PARTITION_H_INCLUDED
|
||||||
|
|
||||||
|
int
|
||||||
|
partition_unif_idx(int ndims, int nc,
|
||||||
|
const int *fine_d,
|
||||||
|
const int *coarse_d,
|
||||||
|
const int *idx,
|
||||||
|
int *p);
|
||||||
|
|
||||||
|
int
|
||||||
|
partition_compress(int n, int *p);
|
||||||
|
|
||||||
|
#endif /* PARTITION_H_INLCUDED */
|
||||||
Loading…
Reference in New Issue
Block a user