[MOAB-dev] r1933 - MOAB/trunk/tools/mbcoupler
tautges at mcs.anl.gov
tautges at mcs.anl.gov
Thu Jun 26 14:09:38 CDT 2008
Author: tautges
Date: 2008-06-26 14:09:38 -0500 (Thu, 26 Jun 2008)
New Revision: 1933
Modified:
MOAB/trunk/tools/mbcoupler/MBCoupler.cpp
MOAB/trunk/tools/mbcoupler/MBElemUtil.cpp
MOAB/trunk/tools/mbcoupler/MBElemUtil.hpp
MOAB/trunk/tools/mbcoupler/findpt_test.cpp
MOAB/trunk/tools/mbcoupler/types.h
Log:
Hooking up point interpolation from James Lottes' higher-order point location, from code written by Alvaro (thanks Alvaro!). Point location looks to be working, at least in the case where all points are local.
Modified: MOAB/trunk/tools/mbcoupler/MBCoupler.cpp
===================================================================
--- MOAB/trunk/tools/mbcoupler/MBCoupler.cpp 2008-06-26 18:39:15 UTC (rev 1932)
+++ MOAB/trunk/tools/mbcoupler/MBCoupler.cpp 2008-06-26 19:09:38 UTC (rev 1933)
@@ -394,7 +394,7 @@
//test to find out in which hex the point is
// get natural coordinates
- MBElemUtil::nat_coords_trilinear_hex(&coords_vert[0], MBCartVect(xyz),
+ MBElemUtil::nat_coords_trilinear_hex2(&coords_vert[0], MBCartVect(xyz),
tmp_nat_coords, 1e-10);
if (-1.0 <= tmp_nat_coords[0] && tmp_nat_coords[0] <= 1.0 &&
-1.0 <= tmp_nat_coords[1] && tmp_nat_coords[1] <= 1.0 &&
Modified: MOAB/trunk/tools/mbcoupler/MBElemUtil.cpp
===================================================================
--- MOAB/trunk/tools/mbcoupler/MBElemUtil.cpp 2008-06-26 18:39:15 UTC (rev 1932)
+++ MOAB/trunk/tools/mbcoupler/MBElemUtil.cpp 2008-06-26 19:09:38 UTC (rev 1933)
@@ -2,6 +2,7 @@
#include <assert.h>
#include "MBElemUtil.hpp"
+#include "types.h"
//
// nat_coords_trilinear_hex
@@ -140,6 +141,49 @@
}
+//
+// nat_coords_trilinear_hex2
+// Duplicate functionality of nat_coords_trilinear_hex using hex_findpt
+//
+void MBElemUtil::nat_coords_trilinear_hex2(MBCartVect hex[8],
+ MBCartVect xyz,
+ MBCartVect &ncoords,
+ double etol)
+
+{
+ const int ndim = 3;
+ const int nverts = 8;
+ const int vertMap[nverts] = {0,1,3,2, 4,5,7,6}; //Map from nat to lex ordering
+
+ const int n = 2; //linear
+ real coords[ndim*nverts]; //buffer
+
+ real *xm[ndim];
+ for(int i=0; i<ndim; i++)
+ xm[i] = coords + i*nverts;
+
+ //stuff hex into coords
+ for(int i=0; i<nverts; i++){
+ real vcoord[ndim];
+ hex[i].get(vcoord);
+
+ for(int d=0; d<ndim; d++)
+ coords[d*nverts + vertMap[i]] = vcoord[d];
+
+ }
+
+ double dist = 0.0;
+ MBElemUtil::hex_findpt(xm, n, xyz, ncoords, dist);
+ if (3*EPS < dist) {
+ // outside element, set extremal values to something outside range
+ for (int j = 0; j < 3; j++) {
+ if (ncoords[j] < (-1.0-etol) || ncoords[j] > (1.0+etol))
+ ncoords[j] *= 10;
+ }
+ }
+
+}
+
bool MBElemUtil::point_in_trilinear_hex(MBCartVect *hex,
MBCartVect xyz,
double etol)
Modified: MOAB/trunk/tools/mbcoupler/MBElemUtil.hpp
===================================================================
--- MOAB/trunk/tools/mbcoupler/MBElemUtil.hpp 2008-06-26 18:39:15 UTC (rev 1932)
+++ MOAB/trunk/tools/mbcoupler/MBElemUtil.hpp 2008-06-26 19:09:38 UTC (rev 1933)
@@ -19,7 +19,14 @@
MBCartVect box_max,
double etol);
+ //wrapper to hex_findpt
+ static void nat_coords_trilinear_hex2(MBCartVect*,
+ MBCartVect,
+ MBCartVect&,
+ double);
+
+
static void hex_findpt(double *xm[3],
int n,
MBCartVect xyz,
Modified: MOAB/trunk/tools/mbcoupler/findpt_test.cpp
===================================================================
--- MOAB/trunk/tools/mbcoupler/findpt_test.cpp 2008-06-26 18:39:15 UTC (rev 1932)
+++ MOAB/trunk/tools/mbcoupler/findpt_test.cpp 2008-06-26 19:09:38 UTC (rev 1933)
@@ -6,9 +6,10 @@
extern "C"{
#include "errmem.h" //for tmalloc, convenient but not C++
+#include "types.h"
}
-int main()
+void test_hex_findpt()
{
MBCartVect xyz(.5,.3,.4);
@@ -29,14 +30,14 @@
//Stuff xm with sample data
for(int k=0; k<n; k++){
for(int j=0; j<n; j++){
- for(int i=0; i<n; i++){
+ for(int i=0; i<n; i++){
- xm[0][node] = i*scale;
- xm[1][node] = j*scale;
- xm[2][node] = k*scale;
-
- node++;
- }
+ xm[0][node] = i*scale;
+ xm[1][node] = j*scale;
+ xm[2][node] = k*scale;
+
+ node++;
+ }
}
}
@@ -47,3 +48,40 @@
" distance: "<< dist << endl;
}
+
+
+
+void test_nat_coords_trilinear_hex2()
+{
+ MBCartVect hex[8];
+ MBCartVect xyz(.5,.3,.4);
+ MBCartVect ncoords;;
+ double etol;
+
+ MBElemUtil u;
+
+ //Make our sample hex the unit cube [0,1]**3
+ hex[0] = MBCartVect(0,0,0);
+ hex[1] = MBCartVect(1,0,0);
+ hex[2] = MBCartVect(1,1,0);
+ hex[3] = MBCartVect(0,1,0);
+ hex[4] = MBCartVect(0,0,1);
+ hex[5] = MBCartVect(1,0,1);
+ hex[6] = MBCartVect(1,1,1);
+ hex[7] = MBCartVect(0,1,1);
+
+ etol = .1 ; //ignored by nat_coords
+
+ u.nat_coords_trilinear_hex2(hex, xyz, ncoords, etol);
+
+ cout << "Coords of " << xyz << " are: "<< ncoords << endl;
+
+}
+
+
+
+int main(){
+
+ test_nat_coords_trilinear_hex2();
+ //test_hex_findpt();
+}
Modified: MOAB/trunk/tools/mbcoupler/types.h
===================================================================
--- MOAB/trunk/tools/mbcoupler/types.h 2008-06-26 18:39:15 UTC (rev 1932)
+++ MOAB/trunk/tools/mbcoupler/types.h 2008-06-26 19:09:38 UTC (rev 1933)
@@ -1,5 +1,6 @@
#ifndef TYPES_H
#define TYPES_H
+#include "float.h"
/* integer type to use for everything */
#if defined(USE_LONG)
More information about the moab-dev
mailing list