From 691c3a3b7b46f1d34743fbe4517b08dcc70ebe17 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?B=C3=A5rd=20Skaflestad?= Date: Sun, 9 Jan 2011 15:15:43 +0000 Subject: [PATCH] By jrn's blessing, enable the alternative implementation of 'processgrid'. Also remove the existing implementation. MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit While here, assert copyright for 2011. Signed-off-by: Bård Skaflestad --- processgrid.c | 239 +------------------------------------------------- 1 file changed, 2 insertions(+), 237 deletions(-) diff --git a/processgrid.c b/processgrid.c index c03233ed..276c3fa5 100644 --- a/processgrid.c +++ b/processgrid.c @@ -13,8 +13,8 @@ //===========================================================================*/ /* -Copyright 2009, 2010 SINTEF ICT, Applied Mathematics. -Copyright 2009, 2010 Statoil ASA. +Copyright 2009, 2010, 2011 SINTEF ICT, Applied Mathematics. +Copyright 2009, 2010, 2011 Statoil ASA. This file is part of The Open Reservoir Simulator Project (OpenRS). @@ -46,7 +46,6 @@ along with OpenRS. If not, see . #include "mxgrdecl.h" -#if 0 /* ---------------------------------------------------------------------- */ static mxArray * allocate_nodes(size_t nnodes) @@ -454,237 +453,3 @@ mexFunction(int nlhs, mxArray *plhs[], mexErrMsgTxt(errmsg); } } - -#else - -void fill_grid(mxArray **out, struct processed_grid *grid) -{ - const char *names[] = {"nodes", "faces", "cells", "cartDims", "type"}; - mxArray *G = mxCreateStructMatrix(1,1,sizeof names / sizeof names[0],names); - - int i,j; - double *ptr; - - - /* nodes */ - const char *n2[] = {"num", "coords"}; - mxArray *nodes = mxCreateStructMatrix(1,1,sizeof n2 / sizeof n2[0],n2); - mxSetField(nodes, 0, "num", mxCreateDoubleScalar(grid->number_of_nodes)); - - mxArray *nodecoords = mxCreateDoubleMatrix(grid->number_of_nodes, 3, mxREAL); - ptr = mxGetPr(nodecoords); - - for (j=0;j<3;++j) - { - for(i=0; inumber_of_nodes; ++i) - { - ptr[i+grid->number_of_nodes*j] = grid->node_coordinates[3*i+j]; - } - } - mxSetField(nodes, 0, "coords", nodecoords); - mxSetField(G, 0, "nodes", nodes); - - - /* faces */ - const char *n3[] = {"num", "neighbors", "nodes", "numNodes", "nodePos", "tag"}; - mxArray *faces = mxCreateStructMatrix(1,1,sizeof n3 / sizeof n3[0], n3); - - - mxSetField(faces, 0, "num", mxCreateDoubleScalar(grid->number_of_faces)); - - mxArray *faceneighbors = mxCreateNumericMatrix(grid->number_of_faces, 2, - mxINT32_CLASS, mxREAL); - { - int *iptr = mxGetData(faceneighbors); - for(j=0; j<2; ++j) - { - for (i=0; inumber_of_faces; ++i) - { - int ix = grid->face_neighbors[2*i+j]; - if (ix == -1) - { - iptr[i+grid->number_of_faces*j] = 0; - } - else - { - iptr[i+grid->number_of_faces*j] = ix+1; - } - } - } - } - mxSetField(faces, 0, "neighbors", faceneighbors); - mxArray *numnodes = mxCreateNumericMatrix(grid->number_of_faces, 1, - mxINT32_CLASS, mxREAL); - mxArray *nodepos = mxCreateNumericMatrix(grid->number_of_faces+1, 1, - mxINT32_CLASS, mxREAL); - { - int *iptr1 = mxGetData(numnodes); - int *iptr2 = mxGetData(nodepos); - iptr2[0] = 1; - for (i=0; inumber_of_faces; ++i) - { - iptr1[i] = grid->face_ptr[i+1] - grid->face_ptr[i]; - iptr2[i+1] = iptr2[i] + iptr1[i]; - } - } - mxSetField(faces, 0, "numNodes", numnodes); - mxSetField(faces, 0, "nodePos", nodepos); - - mxArray *tags = mxCreateNumericMatrix(grid->number_of_faces, 1, - mxINT32_CLASS, mxREAL); - { - int *iptr = mxGetData(tags); - for (i = 0; i < grid->number_of_faces; ++i) - { - iptr[i] = grid->face_tag[i] + 1; - } - } - mxSetField(faces, 0, "tag", tags); - - - const char *n4[] = {"num", "faces", "facePos", "indexMap"}; - mxArray *cells = mxCreateStructMatrix(1,1,sizeof n4 / sizeof n4[0], n4); - - mxSetField(cells, 0, "num", mxCreateDoubleScalar(grid->number_of_cells)); - - mxArray *map = mxCreateNumericMatrix(grid->number_of_cells, 1, - mxINT32_CLASS, mxREAL); - { - int *iptr = mxGetData(map); - for(i=0; inumber_of_cells; ++i) - { - iptr[i] = grid->local_cell_index[i]+1; - } - } - mxSetField(cells, 0, "indexMap", map); - - mxArray *facepos = mxCreateNumericMatrix(grid->number_of_cells+1, 1, - mxINT32_CLASS, mxREAL); - { - int *iptr = mxGetData(facepos); - for(i=0; inumber_of_cells+1; ++i) - { - iptr[i] = 0; - } - for (i=0; i<2*grid->number_of_faces; ++i) - { - int c=grid->face_neighbors[i]; - if(c != -1) - { - iptr[c+1]++; - } - } - iptr[0] = 1; - for(i=0; inumber_of_cells; ++i) - { - iptr[i+1] += iptr[i]; - } - } - mxSetField(cells, 0, "facePos", facepos); - - - - int *counter = calloc(grid->number_of_cells, sizeof(*counter)); - int num_half_faces = 0; - { - int *iptr = mxGetData(facepos); - for(i=0; inumber_of_cells; ++i) - { - counter[i] = num_half_faces; - num_half_faces += iptr[i+1]-iptr[i]; - } - } - - mxArray *cellfaces = mxCreateNumericMatrix(num_half_faces, 2, - mxINT32_CLASS, mxREAL); - { - int *iptr = mxGetData(cellfaces); - int c1, c2, cf_tag; - for (i=0; inumber_of_faces; ++i) - { - cf_tag = 2 * grid->face_tag[i] + 1; /* [1, 3, 5] */ - c1 = grid->face_neighbors[2*i+0]; - c2 = grid->face_neighbors[2*i+1]; - if(c1 != -1) { - iptr[counter[c1] + 0*num_half_faces] = i+1; - iptr[counter[c1] + 1*num_half_faces] = cf_tag + 1; /* out */ - counter[c1] += 1; - } - if(c2 != -1) { - iptr[counter[c2] + 0*num_half_faces] = i+1; - iptr[counter[c2] + 1*num_half_faces] = cf_tag + 0; /* in */ - counter[c2] += 1; - } - } - } - mxSetField(cells, 0, "faces", cellfaces); - - mxSetField(G, 0, "cells", cells); - - int n = grid->face_ptr[grid->number_of_faces]; - - mxArray *facenodes = mxCreateNumericMatrix(n, 1, - mxINT32_CLASS, mxREAL); - { - int *iptr = mxGetData(facenodes); - for (i=0; iface_nodes[i]+1; - } - } - mxSetField(faces, 0, "nodes", facenodes); - - mxSetField(G, 0, "faces", faces); - - - free(counter); - - /* cartDims */ - mxArray *cartDims = mxCreateDoubleMatrix(1, 3, mxREAL); - ptr = mxGetPr(cartDims); - ptr[0] = grid->dimensions[0]; - ptr[1] = grid->dimensions[1]; - ptr[2] = grid->dimensions[2]; - - mxSetField(G, 0, "cartDims", cartDims); - - /* type */ - mxArray *type = mxCreateCellMatrix(1, 1); - mxSetCell(type, 0, mxCreateString("processgrid")); - - mxSetField(G, 0, "type", type); - - out[0] = G; -} - - - -/* Gateway routine for Matlab mex function. */ -/*-------------------------------------------------------*/ -void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) -{ - /* Set up data passed from Matlab */ - struct grdecl g; - struct processed_grid out; - double tolerance = 0.0; - - mx_init_grdecl(&g, prhs[0]); - if (nrhs == 2) - { - tolerance = mxGetScalar (prhs[1]); - } - - process_grdecl(&g, tolerance, &out); - - - if (plhs >0) - { - /* write to matlab struct */ - fill_grid(plhs, &out); - } - - - /* Free whatever was allocated in initGrdecl. */ - free_processed_grid(&out); -} -#endif