1 //============================================================================ 2 // Copyright (c) Kitware, Inc. 3 // All rights reserved. 4 // See LICENSE.txt for details. 5 // 6 // This software is distributed WITHOUT ANY WARRANTY; without even 7 // the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR 8 // PURPOSE. See the above copyright notice for more information. 9 //============================================================================ 10 // Copyright (c) 2016, Los Alamos National Security, LLC 11 // All rights reserved. 12 // 13 // Copyright 2016. Los Alamos National Security, LLC. 14 // This software was produced under U.S. Government contract DE-AC52-06NA25396 15 // for Los Alamos National Laboratory (LANL), which is operated by 16 // Los Alamos National Security, LLC for the U.S. Department of Energy. 17 // The U.S. Government has rights to use, reproduce, and distribute this 18 // software. NEITHER THE GOVERNMENT NOR LOS ALAMOS NATIONAL SECURITY, LLC 19 // MAKES ANY WARRANTY, EXPRESS OR IMPLIED, OR ASSUMES ANY LIABILITY FOR THE 20 // USE OF THIS SOFTWARE. If software is modified to produce derivative works, 21 // such modified software should be clearly marked, so as not to confuse it 22 // with the version available from LANL. 23 // 24 // Additionally, redistribution and use in source and binary forms, with or 25 // without modification, are permitted provided that the following conditions 26 // are met: 27 // 28 // 1. Redistributions of source code must retain the above copyright notice, 29 // this list of conditions and the following disclaimer. 30 // 2. Redistributions in binary form must reproduce the above copyright notice, 31 // this list of conditions and the following disclaimer in the documentation 32 // and/or other materials provided with the distribution. 33 // 3. Neither the name of Los Alamos National Security, LLC, Los Alamos 34 // National Laboratory, LANL, the U.S. Government, nor the names of its 35 // contributors may be used to endorse or promote products derived from 36 // this software without specific prior written permission. 37 // 38 // THIS SOFTWARE IS PROVIDED BY LOS ALAMOS NATIONAL SECURITY, LLC AND 39 // CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, 40 // BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS 41 // FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL LOS ALAMOS 42 // NATIONAL SECURITY, LLC OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, 43 // INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, 44 // BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF 45 // USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY 46 // THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT 47 // (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF 48 // THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 49 //============================================================================ 50 51 // This code is based on the algorithm presented in the paper: 52 // “Parallel Peak Pruning for Scalable SMP Contour Tree Computation.” 53 // Hamish Carr, Gunther Weber, Christopher Sewell, and James Ahrens. 54 // Proceedings of the IEEE Symposium on Large Data Analysis and Visualization 55 // (LDAV), October 2016, Baltimore, Maryland. 56 57 //======================================================================================= 58 // 59 // COMMENTS: 60 // 61 // This functor replaces a parallel loop examining neighbours - again, for arbitrary 62 // meshes, it needs to be a reduction, but for regular meshes, it's faster this way. 63 // 64 // Any vector needed by the functor for lookup purposes will be passed as a parameter to 65 // the constructor and saved, with the actual function call being the operator () 66 // 67 // Vectors marked I/O are intrinsically risky unless there is an algorithmic guarantee 68 // that the read/writes are completely independent - which for our case actually occurs 69 // The I/O vectors should therefore be justified in comments both here & in caller 70 // 71 //======================================================================================= 72 73 #ifndef vtkm_worklet_contourtree_mesh3d_dem_vertex_starter_h 74 #define vtkm_worklet_contourtree_mesh3d_dem_vertex_starter_h 75 76 #include <iostream> 77 #include <vtkm/worklet/WorkletMapField.h> 78 #include <vtkm/worklet/contourtree/Mesh3D_DEM_Triangulation_Macros.h> 79 #include <vtkm/worklet/contourtree/VertexValueComparator.h> 80 81 namespace vtkm 82 { 83 namespace worklet 84 { 85 namespace contourtree 86 { 87 88 // Worklet for setting initial chain maximum value 89 template <typename T> 90 class Mesh3D_DEM_VertexStarter : public vtkm::worklet::WorkletMapField 91 { 92 public: 93 using TagType = vtkm::List<T>; 94 95 using ControlSignature = void(FieldIn vertex, // (input) index of vertex 96 WholeArrayIn values, // (input) values within mesh 97 FieldOut chain, // (output) modify the chains 98 FieldOut linkMask); // (output) modify the mask 99 using ExecutionSignature = void(_1, _2, _3, _4); 100 using InputDomain = _1; 101 102 vtkm::Id nRows; // (input) number of rows in 3D 103 vtkm::Id nCols; // (input) number of cols in 3D 104 vtkm::Id nSlices; // (input) number of cols in 3D 105 bool ascending; // ascending or descending (join or split tree) 106 107 // Constructor 108 VTKM_EXEC_CONT Mesh3D_DEM_VertexStarter(vtkm::Id NRows,vtkm::Id NCols,vtkm::Id NSlices,bool Ascending)109 Mesh3D_DEM_VertexStarter(vtkm::Id NRows, vtkm::Id NCols, vtkm::Id NSlices, bool Ascending) 110 : nRows(NRows) 111 , nCols(NCols) 112 , nSlices(NSlices) 113 , ascending(Ascending) 114 { 115 } 116 117 // Locate the next vertex in direction indicated 118 template <typename InFieldPortalType> operator()119 VTKM_EXEC void operator()(const vtkm::Id& vertex, 120 const InFieldPortalType& values, 121 vtkm::Id& chain, 122 vtkm::Id& linkMask) const 123 { 124 VertexValueComparator<InFieldPortalType> lessThan(values); 125 vtkm::Id row = VERTEX_ROW_3D(vertex, nRows, nCols); 126 vtkm::Id col = VERTEX_COL_3D(vertex, nRows, nCols); 127 vtkm::Id slice = VERTEX_SLICE_3D(vertex, nRows, nCols); 128 129 vtkm::Id destination = vertex; 130 vtkm::Id mask = 0; 131 132 bool isLeft = (col == 0); 133 bool isRight = (col == nCols - 1); 134 bool isTop = (row == 0); 135 bool isBottom = (row == nRows - 1); 136 bool isFront = (slice == 0); 137 bool isBack = (slice == nSlices - 1); 138 139 // This order of processing must be maintained to match the LinkComponentCaseTables 140 // and to return the correct destination extremum 141 for (vtkm::Id edgeNo = (N_INCIDENT_EDGES_3D - 1); edgeNo >= 0; edgeNo--) 142 { 143 vtkm::Id nbr; 144 145 switch (edgeNo) 146 { 147 //////////////////////////////////////////////////////// 148 case 13: // down right back 149 if (isBack || isRight || isBottom) 150 break; 151 nbr = vertex + (nRows * nCols) + nCols + 1; 152 if (lessThan(vertex, nbr, ascending)) 153 break; 154 mask |= 0x2000; 155 destination = nbr; 156 break; 157 158 case 12: // down back 159 if (isBack || isBottom) 160 break; 161 nbr = vertex + (nRows * nCols) + nCols; 162 if (lessThan(vertex, nbr, ascending)) 163 break; 164 mask |= 0x1000; 165 destination = nbr; 166 break; 167 168 case 11: // right back 169 if (isBack || isRight) 170 break; 171 nbr = vertex + (nRows * nCols) + 1; 172 if (lessThan(vertex, nbr, ascending)) 173 break; 174 mask |= 0x800; 175 destination = nbr; 176 break; 177 178 case 10: // back 179 if (isBack) 180 break; 181 nbr = vertex + (nRows * nCols); 182 if (lessThan(vertex, nbr, ascending)) 183 break; 184 mask |= 0x400; 185 destination = nbr; 186 break; 187 188 case 9: // down right 189 if (isBottom || isRight) 190 break; 191 nbr = vertex + nCols + 1; 192 if (lessThan(vertex, nbr, ascending)) 193 break; 194 mask |= 0x200; 195 destination = nbr; 196 break; 197 198 case 8: // down 199 if (isBottom) 200 break; 201 nbr = vertex + nCols; 202 if (lessThan(vertex, nbr, ascending)) 203 break; 204 mask |= 0x100; 205 destination = nbr; 206 break; 207 208 case 7: // right 209 if (isRight) 210 break; 211 nbr = vertex + 1; 212 if (lessThan(vertex, nbr, ascending)) 213 break; 214 mask |= 0x80; 215 destination = nbr; 216 break; 217 218 case 6: // up left 219 if (isLeft || isTop) 220 break; 221 nbr = vertex - nCols - 1; 222 if (lessThan(vertex, nbr, ascending)) 223 break; 224 mask |= 0x40; 225 destination = nbr; 226 break; 227 228 case 5: // left 229 if (isLeft) 230 break; 231 nbr = vertex - 1; 232 if (lessThan(vertex, nbr, ascending)) 233 break; 234 mask |= 0x20; 235 destination = nbr; 236 break; 237 238 case 4: // left front 239 if (isLeft || isFront) 240 break; 241 nbr = vertex - (nRows * nCols) - 1; 242 if (lessThan(vertex, nbr, ascending)) 243 break; 244 mask |= 0x10; 245 destination = nbr; 246 break; 247 248 case 3: // front 249 if (isFront) 250 break; 251 nbr = vertex - (nRows * nCols); 252 if (lessThan(vertex, nbr, ascending)) 253 break; 254 mask |= 0x08; 255 destination = nbr; 256 break; 257 258 case 2: // up front 259 if (isTop || isFront) 260 break; 261 nbr = vertex - (nRows * nCols) - nCols; 262 if (lessThan(vertex, nbr, ascending)) 263 break; 264 mask |= 0x04; 265 destination = nbr; 266 break; 267 268 case 1: // up 269 if (isTop) 270 break; 271 nbr = vertex - nCols; 272 if (lessThan(vertex, nbr, ascending)) 273 break; 274 mask |= 0x02; 275 destination = nbr; 276 break; 277 278 case 0: // up left front 279 if (isTop || isLeft || isFront) 280 break; 281 nbr = vertex - (nRows * nCols) - nCols - 1; 282 if (lessThan(vertex, nbr, ascending)) 283 break; 284 mask |= 0x01; 285 destination = nbr; 286 break; 287 } // switch on edgeNo 288 } // per edge 289 290 linkMask = mask; 291 chain = destination; 292 } // operator() 293 }; // Mesh3D_DEM_VertexStarter 294 } 295 } 296 } 297 298 #endif 299