1 /**
2  * MOAB, a Mesh-Oriented datABase, is a software component for creating,
3  * storing and accessing finite element mesh data.
4  *
5  * Copyright 2004 Sandia Corporation.  Under the terms of Contract
6  * DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government
7  * retains certain rights in this software.
8  *
9  * This library is free software; you can redistribute it and/or
10  * modify it under the terms of the GNU Lesser General Public
11  * License as published by the Free Software Foundation; either
12  * version 2.1 of the License, or (at your option) any later version.
13  *
14  */
15 
16 //-------------------------------------------------------------------------
17 // Purpose       : This file contains utility functions for use in MOAB
18 //
19 // Special Notes : This is a pure virtual class, to prevent instantiation.
20 //                 All functions are static, called like this:
21 //                 Util::function_name();
22 //-------------------------------------------------------------------------
23 #include "moab/Util.hpp"
24 #include "moab/Interface.hpp"
25 #include <assert.h>
26 #include <algorithm>
27 #include <limits>
28 #if defined(_MSC_VER) || defined(__MINGW32__)
29 #  include <float.h>
30 #  define finite(A) _finite(A)
31 #ifndef MOAB_HAVE_FINITE
32 #define MOAB_HAVE_FINITE
33 #endif
34 #endif
35 
36 namespace moab {
37 
38 //! temporary normal function for MBEntities.  This should be moved to
39 //! an appropriate MB algorithms file
40 
normal(Interface * MB,EntityHandle handle,double & x,double & y,double & z)41 void Util::normal(Interface* MB, EntityHandle handle, double& x, double& y, double& z)
42 {
43    // get connectivity
44    const EntityHandle *connectivity = NULL;
45    int number_nodes = 0;
46    // TODO make the return value nonvoid
47    ErrorCode rval = MB->get_connectivity(handle, connectivity, number_nodes, true);
48    MB_CHK_SET_ERR_RET(rval, "can't get_connectivity");
49    assert(number_nodes >= 3);
50 
51    // get_coordinates
52    double coords[3][3];
53    MB->get_coords(&(connectivity[0]), 1, coords[0]);
54    MB->get_coords(&(connectivity[1]), 1, coords[1]);
55    MB->get_coords(&(connectivity[2]), 1, coords[2]);
56 
57    double vecs[2][3];
58    vecs[0][0] = coords[1][0] - coords[0][0];
59    vecs[0][1] = coords[1][1] - coords[0][1];
60    vecs[0][2] = coords[1][2] - coords[0][2];
61    vecs[1][0] = coords[2][0] - coords[0][0];
62    vecs[1][1] = coords[2][1] - coords[0][1];
63    vecs[1][2] = coords[2][2] - coords[0][2];
64 
65    x = vecs[0][1] * vecs[1][2] - vecs[0][2] * vecs[1][1];
66    y = vecs[0][2] * vecs[1][0] - vecs[0][0] * vecs[1][2];
67    z = vecs[0][0] * vecs[1][1] - vecs[0][1] * vecs[1][0];
68 
69    double mag = sqrt(x*x + y*y + z*z);
70    if(mag > std::numeric_limits<double>::epsilon())
71    {
72      x /= mag;
73      y /= mag;
74      z /= mag;
75    }
76 }
77 
centroid(Interface * MB,EntityHandle handle,CartVect & coord)78 void Util::centroid(Interface *MB, EntityHandle handle, CartVect &coord)
79 {
80    const EntityHandle *connectivity = NULL;
81    int number_nodes = 0;
82    // TODO make the return value nonvoid
83    ErrorCode rval = MB->get_connectivity(handle, connectivity, number_nodes, true);
84    MB_CHK_SET_ERR_RET(rval, "can't get_connectivity");
85 
86    coord[0]=coord[1]=coord[2]=0.0;
87 
88    for(int i = 0; i< number_nodes; i++)
89    {
90       double node_coords[3];
91       MB->get_coords(&(connectivity[i]), 1, node_coords);
92 
93       coord[0]+=node_coords[0];
94       coord[1]+=node_coords[1];
95       coord[2]+=node_coords[2];
96    }
97 
98    coord[0]/=(double)number_nodes;
99    coord[1]/=(double)number_nodes;
100    coord[2]/=(double)number_nodes;
101 }
102 
103 /*//This function calculates the coordinates for the centers of each edges of the entity specified by handle. The coordinates are returned in the list coords_list
104 void Util::edge_centers(Interface *MB, EntityHandle handle, std::vector<Coord> &coords_list)
105 {
106   MB canon_tool(MB);
107   EntityType type;
108   int i = 0;
109   int number_nodes = 0;
110   double coords[2][3];
111   const EntityHandle *connectivity;
112 
113   MB->get_connectivity(handle, connectivity, number_nodes,true);
114 
115   MB->type_from_handle(handle,type);
116 
117   const struct MBCN::ConnMap* conn_map = &(canon_tool.mConnectivityMap[type][0]); //get edge sub_elements
118 
119   coords_list.resize(conn_map->num_sub_elements);
120 
121   for(i = 0; i<conn_map->num_sub_elements; i++)
122   {
123 
124     MB->get_coords(connectivity[conn_map->conn[i][0]], coords[0]);
125     MB->get_coords(connectivity[conn_map->conn[i][1]], coords[1]);
126 
127     coords_list[i].x = (coords[0][0] + coords[1][0])/2.0;
128     coords_list[i].y = (coords[0][1] + coords[1][1])/2.0;
129     coords_list[i].z = (coords[0][2] + coords[1][2])/2.0;
130   }
131 }
132 */
133 
134 /*
135 void Util::face_centers(Interface *MB, EntityHandle handle, std::vector<Coord> &coords_list)
136 {
137   MB canon_tool(MB);
138   EntityType type;
139   int i = 0;
140   int number_nodes = 0;
141   double node_coords[3];
142   const EntityHandle *connectivity;
143 
144   MB->get_connectivity(handle, connectivity, number_nodes,true);
145 
146   MB->type_from_handle(handle,type);
147 
148   const struct MBCN::ConnMap* conn_map = &(canon_tool.mConnectivityMap[type][1]); //get face sub_elements
149 
150   coords_list.resize(conn_map->num_sub_elements);
151 
152   for(i = 0; i<conn_map->num_sub_elements;i++)
153   {
154     int number_nodes_per_element = conn_map->num_nodes_per_sub_element[i];
155 
156     for(int j = 0; j<number_nodes_per_element; j++)
157     {
158       MB->get_coords(connectivity[conn_map->conn[i][j]], node_coords);
159 
160       coords_list[i].x+=node_coords[0];
161       coords_list[i].y+=node_coords[1];
162       coords_list[i].z+=node_coords[2];
163     }
164 
165     coords_list[i].x/=(double)number_nodes_per_element;
166     coords_list[i].y/=(double)number_nodes_per_element;
167     coords_list[i].z/=(double)number_nodes_per_element;
168   }
169 }
170 */
171 
172 // Explicit template specializations
173 template bool Util::is_finite<double>(double value);
174 template bool Util::is_finite<int>(int value);
175 template bool Util::is_finite<unsigned int>(unsigned int value);
176 template bool Util::is_finite<long>(long value);
177 
178 } // namespace moab
179 
180