1 //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~//
2 // Copyright (c) Lawrence Livermore National Security, LLC and other Ascent
3 // Project developers. See top-level LICENSE AND COPYRIGHT files for dates and
4 // other details. No copyright assignment is required to contribute to Ascent.
5 //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~//
6
7
8 //-----------------------------------------------------------------------------
9 ///
10 /// file: ascent_data_adapter.cpp
11 ///
12 //-----------------------------------------------------------------------------
13 #include "ascent_vtkh_data_adapter.hpp"
14
15 // standard lib includes
16 #include <iostream>
17 #include <string.h>
18 #include <limits.h>
19 #include <cstdlib>
20 #include <sstream>
21 #include <type_traits>
22
23 // third party includes
24
25 // mpi
26 #ifdef ASCENT_MPI_ENABLED
27 #include <mpi.h>
28 #endif
29
30 // VTKm includes
31 #define VTKM_USE_DOUBLE_PRECISION
32 #include <vtkm/cont/DataSet.h>
33 #include <vtkm/cont/ArrayCopy.h>
34 #include <vtkm/cont/ArrayHandleExtractComponent.h>
35 #include <vtkm/cont/CoordinateSystem.h>
36 #include <vtkh/DataSet.hpp>
37
38 // other ascent includes
39 #include <ascent_logging.hpp>
40 #include <ascent_block_timer.hpp>
41 #include <ascent_mpi_utils.hpp>
42 #include <vtkh/utils/vtkm_array_utils.hpp>
43 #include <vtkh/utils/vtkm_dataset_info.hpp>
44
45 #include <conduit_blueprint.hpp>
46
47 using namespace std;
48 using namespace conduit;
49
50 //-----------------------------------------------------------------------------
51 // -- begin ascent:: --
52 //-----------------------------------------------------------------------------
53 namespace ascent
54 {
55
56 //-----------------------------------------------------------------------------
57 // -- begin detail:: --
58 //-----------------------------------------------------------------------------
59 namespace detail
60 {
61
topo_origin(const conduit::Node & n_topo)62 vtkm::Id3 topo_origin(const conduit::Node &n_topo)
63 {
64 vtkm::Id3 topo_origin(0,0,0);
65 // maintain backwards compatibility between
66 // i and i0 versions
67 if(n_topo.has_path("elements/origin"))
68 {
69 const conduit::Node &origin = n_topo["elements/origin"];
70
71 if(origin.has_path("i"))
72 {
73 topo_origin[0] = n_topo["elements/origin/i"].to_int32();
74 }
75 if(origin.has_path("i0"))
76 {
77 topo_origin[0] = n_topo["elements/origin/i0"].to_int32();
78 }
79
80 if(origin.has_path("j"))
81 {
82 topo_origin[1] = n_topo["elements/origin/j"].to_int32();
83 }
84 if(origin.has_path("j0"))
85 {
86 topo_origin[1] = n_topo["elements/origin/j0"].to_int32();
87 }
88
89 if(origin.has_path("k"))
90 {
91 topo_origin[2] = n_topo["elements/origin/k"].to_int32();
92 }
93 if(origin.has_path("k0"))
94 {
95 topo_origin[2] = n_topo["elements/origin/k0"].to_int32();
96 }
97 }
98
99 return topo_origin;
100 }
101
102 template<typename T>
103 const T* GetNodePointer(const conduit::Node &node);
104
105 template<>
GetNodePointer(const conduit::Node & node)106 const float64* GetNodePointer<float64>(const conduit::Node &node)
107 {
108 return node.as_float64_ptr();
109 }
110
111 template<>
GetNodePointer(const conduit::Node & node)112 const float32* GetNodePointer<float32>(const conduit::Node &node)
113 {
114 return node.as_float32_ptr();
115 }
116
117 template<typename T>
CopyArray(vtkm::cont::ArrayHandle<T> & vtkm_handle,const T * vals_ptr,const int size,bool zero_copy)118 void CopyArray(vtkm::cont::ArrayHandle<T> &vtkm_handle, const T* vals_ptr, const int size, bool zero_copy)
119 {
120 vtkm::CopyFlag copy = vtkm::CopyFlag::On;
121 if(zero_copy)
122 {
123 copy = vtkm::CopyFlag::Off;
124 }
125
126 vtkm_handle = vtkm::cont::make_ArrayHandle(vals_ptr, size, copy);
127 }
128
129 template<typename T>
130 vtkm::cont::CoordinateSystem
GetExplicitCoordinateSystem(const conduit::Node & n_coords,const std::string name,int & ndims,bool zero_copy)131 GetExplicitCoordinateSystem(const conduit::Node &n_coords,
132 const std::string name,
133 int &ndims,
134 bool zero_copy)
135 {
136 int nverts = n_coords["values/x"].dtype().number_of_elements();
137 bool is_interleaved = blueprint::mcarray::is_interleaved(n_coords["values"]);
138
139 // some interleaved cases aren't working
140 // disabling this path until we find out what is going wrong.
141 is_interleaved = false;
142
143 ndims = 2;
144
145 // n_coords_conv holds contig data if we have stride-ed but
146 // non-interleaved values
147 Node n_coords_conv;
148
149 const T* x_coords_ptr = NULL;
150 const T* y_coords_ptr = NULL;
151 const T *z_coords_ptr = NULL;
152
153 // if we are an interleaved mcarray, or compact we can
154 // directly use the pointer with vtk-m.
155 // otherwise, we need to compact.
156
157 if(is_interleaved || n_coords["values/x"].is_compact())
158 {
159 x_coords_ptr = GetNodePointer<T>(n_coords["values/x"]);
160 }
161 else
162 {
163 n_coords["values/x"].compact_to(n_coords_conv["x"]);
164 x_coords_ptr = GetNodePointer<T>(n_coords_conv["x"]);
165 // since we had to copy and compact the data, we can't zero copy
166 zero_copy = false;
167 }
168
169 if(is_interleaved || n_coords["values/y"].is_compact())
170 {
171 y_coords_ptr = GetNodePointer<T>(n_coords["values/y"]);
172 }
173 else
174 {
175 n_coords["values/y"].compact_to(n_coords_conv["y"]);
176 y_coords_ptr = GetNodePointer<T>(n_coords_conv["y"]);
177 // since we had to copy and compact the data, we can't zero copy
178 zero_copy = false;
179 }
180
181 if(n_coords.has_path("values/z"))
182 {
183 ndims = 3;
184 if(is_interleaved || n_coords["values/z"].is_compact())
185 {
186 z_coords_ptr = GetNodePointer<T>(n_coords["values/z"]);
187 }
188 else
189 {
190 n_coords["values/z"].compact_to(n_coords_conv["z"]);
191 z_coords_ptr = GetNodePointer<T>(n_coords_conv["z"]);
192 // since we had to copy and compact the data, we can't zero copy
193 zero_copy = false;
194 }
195 }
196
197 if(!is_interleaved)
198 {
199 vtkm::cont::ArrayHandle<T> x_coords_handle;
200 vtkm::cont::ArrayHandle<T> y_coords_handle;
201 vtkm::cont::ArrayHandle<T> z_coords_handle;
202
203 detail::CopyArray(x_coords_handle, x_coords_ptr, nverts, zero_copy);
204 detail::CopyArray(y_coords_handle, y_coords_ptr, nverts, zero_copy);
205
206 if(ndims == 3)
207 {
208 detail::CopyArray(z_coords_handle, z_coords_ptr, nverts, zero_copy);
209 }
210 else
211 {
212 z_coords_handle.Allocate(nverts);
213 // This does not get initialized to zero
214 T *z = vtkh::GetVTKMPointer(z_coords_handle);
215 memset(z, 0.0, nverts * sizeof(T));
216 }
217
218 return vtkm::cont::CoordinateSystem(name,
219 make_ArrayHandleSOA(x_coords_handle,
220 y_coords_handle,
221 z_coords_handle));
222 }
223 else // NOTE: This case is disabled.
224 {
225 // we have interleaved coordinates x0,y0,z0,x1,y1,z1...
226 const T* coords_ptr = GetNodePointer<T>(n_coords["values/x"]);
227 vtkm::cont::ArrayHandle<vtkm::Vec<T, 3>> coords;
228 // we cannot zero copy 2D interleaved arrays into vtkm
229 if(ndims == 3 || true) // TODO: need way to detect 3d interleaved components that has
230 // only has xy in conduit
231 {
232 // this case was failing from Nyx + AMReX
233 // still haven't been able to reproduce with a simpler test
234 detail::CopyArray(coords, (vtkm::Vec<T,3>*)coords_ptr, nverts, zero_copy);
235 }
236 else
237 {
238 // 2D interleaved array case
239 vtkm::cont::ArrayHandle<T> x_coords_handle;
240 vtkm::cont::ArrayHandle<T> y_coords_handle;
241 vtkm::cont::ArrayHandle<T> z_coords_handle;
242
243 x_coords_handle.Allocate(nverts);
244 y_coords_handle.Allocate(nverts);
245 z_coords_handle.Allocate(nverts);
246
247 auto x_portal = x_coords_handle.WritePortal();
248 auto y_portal = y_coords_handle.WritePortal();
249
250 const T* coords_ptr = GetNodePointer<T>(n_coords["values/x"]);
251
252 T *z = (T*) vtkh::GetVTKMPointer(z_coords_handle);
253 memset(z, 0.0, nverts * sizeof(T));
254
255 for(int i = 0; i < nverts; ++i)
256 {
257 x_portal.Set(i, coords_ptr[i*2+0]);
258 y_portal.Set(i, coords_ptr[i*2+1]);
259 }
260
261 return vtkm::cont::CoordinateSystem(name,
262 make_ArrayHandleSOA(x_coords_handle,
263 y_coords_handle,
264 z_coords_handle));
265 }
266
267 return vtkm::cont::CoordinateSystem(name, coords);
268 }
269
270 }
271
272 template<typename T>
GetField(const conduit::Node & node,const std::string field_name,const std::string assoc_str,const std::string topo_str,bool zero_copy)273 vtkm::cont::Field GetField(const conduit::Node &node,
274 const std::string field_name,
275 const std::string assoc_str,
276 const std::string topo_str,
277 bool zero_copy)
278 {
279 vtkm::CopyFlag copy = vtkm::CopyFlag::On;
280 if(zero_copy)
281 {
282 copy = vtkm::CopyFlag::Off;
283 }
284 vtkm::cont::Field::Association vtkm_assoc = vtkm::cont::Field::Association::ANY;
285 if(assoc_str == "vertex")
286 {
287 vtkm_assoc = vtkm::cont::Field::Association::POINTS;
288 }
289 else if(assoc_str == "element")
290 {
291 vtkm_assoc = vtkm::cont::Field::Association::CELL_SET;
292 }
293 else
294 {
295 ASCENT_ERROR("Cannot add field association "<<assoc_str<<" from field "<<field_name);
296 }
297
298 int num_vals = node.dtype().number_of_elements();
299
300 const T *values_ptr = node.value();
301
302 vtkm::cont::Field field;
303 field = vtkm::cont::make_Field(field_name,
304 vtkm_assoc,
305 values_ptr,
306 num_vals,
307 copy);
308 return field;
309 }
310
311 template<typename T>
GetVectorField(T * values_ptr,const int num_vals,const std::string field_name,const std::string assoc_str,const std::string topo_str,bool zero_copy)312 vtkm::cont::Field GetVectorField(T *values_ptr,
313 const int num_vals,
314 const std::string field_name,
315 const std::string assoc_str,
316 const std::string topo_str,
317 bool zero_copy)
318 {
319 vtkm::CopyFlag copy = vtkm::CopyFlag::On;
320 if(zero_copy)
321 {
322 copy = vtkm::CopyFlag::Off;
323 }
324 vtkm::cont::Field::Association vtkm_assoc = vtkm::cont::Field::Association::ANY;
325 if(assoc_str == "vertex")
326 {
327 vtkm_assoc = vtkm::cont::Field::Association::POINTS;
328 }
329 else if(assoc_str == "element")
330 {
331 vtkm_assoc = vtkm::cont::Field::Association::CELL_SET;
332 }
333 else
334 {
335 ASCENT_ERROR("Cannot add vector field with association "
336 <<assoc_str<<" field_name "<<field_name);
337 }
338
339 vtkm::cont::Field field;
340 field = vtkm::cont::make_Field(field_name,
341 vtkm_assoc,
342 values_ptr,
343 num_vals,
344 copy);
345
346 return field;
347 }
348
349 //
350 // extract a vector from 3 separate arrays
351 //
352 template<typename T>
ExtractVector(vtkm::cont::DataSet * dset,const conduit::Node & u,const conduit::Node & v,const conduit::Node & w,const int num_vals,const int dims,const std::string field_name,const std::string assoc_str,const std::string topo_name,bool zero_copy)353 void ExtractVector(vtkm::cont::DataSet *dset,
354 const conduit::Node &u,
355 const conduit::Node &v,
356 const conduit::Node &w,
357 const int num_vals,
358 const int dims,
359 const std::string field_name,
360 const std::string assoc_str,
361 const std::string topo_name,
362 bool zero_copy)
363 {
364 // TODO: Do we need to fix this for striding?
365 // GetField<T> expects compact
366 if(dims != 2 && dims != 3)
367 {
368 ASCENT_ERROR("Extract vector: only 2 and 3 dims supported given "<<dims);
369 }
370
371 vtkm::cont::Field::Association vtkm_assoc = vtkm::cont::Field::Association::ANY;
372 if(assoc_str == "vertex")
373 {
374 vtkm_assoc = vtkm::cont::Field::Association::POINTS;
375 }
376 else if (assoc_str == "element")
377 {
378 vtkm_assoc = vtkm::cont::Field::Association::CELL_SET;
379 }
380 else
381 {
382 ASCENT_ERROR("Cannot add vector field with association "
383 <<assoc_str<<" field_name "<<field_name);
384 }
385
386 if(dims == 2)
387 {
388 const T *x_ptr = GetNodePointer<T>(u);
389 const T *y_ptr = GetNodePointer<T>(v);
390
391 vtkm::cont::ArrayHandle<T> x_handle;
392 vtkm::cont::ArrayHandle<T> y_handle;
393
394 // always zero copy because we are about to make a copy
395 detail::CopyArray(x_handle, x_ptr, num_vals, true);
396 detail::CopyArray(y_handle, y_ptr, num_vals, true);
397
398
399 auto composite = make_ArrayHandleSOA(x_handle,
400 y_handle);
401
402 vtkm::cont::ArrayHandle<vtkm::Vec<T,2>> interleaved_handle;
403 interleaved_handle.Allocate(num_vals);
404 // Calling this without forcing serial could cause serious problems
405 {
406 vtkm::cont::ScopedRuntimeDeviceTracker tracker(vtkm::cont::DeviceAdapterTagSerial{});
407 vtkm::cont::ArrayCopy(composite, interleaved_handle);
408 }
409
410 vtkm::cont::Field field(field_name, vtkm_assoc, interleaved_handle);
411 dset->AddField(field);
412 }
413
414 if(dims == 3)
415 {
416 const T *x_ptr = GetNodePointer<T>(u);
417 const T *y_ptr = GetNodePointer<T>(v);
418 const T *z_ptr = GetNodePointer<T>(w);
419
420 vtkm::cont::ArrayHandle<T> x_handle;
421 vtkm::cont::ArrayHandle<T> y_handle;
422 vtkm::cont::ArrayHandle<T> z_handle;
423
424 // always zero copy because we are about to make a copy
425 detail::CopyArray(x_handle, x_ptr, num_vals, true);
426 detail::CopyArray(y_handle, y_ptr, num_vals, true);
427 detail::CopyArray(z_handle, z_ptr, num_vals, true);
428
429 auto composite = make_ArrayHandleSOA(x_handle,
430 y_handle,
431 z_handle);
432
433 vtkm::cont::ArrayHandle<vtkm::Vec<T,3>> interleaved_handle;
434 interleaved_handle.Allocate(num_vals);
435 // Calling this without forcing serial could cause serious problems
436 {
437 vtkm::cont::ScopedRuntimeDeviceTracker tracker(vtkm::cont::DeviceAdapterTagSerial{});
438 vtkm::cont::ArrayCopy(composite, interleaved_handle);
439 }
440
441 vtkm::cont::Field field(field_name, vtkm_assoc, interleaved_handle);
442 dset->AddField(field);
443 }
444 }
445
446
VTKmCellShape(const std::string shape_type,vtkm::UInt8 & shape_id,vtkm::IdComponent & num_indices)447 void VTKmCellShape(const std::string shape_type,
448 vtkm::UInt8 &shape_id,
449 vtkm::IdComponent &num_indices)
450 {
451 shape_id = 0;
452 num_indices = 0;
453 if(shape_type == "tri")
454 {
455 shape_id = 5;
456 num_indices = 3;
457 }
458 else if(shape_type == "quad")
459 {
460 shape_id = 9;
461 num_indices = 4;
462 }
463 else if(shape_type == "tet")
464 {
465 shape_id = 10;
466 num_indices = 4;
467 }
468 else if(shape_type == "hex")
469 {
470 shape_id = 12;
471 num_indices = 8;
472 }
473 else if(shape_type == "point")
474 {
475 shape_id = 1;
476 num_indices = 1;
477 }
478 else if(shape_type == "line")
479 {
480 shape_id = 3;
481 num_indices = 2;
482 }
483 else
484 {
485 ASCENT_ERROR("Unsupported cell type "<<shape_type);
486 }
487 }
488
489 };
490 //-----------------------------------------------------------------------------
491 // -- end detail:: --
492 //-----------------------------------------------------------------------------
493
494 //-----------------------------------------------------------------------------
495 // VTKHDataAdapter public methods
496 //-----------------------------------------------------------------------------
497
498 VTKHCollection*
BlueprintToVTKHCollection(const conduit::Node & n,bool zero_copy)499 VTKHDataAdapter::BlueprintToVTKHCollection(const conduit::Node &n,
500 bool zero_copy)
501 {
502 // We must separate different topologies into
503 // different vtkh data sets
504
505 const int num_domains = n.number_of_children();
506
507 VTKHCollection *res = new VTKHCollection();
508 std::map<std::string, vtkh::DataSet> datasets;
509 vtkm::UInt64 cycle = 0;
510 double time = 0;
511
512 for(int i = 0; i < num_domains; ++i)
513 {
514 const conduit::Node &dom = n.child(i);
515 const std::vector<std::string> topo_names = dom["topologies"].child_names();
516
517 if(!dom.has_path("state/domain_id"))
518 {
519 ASCENT_ERROR("Must have a domain_id to convert blueprint to vtkh");
520 }
521
522 int domain_id = dom["state/domain_id"].to_int();
523
524 if(dom.has_path("state/cycle"))
525 {
526 cycle = dom["state/cycle"].to_uint64();
527 }
528
529 if(dom.has_path("state/time"))
530 {
531 time = dom["state/time"].to_float64();
532 }
533
534 for(int t = 0; t < topo_names.size(); ++t)
535 {
536 const std::string topo_name = topo_names[t];
537 vtkm::cont::DataSet *dset = BlueprintToVTKmDataSet(dom, zero_copy, topo_name);
538 datasets[topo_name].AddDomain(*dset,domain_id);
539 delete dset;
540 }
541
542 }
543
544 for(auto dset_it : datasets)
545 {
546 res->add(dset_it.second, dset_it.first);
547 }
548
549 return res;
550 }
551
552 //-----------------------------------------------------------------------------
553 vtkh::DataSet *
BlueprintToVTKHDataSet(const Node & node,const std::string & topo_name,bool zero_copy)554 VTKHDataAdapter::BlueprintToVTKHDataSet(const Node &node,
555 const std::string &topo_name,
556 bool zero_copy)
557 {
558
559 // treat everything as a multi-domain data set
560
561 vtkh::DataSet *res = new vtkh::DataSet;
562
563 int num_domains = 0;
564
565 // get the number of domains and check for id consistency
566 num_domains = node.number_of_children();
567
568 for(int i = 0; i < num_domains; ++i)
569 {
570 const conduit::Node &dom = node.child(i);
571 vtkm::cont::DataSet *dset = VTKHDataAdapter::BlueprintToVTKmDataSet(dom,
572 zero_copy,
573 topo_name);
574 int domain_id = dom["state/domain_id"].to_int();
575
576 if(dom.has_path("state/cycle"))
577 {
578 vtkm::UInt64 cycle = dom["state/cycle"].to_uint64();
579 res->SetCycle(cycle);
580 }
581
582 res->AddDomain(*dset,domain_id);
583 // vtk-m will shallow copy the data assoced with dset
584 // clean up our copy
585 delete dset;
586
587 }
588 return res;
589 }
590
591 //-----------------------------------------------------------------------------
592 vtkh::DataSet *
VTKmDataSetToVTKHDataSet(vtkm::cont::DataSet * dset)593 VTKHDataAdapter::VTKmDataSetToVTKHDataSet(vtkm::cont::DataSet *dset)
594 {
595 // wrap a single VTKm data set into a VTKH dataset
596 vtkh::DataSet *res = new vtkh::DataSet;
597 int domain_id = 0; // TODO, MPI_TASK_ID ?
598 res->AddDomain(*dset,domain_id);
599 return res;
600 }
601
602 //-----------------------------------------------------------------------------
603 vtkm::cont::DataSet *
BlueprintToVTKmDataSet(const Node & node,bool zero_copy,const std::string & topo_name_str)604 VTKHDataAdapter::BlueprintToVTKmDataSet(const Node &node,
605 bool zero_copy,
606 const std::string &topo_name_str)
607 {
608 vtkm::cont::DataSet * result = NULL;
609
610 std::string topo_name = topo_name_str;
611
612 // we must find the topolgy they asked for
613 if(!node["topologies"].has_child(topo_name))
614 {
615 ASCENT_ERROR("Invalid topology name: " << topo_name);
616 }
617
618 // as long as mesh blueprint verify true, we access data without fear.
619
620 const Node &n_topo = node["topologies"][topo_name];
621 string mesh_type = n_topo["type"].as_string();
622
623 string coords_name = n_topo["coordset"].as_string();
624 const Node &n_coords = node["coordsets"][coords_name];
625
626 int neles = 0;
627 int nverts = 0;
628
629 if( mesh_type == "uniform")
630 {
631 result = UniformBlueprintToVTKmDataSet(coords_name,
632 n_coords,
633 topo_name,
634 n_topo,
635 neles,
636 nverts);
637 }
638 else if(mesh_type == "rectilinear")
639 {
640 result = RectilinearBlueprintToVTKmDataSet(coords_name,
641 n_coords,
642 topo_name,
643 n_topo,
644 neles,
645 nverts,
646 zero_copy);
647
648 }
649 else if(mesh_type == "structured")
650 {
651 result = StructuredBlueprintToVTKmDataSet(coords_name,
652 n_coords,
653 topo_name,
654 n_topo,
655 neles,
656 nverts,
657 zero_copy);
658 }
659 else if( mesh_type == "unstructured")
660 {
661 result = UnstructuredBlueprintToVTKmDataSet(coords_name,
662 n_coords,
663 topo_name,
664 n_topo,
665 neles,
666 nverts,
667 zero_copy);
668 }
669 else
670 {
671 ASCENT_ERROR("Unsupported topology/type:" << mesh_type);
672 }
673
674
675 if(node.has_child("fields"))
676 {
677 // add all of the fields:
678 NodeConstIterator itr = node["fields"].children();
679 while(itr.has_next())
680 {
681
682 const Node &n_field = itr.next();
683 std::string field_name = itr.name();
684 if(n_field["topology"].as_string() != topo_name)
685 {
686 // these are not the fields we are looking for
687 continue;
688 }
689
690 // skip vector fields for now, we need to add
691 // more logic to AddField
692 const int num_children = n_field["values"].number_of_children();
693
694 if(num_children == 0 || num_children == 1)
695 {
696
697 AddField(field_name,
698 n_field,
699 topo_name,
700 neles,
701 nverts,
702 result,
703 zero_copy);
704 }
705 else if(num_children == 2 )
706 {
707 AddVectorField(field_name,
708 n_field,
709 topo_name,
710 neles,
711 nverts,
712 result,
713 2,
714 zero_copy);
715 }
716 else if(num_children == 3 )
717 {
718 AddVectorField(field_name,
719 n_field,
720 topo_name,
721 neles,
722 nverts,
723 result,
724 3,
725 zero_copy);
726 }
727 else
728 {
729 ASCENT_INFO("skipping field "<<field_name<<" with "<<num_children<<" comps");
730 }
731 }
732 }
733 return result;
734 }
735
736
737 //-----------------------------------------------------------------------------
738 class ExplicitArrayHelper
739 {
740 public:
741 // Helper function to create explicit coordinate arrays for vtkm data sets
CreateExplicitArrays(vtkm::cont::ArrayHandle<vtkm::UInt8> & shapes,vtkm::cont::ArrayHandle<vtkm::IdComponent> & num_indices,const std::string & shape_type,const vtkm::Id & conn_size,vtkm::IdComponent & dimensionality,int & neles)742 void CreateExplicitArrays(vtkm::cont::ArrayHandle<vtkm::UInt8> &shapes,
743 vtkm::cont::ArrayHandle<vtkm::IdComponent> &num_indices,
744 const std::string &shape_type,
745 const vtkm::Id &conn_size,
746 vtkm::IdComponent &dimensionality,
747 int &neles)
748 {
749 vtkm::UInt8 shape_id = 0;
750 vtkm::IdComponent indices= 0;
751 if(shape_type == "tri")
752 {
753 shape_id = 3;
754 indices = 3;
755 // note: vtkm cell dimensions are topological
756 dimensionality = 2;
757 }
758 else if(shape_type == "quad")
759 {
760 shape_id = 9;
761 indices = 4;
762 // note: vtkm cell dimensions are topological
763 dimensionality = 2;
764 }
765 else if(shape_type == "tet")
766 {
767 shape_id = 10;
768 indices = 4;
769 dimensionality = 3;
770 }
771 else if(shape_type == "hex")
772 {
773 shape_id = 12;
774 indices = 8;
775 dimensionality = 3;
776 }
777 else if(shape_type == "points")
778 {
779 shape_id = 1;
780 indices = 1;
781 dimensionality = 1;
782 }
783 // TODO: Not supported in blueprint yet ...
784 // else if(shape_type == "wedge")
785 // {
786 // shape_id = 13;
787 // indices = 6;
788 // dimensionality = 3;
789 // }
790 // else if(shape_type == "pyramid")
791 // {
792 // shape_id = 14;
793 // indices = 5;
794 // dimensionality = 3;
795 // }
796 else
797 {
798 ASCENT_ERROR("Unsupported element shape " << shape_type);
799 }
800
801 if(conn_size < indices)
802 ASCENT_ERROR("Connectivity array size " <<conn_size << " must be at least size " << indices);
803 if(conn_size % indices != 0)
804 ASCENT_ERROR("Connectivity array size " <<conn_size << " be evenly divided by indices size" << indices);
805
806 const vtkm::Id num_shapes = conn_size / indices;
807
808 neles = num_shapes;
809
810 shapes.Allocate(num_shapes);
811 num_indices.Allocate(num_shapes);
812
813 // We could memset these and then zero copy them but that
814 // would make us responsible for the data. If we just create
815 // them, smart pointers will automatically delete them.
816 // Hopefull the compiler turns this into a memset.
817
818 const vtkm::UInt8 shape_value = shape_id;
819 const vtkm::IdComponent indices_value = indices;
820 auto shapes_portal = shapes.WritePortal();
821 auto num_indices_portal = num_indices.WritePortal();
822 #ifdef ASCENT_USE_OPENMP
823 #pragma omp parallel for
824 #endif
825 for (int i = 0; i < num_shapes; ++i)
826 {
827 shapes_portal.Set(i, shape_value);
828 num_indices_portal.Set(i, indices_value);
829 }
830 }
831 };
832 //-----------------------------------------------------------------------------
833
834 vtkm::cont::DataSet *
UniformBlueprintToVTKmDataSet(const std::string & coords_name,const Node & n_coords,const std::string & topo_name,const Node & n_topo,int & neles,int & nverts)835 VTKHDataAdapter::UniformBlueprintToVTKmDataSet
836 (const std::string &coords_name, // input string with coordset name
837 const Node &n_coords, // input mesh bp coordset (assumed uniform)
838 const std::string &topo_name, // input string with topo name
839 const Node &n_topo, // input mesh bp topo
840 int &neles, // output, number of eles
841 int &nverts) // output, number of verts
842 {
843 //
844 // blueprint uniform coord set provides:
845 //
846 // dims/{i,j,k}
847 // origin/{x,y,z} (optional)
848 // spacing/{dx,dy,dz} (optional)
849
850 //Create implicit vtkm coordinate system
851 vtkm::cont::DataSet *result = new vtkm::cont::DataSet();
852
853 const Node &n_dims = n_coords["dims"];
854
855 int dims_i = n_dims["i"].to_int();
856 int dims_j = n_dims["j"].to_int();
857 int dims_k = 1;
858
859 bool is_2d = true;
860
861 // check for 3d
862 if(n_dims.has_path("k"))
863 {
864 dims_k = n_dims["k"].to_int();
865 is_2d = false;
866 }
867
868
869
870 float64 origin_x = 0.0;
871 float64 origin_y = 0.0;
872 float64 origin_z = 0.0;
873
874
875 float64 spacing_x = 1.0;
876 float64 spacing_y = 1.0;
877 float64 spacing_z = 1.0;
878
879
880 if(n_coords.has_child("origin"))
881 {
882 const Node &n_origin = n_coords["origin"];
883
884 if(n_origin.has_child("x"))
885 {
886 origin_x = n_origin["x"].to_float64();
887 }
888
889 if(n_origin.has_child("y"))
890 {
891 origin_y = n_origin["y"].to_float64();
892 }
893
894 if(n_origin.has_child("z"))
895 {
896 origin_z = n_origin["z"].to_float64();
897 }
898 }
899
900 if(n_coords.has_path("spacing"))
901 {
902 const Node &n_spacing = n_coords["spacing"];
903
904 if(n_spacing.has_path("dx"))
905 {
906 spacing_x = n_spacing["dx"].to_float64();
907 }
908
909 if(n_spacing.has_path("dy"))
910 {
911 spacing_y = n_spacing["dy"].to_float64();
912 }
913
914 if(n_spacing.has_path("dz"))
915 {
916 spacing_z = n_spacing["dz"].to_float64();
917 }
918 }
919
920 // todo, should this be float64 -- or should we read float32 above?
921
922 vtkm::Vec<vtkm::Float32,3> origin(origin_x,
923 origin_y,
924 origin_z);
925
926 vtkm::Vec<vtkm::Float32,3> spacing(spacing_x,
927 spacing_y,
928 spacing_z);
929
930 vtkm::Id3 dims(dims_i,
931 dims_j,
932 dims_k);
933
934 // todo, use actually coordset and topo names?
935 result->AddCoordinateSystem( vtkm::cont::CoordinateSystem(coords_name.c_str(),
936 dims,
937 origin,
938 spacing));
939 vtkm::Id3 topo_origin = detail::topo_origin(n_topo);
940 if(is_2d)
941 {
942 vtkm::Id2 dims2(dims[0], dims[1]);
943 vtkm::cont::CellSetStructured<2> cell_set;
944 cell_set.SetPointDimensions(dims2);
945 vtkm::Id2 origin2(topo_origin[0], topo_origin[1]);
946 cell_set.SetGlobalPointIndexStart(origin2);
947 result->SetCellSet(cell_set);
948 }
949 else
950 {
951 vtkm::cont::CellSetStructured<3> cell_set;
952 cell_set.SetPointDimensions(dims);
953 cell_set.SetGlobalPointIndexStart(topo_origin);
954 result->SetCellSet(cell_set);
955 }
956
957 neles = (dims_i - 1) * (dims_j - 1);
958 if(dims_k > 1)
959 {
960 neles *= (dims_k - 1);
961 }
962
963 nverts = dims_i * dims_j;
964 if(dims_k > 1)
965 {
966 nverts *= dims_k;
967 }
968
969 return result;
970 }
971
972
973 //-----------------------------------------------------------------------------
974
975 vtkm::cont::DataSet *
RectilinearBlueprintToVTKmDataSet(const std::string & coords_name,const Node & n_coords,const std::string & topo_name,const Node & n_topo,int & neles,int & nverts,bool zero_copy)976 VTKHDataAdapter::RectilinearBlueprintToVTKmDataSet
977 (const std::string &coords_name, // input string with coordset name
978 const Node &n_coords, // input mesh bp coordset (assumed rectilinear)
979 const std::string &topo_name, // input string with topo name
980 const Node &n_topo, // input mesh bp topo
981 int &neles, // output, number of eles
982 int &nverts, // output, number of verts
983 bool zero_copy) // attempt to zero copy
984 {
985 vtkm::cont::DataSet *result = new vtkm::cont::DataSet();
986
987 int x_npts = n_coords["values/x"].dtype().number_of_elements();
988 int y_npts = n_coords["values/y"].dtype().number_of_elements();
989 int z_npts = 0;
990
991 int32 ndims = 2;
992
993 // todo assumes float64
994 const float64 *x_coords_ptr = n_coords["values/x"].as_float64_ptr();
995 const float64 *y_coords_ptr = n_coords["values/y"].as_float64_ptr();
996 const float64 *z_coords_ptr = NULL;
997
998 if(n_coords.has_path("values/z"))
999 {
1000 ndims = 3;
1001 z_npts = n_coords["values/z"].dtype().number_of_elements();
1002 z_coords_ptr = n_coords["values/z"].as_float64_ptr();
1003 }
1004
1005 vtkm::cont::ArrayHandle<vtkm::Float64> x_coords_handle;
1006 vtkm::cont::ArrayHandle<vtkm::Float64> y_coords_handle;
1007 vtkm::cont::ArrayHandle<vtkm::Float64> z_coords_handle;
1008
1009 if(zero_copy)
1010 {
1011 x_coords_handle = vtkm::cont::make_ArrayHandle(x_coords_ptr, x_npts, vtkm::CopyFlag::Off);
1012 y_coords_handle = vtkm::cont::make_ArrayHandle(y_coords_ptr, y_npts, vtkm::CopyFlag::Off);
1013 }
1014 else
1015 {
1016 x_coords_handle.Allocate(x_npts);
1017 y_coords_handle.Allocate(y_npts);
1018
1019 vtkm::Float64 *x = vtkh::GetVTKMPointer(x_coords_handle);
1020 memcpy(x, x_coords_ptr, sizeof(float64) * x_npts);
1021 vtkm::Float64 *y = vtkh::GetVTKMPointer(y_coords_handle);
1022 memcpy(y, y_coords_ptr, sizeof(float64) * y_npts);
1023 }
1024
1025 if(ndims == 3)
1026 {
1027 if(zero_copy)
1028 {
1029 z_coords_handle = vtkm::cont::make_ArrayHandle(z_coords_ptr, z_npts, vtkm::CopyFlag::Off);
1030 }
1031 else
1032 {
1033 z_coords_handle.Allocate(z_npts);
1034 vtkm::Float64 *z = vtkh::GetVTKMPointer(z_coords_handle);
1035 memcpy(z, z_coords_ptr, sizeof(float64) * z_npts);
1036 }
1037 }
1038 else
1039 {
1040 z_coords_handle.Allocate(1);
1041 z_coords_handle.WritePortal().Set(0, 0.0);
1042 }
1043
1044 static_assert(std::is_same<vtkm::FloatDefault, double>::value,
1045 "VTK-m needs to be configured with 'VTKm_USE_DOUBLE_PRECISION=ON'");
1046 vtkm::cont::ArrayHandleCartesianProduct<
1047 vtkm::cont::ArrayHandle<vtkm::FloatDefault>,
1048 vtkm::cont::ArrayHandle<vtkm::FloatDefault>,
1049 vtkm::cont::ArrayHandle<vtkm::FloatDefault> > coords;
1050
1051 coords = vtkm::cont::make_ArrayHandleCartesianProduct(x_coords_handle,
1052 y_coords_handle,
1053 z_coords_handle);
1054
1055 vtkm::cont::CoordinateSystem coordinate_system(coords_name.c_str(),
1056 coords);
1057 result->AddCoordinateSystem(coordinate_system);
1058
1059 vtkm::Id3 topo_origin = detail::topo_origin(n_topo);
1060
1061 if (ndims == 2)
1062 {
1063 vtkm::cont::CellSetStructured<2> cell_set;
1064 cell_set.SetPointDimensions(vtkm::make_Vec(x_npts,
1065 y_npts));
1066 vtkm::Id2 origin2(topo_origin[0], topo_origin[1]);
1067 cell_set.SetGlobalPointIndexStart(origin2);
1068 result->SetCellSet(cell_set);
1069 }
1070 else
1071 {
1072 vtkm::cont::CellSetStructured<3> cell_set;
1073 cell_set.SetPointDimensions(vtkm::make_Vec(x_npts,
1074 y_npts,
1075 z_npts));
1076 cell_set.SetGlobalPointIndexStart(topo_origin);
1077 result->SetCellSet(cell_set);
1078 }
1079
1080 nverts = x_npts * y_npts;
1081 neles = (x_npts - 1) * (y_npts - 1);
1082 if(ndims > 2)
1083 {
1084 nverts *= z_npts;
1085 neles *= (z_npts - 1);
1086 }
1087
1088 return result;
1089 }
1090
1091 //-----------------------------------------------------------------------------
1092
1093 vtkm::cont::DataSet *
StructuredBlueprintToVTKmDataSet(const std::string & coords_name,const Node & n_coords,const std::string & topo_name,const Node & n_topo,int & neles,int & nverts,bool zero_copy)1094 VTKHDataAdapter::StructuredBlueprintToVTKmDataSet
1095 (const std::string &coords_name, // input string with coordset name
1096 const Node &n_coords, // input mesh bp coordset (assumed rectilinear)
1097 const std::string &topo_name, // input string with topo name
1098 const Node &n_topo, // input mesh bp topo
1099 int &neles, // output, number of eles
1100 int &nverts, // output, number of verts
1101 bool zero_copy) // attempt to zero copy
1102 {
1103 vtkm::cont::DataSet *result = new vtkm::cont::DataSet();
1104
1105 nverts = n_coords["values/x"].dtype().number_of_elements();
1106
1107 int ndims = 0;
1108
1109 vtkm::cont::CoordinateSystem coords;
1110 if(n_coords["values/x"].dtype().is_float64())
1111 {
1112 coords = detail::GetExplicitCoordinateSystem<float64>(n_coords,
1113 coords_name,
1114 ndims,
1115 zero_copy);
1116 }
1117 else if(n_coords["values/x"].dtype().is_float32())
1118 {
1119 coords = detail::GetExplicitCoordinateSystem<float32>(n_coords,
1120 coords_name,
1121 ndims,
1122 zero_copy);
1123 }
1124 else
1125 {
1126 ASCENT_ERROR("Coordinate system must be floating point values");
1127 }
1128
1129 result->AddCoordinateSystem(coords);
1130
1131 int32 x_elems = n_topo["elements/dims/i"].as_int32();
1132 int32 y_elems = n_topo["elements/dims/j"].as_int32();
1133
1134 vtkm::Id3 topo_origin = detail::topo_origin(n_topo);
1135
1136 if (ndims == 2)
1137 {
1138 vtkm::cont::CellSetStructured<2> cell_set;
1139 cell_set.SetPointDimensions(vtkm::make_Vec(x_elems+1,
1140 y_elems+1));
1141 vtkm::Id2 origin2(topo_origin[0], topo_origin[1]);
1142 cell_set.SetGlobalPointIndexStart(origin2);
1143 result->SetCellSet(cell_set);
1144 neles = x_elems * y_elems;
1145 }
1146 else
1147 {
1148 int32 z_elems = n_topo["elements/dims/k"].as_int32();
1149 vtkm::cont::CellSetStructured<3> cell_set;
1150 cell_set.SetPointDimensions(vtkm::make_Vec(x_elems+1,
1151 y_elems+1,
1152 z_elems+1));
1153 cell_set.SetGlobalPointIndexStart(topo_origin);
1154 result->SetCellSet(cell_set);
1155 neles = x_elems * y_elems * z_elems;
1156
1157 }
1158 return result;
1159 }
1160
1161
1162
1163 //-----------------------------------------------------------------------------
1164
1165 vtkm::cont::DataSet *
UnstructuredBlueprintToVTKmDataSet(const std::string & coords_name,const Node & n_coords,const std::string & topo_name,const Node & n_topo,int & neles,int & nverts,bool zero_copy)1166 VTKHDataAdapter::UnstructuredBlueprintToVTKmDataSet
1167 (const std::string &coords_name, // input string with coordset name
1168 const Node &n_coords, // input mesh bp coordset (assumed unstructured)
1169 const std::string &topo_name, // input string with topo name
1170 const Node &n_topo, // input mesh bp topo
1171 int &neles, // output, number of eles
1172 int &nverts, // output, number of verts
1173 bool zero_copy) // attempt to zero copy
1174 {
1175 vtkm::cont::DataSet *result = new vtkm::cont::DataSet();
1176
1177 nverts = n_coords["values/x"].dtype().number_of_elements();
1178
1179 int32 ndims;
1180 vtkm::cont::CoordinateSystem coords;
1181 if(n_coords["values/x"].dtype().is_float64())
1182 {
1183 coords = detail::GetExplicitCoordinateSystem<float64>(n_coords,
1184 coords_name,
1185 ndims,
1186 zero_copy);
1187 }
1188 else if(n_coords["values/x"].dtype().is_float32())
1189 {
1190 coords = detail::GetExplicitCoordinateSystem<float32>(n_coords,
1191 coords_name,
1192 ndims,
1193 zero_copy);
1194 }
1195 else
1196 {
1197 ASCENT_ERROR("Coordinate system must be floating point values");
1198 }
1199
1200 result->AddCoordinateSystem(coords);
1201
1202 // shapes, number of indices, and connectivity.
1203 // Will have to do something different if this is a "zoo"
1204
1205 // TODO: there is a special data set type for single cell types
1206
1207 const Node &n_topo_eles = n_topo["elements"];
1208 std::string ele_shape = n_topo_eles["shape"].as_string();
1209
1210 // TODO: assumes int32, and contiguous
1211
1212 const Node &n_topo_conn = n_topo_eles["connectivity"];
1213
1214 vtkm::cont::ArrayHandle<vtkm::Id> connectivity;
1215
1216 int conn_size = n_topo_conn.dtype().number_of_elements();
1217
1218 if( sizeof(vtkm::Id) == 4)
1219 {
1220 if(n_topo_conn.is_compact() && n_topo_conn.dtype().is_int32())
1221 {
1222 const void *ele_idx_ptr = n_topo_conn.data_ptr();
1223 detail::CopyArray(connectivity, (const vtkm::Id*)ele_idx_ptr, conn_size,zero_copy);
1224 }
1225 else
1226 {
1227 // convert to int32
1228 connectivity.Allocate(conn_size);
1229 void *ptr = (void*) vtkh::GetVTKMPointer(connectivity);
1230 Node n_tmp;
1231 n_tmp.set_external(DataType::int32(conn_size),ptr);
1232 n_topo_conn.to_int32_array(n_tmp);
1233 }
1234 }
1235 else
1236 {
1237 if(n_topo_conn.is_compact() && n_topo_conn.dtype().is_int64())
1238 {
1239 const void *ele_idx_ptr = n_topo_conn.data_ptr();
1240 detail::CopyArray(connectivity, (const vtkm::Id*)ele_idx_ptr, conn_size, zero_copy);
1241 }
1242 else
1243 {
1244 // convert to int64
1245 connectivity.Allocate(conn_size);
1246 void *ptr = (void*) vtkh::GetVTKMPointer(connectivity);
1247 Node n_tmp;
1248 n_tmp.set_external(DataType::int64(conn_size),ptr);
1249 n_topo_conn.to_int64_array(n_tmp);
1250 }
1251 }
1252
1253 vtkm::UInt8 shape_id;
1254 vtkm::IdComponent indices_per;
1255 detail::VTKmCellShape(ele_shape, shape_id, indices_per);
1256 vtkm::cont::CellSetSingleType<> cellset;
1257 cellset.Fill(nverts, shape_id, indices_per, connectivity);
1258 neles = cellset.GetNumberOfCells();
1259 result->SetCellSet(cellset);
1260
1261 return result;
1262 }
1263
1264 //-----------------------------------------------------------------------------
1265
1266 void
AddField(const std::string & field_name,const Node & n_field,const std::string & topo_name,int neles,int nverts,vtkm::cont::DataSet * dset,bool zero_copy)1267 VTKHDataAdapter::AddField(const std::string &field_name,
1268 const Node &n_field,
1269 const std::string &topo_name,
1270 int neles,
1271 int nverts,
1272 vtkm::cont::DataSet *dset,
1273 bool zero_copy) // attempt to zero copy
1274 {
1275 // TODO: how do we deal with vector valued fields?, these will be mcarrays
1276
1277 string assoc_str = n_field["association"].as_string();
1278
1279 vtkm::cont::Field::Association vtkm_assoc = vtkm::cont::Field::Association::ANY;
1280 if(assoc_str == "vertex")
1281 {
1282 vtkm_assoc = vtkm::cont::Field::Association::POINTS;
1283 }
1284 else if(assoc_str == "element")
1285 {
1286 vtkm_assoc = vtkm::cont::Field::Association::CELL_SET;
1287 }
1288 else
1289 {
1290 ASCENT_INFO("VTKm conversion does not support field assoc "<<assoc_str<<". Skipping");
1291 return;
1292 }
1293 if(n_field["values"].number_of_children() > 1)
1294 {
1295 ASCENT_ERROR("Add field can only use zero or one component");
1296 }
1297
1298 bool is_values = n_field["values"].number_of_children() == 0;
1299 const Node &n_vals = is_values ? n_field["values"] : n_field["values"].child(0);
1300 int num_vals = n_vals.dtype().number_of_elements();
1301
1302 if(assoc_str == "vertex" && nverts != num_vals)
1303 {
1304 ASCENT_INFO("Field '"<<field_name<<"' (topology: '" << topo_name <<
1305 "') number of values "<<num_vals<<
1306 " does not match the number of points "<<nverts<<". Skipping");
1307 return;
1308 }
1309
1310 if(assoc_str == "element" && neles != num_vals)
1311 {
1312 if(field_name != "boundary_attribute")
1313 {
1314 ASCENT_INFO("Field '"<<field_name<<"' (topology: '" << topo_name <<
1315 "') number of values "<<num_vals<<
1316 " does not match the number of elements " << neles << ". Skipping");
1317 }
1318 return;
1319 }
1320
1321 try
1322 {
1323 bool supported_type = false;
1324
1325 if(n_vals.is_compact())
1326 {
1327 // we compile vtk-h with fp types
1328 if(n_vals.dtype().is_float32())
1329 {
1330 dset->AddField(detail::GetField<float32>(n_vals, field_name, assoc_str, topo_name, zero_copy));
1331 supported_type = true;
1332 }
1333 else if(n_vals.dtype().is_float64())
1334 {
1335 dset->AddField(detail::GetField<float64>(n_vals, field_name, assoc_str, topo_name, zero_copy));
1336 supported_type = true;
1337 }
1338 }
1339
1340 // vtk-m cant support zero copy for this layout or was not compiled to expose this datatype
1341 // use float64 by default
1342 if(!supported_type)
1343 {
1344 // convert to float64, we use this as a comprise to cover the widest range
1345 vtkm::cont::ArrayHandle<vtkm::Float64> vtkm_arr;
1346 vtkm_arr.Allocate(num_vals);
1347
1348 void *ptr = (void*) vtkh::GetVTKMPointer(vtkm_arr);
1349 Node n_tmp;
1350 n_tmp.set_external(DataType::float64(num_vals),ptr);
1351 n_vals.to_float64_array(n_tmp);
1352
1353 // add field to dataset
1354 if(assoc_str == "vertex")
1355 {
1356 dset->AddField(vtkm::cont::Field(field_name.c_str(),
1357 vtkm::cont::Field::Association::POINTS,
1358 vtkm_arr));
1359 }
1360 else if( assoc_str == "element")
1361 {
1362 dset->AddField(vtkm::cont::Field(field_name.c_str(),
1363 vtkm::cont::Field::Association::CELL_SET,
1364 vtkm_arr));
1365 }
1366 }
1367 }
1368 catch (vtkm::cont::Error error)
1369 {
1370 ASCENT_ERROR("VTKm exception:" << error.GetMessage());
1371 }
1372
1373 }
1374
1375 void
AddVectorField(const std::string & field_name,const Node & n_field,const std::string & topo_name,int neles,int nverts,vtkm::cont::DataSet * dset,const int dims,bool zero_copy)1376 VTKHDataAdapter::AddVectorField(const std::string &field_name,
1377 const Node &n_field,
1378 const std::string &topo_name,
1379 int neles,
1380 int nverts,
1381 vtkm::cont::DataSet *dset,
1382 const int dims,
1383 bool zero_copy) // attempt to zero copy
1384 {
1385 string assoc_str = n_field["association"].as_string();
1386
1387 vtkm::cont::Field::Association vtkm_assoc = vtkm::cont::Field::Association::ANY;
1388 if(assoc_str == "vertex")
1389 {
1390 vtkm_assoc = vtkm::cont::Field::Association::POINTS;
1391 }
1392 else if(assoc_str == "element")
1393 {
1394 vtkm_assoc = vtkm::cont::Field::Association::CELL_SET;
1395 }
1396 else
1397 {
1398 ASCENT_INFO("VTKm conversion does not support field assoc "<<assoc_str<<". Skipping");
1399 return;
1400 }
1401
1402
1403 const Node &n_vals = n_field["values"];
1404 int num_vals = n_vals.child(0).dtype().number_of_elements();
1405 int num_components = n_field["values"].number_of_children();
1406
1407 const conduit::Node &u = n_field["values"].child(0);
1408 bool interleaved = conduit::blueprint::mcarray::is_interleaved(n_vals);
1409 try
1410 {
1411 bool supported_type = false;
1412
1413 if(interleaved)
1414 {
1415 if(dims == 3)
1416 {
1417 // we compile vtk-h with fp types
1418 if(u.dtype().is_float32())
1419 {
1420
1421 using Vec3f32 = vtkm::Vec<vtkm::Float32,3>;
1422 const Vec3f32 *vec_ptr = reinterpret_cast<const Vec3f32*>(u.as_float32_ptr());
1423
1424 dset->AddField(detail::GetVectorField(vec_ptr,
1425 num_vals,
1426 field_name,
1427 assoc_str,
1428 topo_name,
1429 zero_copy));
1430 supported_type = true;
1431 }
1432 else if(u.dtype().is_float64())
1433 {
1434
1435 using Vec3f64 = vtkm::Vec<vtkm::Float64,3>;
1436 const Vec3f64 *vec_ptr = reinterpret_cast<const Vec3f64*>(u.as_float64_ptr());
1437
1438 dset->AddField(detail::GetVectorField(vec_ptr,
1439 num_vals,
1440 field_name,
1441 assoc_str,
1442 topo_name,
1443 zero_copy));
1444 supported_type = true;
1445 }
1446 }
1447 else if(dims == 2)
1448 {
1449 // we compile vtk-h with fp types
1450 if(u.dtype().is_float32())
1451 {
1452
1453 using Vec2f32 = vtkm::Vec<vtkm::Float32,2>;
1454 const Vec2f32 *vec_ptr = reinterpret_cast<const Vec2f32*>(u.as_float32_ptr());
1455
1456 dset->AddField(detail::GetVectorField(vec_ptr,
1457 num_vals,
1458 field_name,
1459 assoc_str,
1460 topo_name,
1461 zero_copy));
1462 supported_type = true;
1463 }
1464 else if(u.dtype().is_float64())
1465 {
1466
1467 using Vec2f64 = vtkm::Vec<vtkm::Float64,2>;
1468 const Vec2f64 *vec_ptr = reinterpret_cast<const Vec2f64*>(u.as_float64_ptr());
1469
1470 dset->AddField(detail::GetVectorField(vec_ptr,
1471 num_vals,
1472 field_name,
1473 assoc_str,
1474 topo_name,
1475 zero_copy));
1476 supported_type = true;
1477 }
1478 }
1479 else
1480 {
1481 ASCENT_ERROR("Vector unsupported dims " << dims);
1482 }
1483 }
1484 else
1485 {
1486 // we have a vector with 2/3 separate arrays
1487 // While vtkm supports ArrayHandleCompositeVectors for
1488 // coordinate systems, it does not support composites
1489 // for fields. Thus we have to copy the data.
1490 if(dims == 3)
1491 {
1492 const conduit::Node &v = n_field["values"].child(1);
1493 const conduit::Node &w = n_field["values"].child(2);
1494
1495 if(u.dtype().is_float32())
1496 {
1497 detail::ExtractVector<float32>(dset,
1498 u,
1499 v,
1500 w,
1501 num_vals,
1502 dims,
1503 field_name,
1504 assoc_str,
1505 topo_name,
1506 zero_copy);
1507 }
1508 else if(u.dtype().is_float64())
1509 {
1510 detail::ExtractVector<float64>(dset,
1511 u,
1512 v,
1513 w,
1514 num_vals,
1515 dims,
1516 field_name,
1517 assoc_str,
1518 topo_name,
1519 zero_copy);
1520 }
1521 }
1522 else if(dims == 2)
1523 {
1524 const conduit::Node &v = n_field["values"].child(1);
1525 conduit::Node fake_w;
1526 if(u.dtype().is_float32())
1527 {
1528 detail::ExtractVector<float32>(dset,
1529 u,
1530 v,
1531 fake_w,
1532 num_vals,
1533 dims,
1534 field_name,
1535 assoc_str,
1536 topo_name,
1537 zero_copy);
1538 }
1539 else if(u.dtype().is_float64())
1540 {
1541 detail::ExtractVector<float64>(dset,
1542 u,
1543 v,
1544 fake_w,
1545 num_vals,
1546 dims,
1547 field_name,
1548 assoc_str,
1549 topo_name,
1550 zero_copy);
1551 }
1552 }
1553 else
1554 {
1555 ASCENT_ERROR("Vector unsupported dims " << dims);
1556 }
1557 }
1558 }
1559 catch (vtkm::cont::Error error)
1560 {
1561 ASCENT_ERROR("VTKm exception:" << error.GetMessage());
1562 }
1563
1564 }
1565
1566 std::string
GetBlueprintCellName(vtkm::UInt8 shape_id)1567 GetBlueprintCellName(vtkm::UInt8 shape_id)
1568 {
1569 std::string name;
1570 if(shape_id == vtkm::CELL_SHAPE_TRIANGLE)
1571 {
1572 name = "tri";
1573 }
1574 else if(shape_id == vtkm::CELL_SHAPE_VERTEX)
1575 {
1576 name = "point";
1577 }
1578 else if(shape_id == vtkm::CELL_SHAPE_LINE)
1579 {
1580 name = "line";
1581 }
1582 else if(shape_id == vtkm::CELL_SHAPE_POLYGON)
1583 {
1584 ASCENT_ERROR("Polygon is not supported in blueprint");
1585 }
1586 else if(shape_id == vtkm::CELL_SHAPE_QUAD)
1587 {
1588 name = "quad";
1589 }
1590 else if(shape_id == vtkm::CELL_SHAPE_TETRA)
1591 {
1592 name = "tet";
1593 }
1594 else if(shape_id == vtkm::CELL_SHAPE_HEXAHEDRON)
1595 {
1596 name = "hex";
1597 }
1598 else if(shape_id == vtkm::CELL_SHAPE_WEDGE)
1599 {
1600 ASCENT_ERROR("Wedge is not supported in blueprint");
1601 }
1602 else if(shape_id == vtkm::CELL_SHAPE_PYRAMID)
1603 {
1604 ASCENT_ERROR("Pyramid is not supported in blueprint");
1605 }
1606 return name;
1607 }
1608
1609 bool
VTKmTopologyToBlueprint(conduit::Node & output,const vtkm::cont::DataSet & data_set,const std::string topo_name,bool zero_copy)1610 VTKHDataAdapter::VTKmTopologyToBlueprint(conduit::Node &output,
1611 const vtkm::cont::DataSet &data_set,
1612 const std::string topo_name,
1613 bool zero_copy)
1614 {
1615
1616 int topo_dims;
1617 bool is_structured = vtkh::VTKMDataSetInfo::IsStructured(data_set, topo_dims);
1618 bool is_uniform = vtkh::VTKMDataSetInfo::IsUniform(data_set);
1619 bool is_rectilinear = vtkh::VTKMDataSetInfo::IsRectilinear(data_set);
1620 vtkm::cont::CoordinateSystem coords = data_set.GetCoordinateSystem();
1621 const std::string coords_name = coords.GetName();
1622 // we cannot access an empty domain
1623 bool is_empty = false;
1624
1625 if(data_set.GetCoordinateSystem().GetData().GetNumberOfValues() == 0 ||
1626 data_set.GetCellSet().GetNumberOfCells() == 0)
1627 {
1628 is_empty = true;
1629 }
1630
1631 if(is_empty)
1632 {
1633 return is_empty;
1634 }
1635
1636 if(is_uniform)
1637 {
1638 auto points = coords.GetData().AsArrayHandle<vtkm::cont::ArrayHandleUniformPointCoordinates>();
1639 auto portal = points.ReadPortal();
1640
1641 auto origin = portal.GetOrigin();
1642 auto spacing = portal.GetSpacing();
1643 auto dims = portal.GetDimensions();
1644
1645 output["topologies/"+topo_name+"/coordset"] = coords_name;
1646 output["topologies/"+topo_name+"/type"] = "uniform";
1647
1648 output["coordsets/"+coords_name+"/type"] = "uniform";
1649 output["coordsets/"+coords_name+"/dims/i"] = (int) dims[0];
1650 output["coordsets/"+coords_name+"/dims/j"] = (int) dims[1];
1651 output["coordsets/"+coords_name+"/dims/k"] = (int) dims[2];
1652 output["coordsets/"+coords_name+"/origin/x"] = (double) origin[0];
1653 output["coordsets/"+coords_name+"/origin/y"] = (double) origin[1];
1654 output["coordsets/"+coords_name+"/origin/z"] = (double) origin[2];
1655 output["coordsets/"+coords_name+"/spacing/dx"] = (double) spacing[0];
1656 output["coordsets/"+coords_name+"/spacing/dy"] = (double) spacing[1];
1657 output["coordsets/"+coords_name+"/spacing/dz"] = (double) spacing[2];
1658 }
1659 else if(is_rectilinear)
1660 {
1661 typedef vtkm::cont::ArrayHandleCartesianProduct<vtkm::cont::ArrayHandle<vtkm::FloatDefault>,
1662 vtkm::cont::ArrayHandle<vtkm::FloatDefault>,
1663 vtkm::cont::ArrayHandle<vtkm::FloatDefault>> Cartesian;
1664
1665 const auto points = coords.GetData().AsArrayHandle<Cartesian>();
1666 auto portal = points.ReadPortal();
1667 auto x_portal = portal.GetFirstPortal();
1668 auto y_portal = portal.GetSecondPortal();
1669 auto z_portal = portal.GetThirdPortal();
1670
1671 // work around for conduit not accepting const pointers
1672 vtkm::FloatDefault *x_ptr = const_cast<vtkm::FloatDefault*>(x_portal.GetArray());
1673 vtkm::FloatDefault *y_ptr = const_cast<vtkm::FloatDefault*>(y_portal.GetArray());
1674 vtkm::FloatDefault *z_ptr = const_cast<vtkm::FloatDefault*>(z_portal.GetArray());
1675
1676 output["topologies/"+topo_name+"/coordset"] = coords_name;
1677 output["topologies/"+topo_name+"/type"] = "rectilinear";
1678
1679 output["coordsets/"+coords_name+"/type"] = "rectilinear";
1680 if(zero_copy)
1681 {
1682 output["coordsets/"+coords_name+"/values/x"].set_external(x_ptr, x_portal.GetNumberOfValues());
1683 output["coordsets/"+coords_name+"/values/y"].set_external(y_ptr, y_portal.GetNumberOfValues());
1684 output["coordsets/"+coords_name+"/values/z"].set_external(z_ptr, z_portal.GetNumberOfValues());
1685 }
1686 else
1687 {
1688 output["coordsets/"+coords_name+"/values/x"].set(x_ptr, x_portal.GetNumberOfValues());
1689 output["coordsets/"+coords_name+"/values/y"].set(y_ptr, y_portal.GetNumberOfValues());
1690 output["coordsets/"+coords_name+"/values/z"].set(z_ptr, z_portal.GetNumberOfValues());
1691 }
1692 }
1693 else
1694 {
1695 int point_dims[3];
1696 //
1697 // This still could be structured, but this will always
1698 // have an explicit coordinate system
1699 output["coordsets/"+coords_name+"/type"] = "explicit";
1700 using Coords32 = vtkm::cont::ArrayHandleSOA<vtkm::Vec<vtkm::Float32, 3>>;
1701 using Coords64 = vtkm::cont::ArrayHandleSOA<vtkm::Vec<vtkm::Float64, 3>>;
1702
1703 using CoordsVec32 = vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Float32,3>>;
1704 using CoordsVec64 = vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Float64,3>>;
1705
1706 vtkm::cont::UnknownArrayHandle coordsHandle(coords.GetData());
1707
1708 if(vtkm::cont::IsType<Coords32>(coordsHandle))
1709 {
1710 Coords32 points;
1711 coordsHandle.AsArrayHandle(points);
1712
1713 auto x_handle = points.GetArray(0);
1714 auto y_handle = points.GetArray(1);
1715 auto z_handle = points.GetArray(2);
1716
1717 point_dims[0] = x_handle.GetNumberOfValues();
1718 point_dims[1] = y_handle.GetNumberOfValues();
1719 point_dims[2] = z_handle.GetNumberOfValues();
1720
1721 if(zero_copy)
1722 {
1723 output["coordsets/"+coords_name+"/values/x"].
1724 set_external(vtkh::GetVTKMPointer(x_handle), point_dims[0]);
1725 output["coordsets/"+coords_name+"/values/y"].
1726 set_external(vtkh::GetVTKMPointer(y_handle), point_dims[1]);
1727 output["coordsets/"+coords_name+"/values/z"].
1728 set_external(vtkh::GetVTKMPointer(z_handle), point_dims[2]);
1729 }
1730 else
1731 {
1732 output["coordsets/"+coords_name+"/values/x"].
1733 set(vtkh::GetVTKMPointer(x_handle), point_dims[0]);
1734 output["coordsets/"+coords_name+"/values/y"].
1735 set(vtkh::GetVTKMPointer(y_handle), point_dims[1]);
1736 output["coordsets/"+coords_name+"/values/z"].
1737 set(vtkh::GetVTKMPointer(z_handle), point_dims[2]);
1738
1739 }
1740
1741 }
1742 else if(coordsHandle.IsType<CoordsVec32>())
1743 {
1744 CoordsVec32 points;
1745 coordsHandle.AsArrayHandle(points);
1746
1747 const int num_vals = points.GetNumberOfValues();
1748 vtkm::Float32 *points_ptr = (vtkm::Float32*)vtkh::GetVTKMPointer(points);
1749 const int byte_size = sizeof(vtkm::Float32);
1750
1751 if(zero_copy)
1752 {
1753 output["coordsets/"+coords_name+"/values/x"].set_external(points_ptr,
1754 num_vals,
1755 byte_size*0, // byte offset
1756 byte_size*3); // stride
1757 output["coordsets/"+coords_name+"/values/y"].set_external(points_ptr,
1758 num_vals,
1759 byte_size*1, // byte offset
1760 sizeof(vtkm::Float32)*3); // stride
1761 output["coordsets/"+coords_name+"/values/z"].set_external(points_ptr,
1762 num_vals,
1763 byte_size*2, // byte offset
1764 byte_size*3); // stride
1765 }
1766 else
1767 {
1768 output["coordsets/"+coords_name+"/values/x"].set(points_ptr,
1769 num_vals,
1770 byte_size*0, // byte offset
1771 byte_size*3); // stride
1772 output["coordsets/"+coords_name+"/values/y"].set(points_ptr,
1773 num_vals,
1774 byte_size*1, // byte offset
1775 sizeof(vtkm::Float32)*3); // stride
1776 output["coordsets/"+coords_name+"/values/z"].set(points_ptr,
1777 num_vals,
1778 byte_size*2, // byte offset
1779 byte_size*3); // stride
1780
1781 }
1782
1783 }
1784 else if(vtkm::cont::IsType<Coords64>(coordsHandle))
1785 {
1786 Coords64 points;
1787 coordsHandle.AsArrayHandle(points);
1788
1789 auto x_handle = points.GetArray(0);
1790 auto y_handle = points.GetArray(1);
1791 auto z_handle = points.GetArray(2);
1792
1793 point_dims[0] = x_handle.GetNumberOfValues();
1794 point_dims[1] = y_handle.GetNumberOfValues();
1795 point_dims[2] = z_handle.GetNumberOfValues();
1796 if(zero_copy)
1797 {
1798 output["coordsets/"+coords_name+"/values/x"].
1799 set_external(vtkh::GetVTKMPointer(x_handle), point_dims[0]);
1800 output["coordsets/"+coords_name+"/values/y"].
1801 set_external(vtkh::GetVTKMPointer(y_handle), point_dims[1]);
1802 output["coordsets/"+coords_name+"/values/z"].
1803 set_external(vtkh::GetVTKMPointer(z_handle), point_dims[2]);
1804 }
1805 else
1806 {
1807 output["coordsets/"+coords_name+"/values/x"].
1808 set(vtkh::GetVTKMPointer(x_handle), point_dims[0]);
1809 output["coordsets/"+coords_name+"/values/y"].
1810 set(vtkh::GetVTKMPointer(y_handle), point_dims[1]);
1811 output["coordsets/"+coords_name+"/values/z"].
1812 set(vtkh::GetVTKMPointer(z_handle), point_dims[2]);
1813
1814 }
1815 }
1816 else if(coordsHandle.IsType<CoordsVec64>())
1817 {
1818 CoordsVec64 points;
1819 coordsHandle.AsArrayHandle(points);
1820
1821 const int num_vals = points.GetNumberOfValues();
1822 vtkm::Float64 *points_ptr = (vtkm::Float64*)vtkh::GetVTKMPointer(points);
1823 const int byte_size = sizeof(vtkm::Float64);
1824
1825 if(zero_copy)
1826 {
1827 output["coordsets/"+coords_name+"/values/x"].set_external(points_ptr,
1828 num_vals,
1829 byte_size*0, // byte offset
1830 byte_size*3); // stride
1831 output["coordsets/"+coords_name+"/values/y"].set_external(points_ptr,
1832 num_vals,
1833 byte_size*1, // byte offset
1834 byte_size*3); // stride
1835 output["coordsets/"+coords_name+"/values/z"].set_external(points_ptr,
1836 num_vals,
1837 byte_size*2, // byte offset
1838 byte_size*3); // stride
1839 }
1840 else
1841 {
1842 output["coordsets/"+coords_name+"/values/x"].set(points_ptr,
1843 num_vals,
1844 byte_size*0, // byte offset
1845 byte_size*3); // stride
1846 output["coordsets/"+coords_name+"/values/y"].set(points_ptr,
1847 num_vals,
1848 byte_size*1, // byte offset
1849 byte_size*3); // stride
1850 output["coordsets/"+coords_name+"/values/z"].set(points_ptr,
1851 num_vals,
1852 byte_size*2, // byte offset
1853 byte_size*3); // stride
1854
1855 }
1856
1857 }
1858 else
1859 {
1860 // Ok vtkm has handed us something we don't know about, and its really
1861 // hard to ask vtkm to tell us what it is. Before we give up, we will
1862 // attempt to copy the data to a known type and copy that copy.
1863 // We can't avoid the double copy since conduit can't take ownership
1864 // and we can't seem to write to a zero copied array
1865
1866 vtkm::cont::ArrayHandle<vtkm::Vec<double,3>> coords_copy;
1867 vtkm::cont::ArrayCopy(coordsHandle, coords_copy);
1868 const int num_vals = coords_copy.GetNumberOfValues();
1869 vtkm::Float64 *points_ptr = (vtkm::Float64*)vtkh::GetVTKMPointer(coords_copy);
1870 const int byte_size = sizeof(vtkm::Float64);
1871
1872
1873 output["coordsets/"+coords_name+"/values/x"].set(points_ptr,
1874 num_vals,
1875 byte_size*0, // byte offset
1876 byte_size*3); // stride
1877 output["coordsets/"+coords_name+"/values/y"].set(points_ptr,
1878 num_vals,
1879 byte_size*1, // byte offset
1880 byte_size*3); // stride
1881 output["coordsets/"+coords_name+"/values/z"].set(points_ptr,
1882 num_vals,
1883 byte_size*2, // byte offset
1884 byte_size*3); // stride
1885 }
1886
1887 vtkm::UInt8 shape_id = 0;
1888 if(is_structured)
1889 {
1890 output["topologies/"+topo_name+"/coordset"] = coords_name;
1891 output["topologies/"+topo_name+"/type"] = "structured";
1892
1893 vtkm::cont::DynamicCellSet dyn_cells = data_set.GetCellSet();
1894 using Structured2D = vtkm::cont::CellSetStructured<2>;
1895 using Structured3D = vtkm::cont::CellSetStructured<3>;
1896 if(dyn_cells.IsSameType(Structured2D()))
1897 {
1898 Structured2D cells = dyn_cells.Cast<Structured2D>();
1899 vtkm::Id2 cell_dims = cells.GetCellDimensions();
1900 output["topologies/"+topo_name+"/elements/dims/i"] = (int) cell_dims[0];
1901 output["topologies/"+topo_name+"/elements/dims/j"] = (int) cell_dims[1];
1902 }
1903 else if(dyn_cells.IsSameType(Structured3D()))
1904 {
1905 Structured3D cells = dyn_cells.Cast<Structured3D>();
1906 vtkm::Id3 cell_dims = cells.GetCellDimensions();
1907 output["topologies/"+topo_name+"/elements/dims/i"] = (int) cell_dims[0];
1908 output["topologies/"+topo_name+"/elements/dims/j"] = (int) cell_dims[1];
1909 output["topologies/"+topo_name+"/elements/dims/k"] = (int) cell_dims[2];
1910 }
1911 else
1912 {
1913 ASCENT_ERROR("Unknown structured cell set");
1914 }
1915
1916 }
1917 else
1918 {
1919 output["topologies/"+topo_name+"/coordset"] = coords_name;
1920 output["topologies/"+topo_name+"/type"] = "unstructured";
1921 vtkm::cont::DynamicCellSet dyn_cells = data_set.GetCellSet();
1922
1923 using SingleType = vtkm::cont::CellSetSingleType<>;
1924 using MixedType = vtkm::cont::CellSetExplicit<>;
1925
1926 if(dyn_cells.IsSameType(SingleType()))
1927 {
1928 SingleType cells = dyn_cells.Cast<SingleType>();
1929 vtkm::UInt8 shape_id = cells.GetCellShape(0);
1930 std::string conduit_name = GetBlueprintCellName(shape_id);
1931 output["topologies/"+topo_name+"/elements/shape"] = conduit_name;
1932
1933 static_assert(sizeof(vtkm::Id) == sizeof(int), "blueprint expects connectivity to be ints");
1934 auto conn = cells.GetConnectivityArray(vtkm::TopologyElementTagCell(),
1935 vtkm::TopologyElementTagPoint());
1936
1937 if(zero_copy)
1938 {
1939 output["topologies/"+topo_name+"/elements/connectivity"].
1940 set_external(vtkh::GetVTKMPointer(conn), conn.GetNumberOfValues());
1941 }
1942 else
1943 {
1944 output["topologies/"+topo_name+"/elements/connectivity"].
1945 set(vtkh::GetVTKMPointer(conn), conn.GetNumberOfValues());
1946 }
1947 }
1948 else if(vtkh::VTKMDataSetInfo::IsSingleCellShape(dyn_cells, shape_id))
1949 {
1950 // If we are here, the we know that the cell set is explicit,
1951 // but only a single cell shape
1952 auto cells = dyn_cells.Cast<vtkm::cont::CellSetExplicit<>>();
1953 auto shapes = cells.GetShapesArray(vtkm::TopologyElementTagCell(),
1954 vtkm::TopologyElementTagPoint());
1955
1956 std::string conduit_name = GetBlueprintCellName(shape_id);
1957 output["topologies/"+topo_name+"/elements/shape"] = conduit_name;
1958
1959 static_assert(sizeof(vtkm::Id) == sizeof(int), "blueprint expects connectivity to be ints");
1960
1961 auto conn = cells.GetConnectivityArray(vtkm::TopologyElementTagCell(),
1962 vtkm::TopologyElementTagPoint());
1963
1964 if(zero_copy)
1965 {
1966 output["topologies/"+topo_name+"/elements/connectivity"].
1967 set_external(vtkh::GetVTKMPointer(conn), conn.GetNumberOfValues());
1968 }
1969 else
1970 {
1971 output["topologies/"+topo_name+"/elements/connectivity"].
1972 set(vtkh::GetVTKMPointer(conn), conn.GetNumberOfValues());
1973 }
1974
1975 }
1976 else
1977 {
1978 ASCENT_ERROR("Mixed explicit types not implemented");
1979 data_set.PrintSummary(std::cout);
1980 MixedType cells = dyn_cells.Cast<MixedType>();
1981 }
1982
1983 }
1984 }
1985 return is_empty;
1986 }
1987
1988 template<typename T, int N>
ConvertVecToNode(conduit::Node & output,std::string path,vtkm::cont::ArrayHandle<vtkm::Vec<T,N>> & handle,bool zero_copy)1989 void ConvertVecToNode(conduit::Node &output,
1990 std::string path,
1991 vtkm::cont::ArrayHandle<vtkm::Vec<T,N>> &handle,
1992 bool zero_copy)
1993 {
1994 static_assert(N > 1 && N < 4, "Vecs must be size 2 or 3");
1995 output[path + "/type"] = "vector";
1996 if(zero_copy)
1997 {
1998 output[path + "/values/u"].set_external((T*) vtkh::GetVTKMPointer(handle),
1999 handle.GetNumberOfValues(),
2000 sizeof(T)*0, // starting offset in bytes
2001 sizeof(T)*N); // stride in bytes
2002 output[path + "/values/v"].set_external((T*) vtkh::GetVTKMPointer(handle),
2003 handle.GetNumberOfValues(),
2004 sizeof(T)*1, // starting offset in bytes
2005 sizeof(T)*N); // stride in bytes
2006 }
2007 else
2008 {
2009 output[path + "/values/u"].set((T*) vtkh::GetVTKMPointer(handle),
2010 handle.GetNumberOfValues(),
2011 sizeof(T)*0, // starting offset in bytes
2012 sizeof(T)*N); // stride in bytes
2013 output[path + "/values/v"].set((T*) vtkh::GetVTKMPointer(handle),
2014 handle.GetNumberOfValues(),
2015 sizeof(T)*1, // starting offset in bytes
2016 sizeof(T)*N); // stride in bytes
2017 }
2018 if(N == 3)
2019 {
2020
2021 if(zero_copy)
2022 {
2023 output[path + "/values/w"].set_external((T*) vtkh::GetVTKMPointer(handle),
2024 handle.GetNumberOfValues(),
2025 sizeof(T)*2, // starting offset in bytes
2026 sizeof(T)*N); // stride in bytes
2027 }
2028 else
2029 {
2030 output[path + "/values/w"].set((T*) vtkh::GetVTKMPointer(handle),
2031 handle.GetNumberOfValues(),
2032 sizeof(T)*2, // starting offset in bytes
2033 sizeof(T)*N); // stride in bytes
2034 }
2035 }
2036 }
2037
2038 void
VTKmFieldToBlueprint(conduit::Node & output,const vtkm::cont::Field & field,const std::string topo_name,bool zero_copy)2039 VTKHDataAdapter::VTKmFieldToBlueprint(conduit::Node &output,
2040 const vtkm::cont::Field &field,
2041 const std::string topo_name,
2042 bool zero_copy)
2043 {
2044 std::string name = field.GetName();
2045 std::string path = "fields/" + name;
2046 bool assoc_points = vtkm::cont::Field::Association::POINTS == field.GetAssociation();
2047 bool assoc_cells = vtkm::cont::Field::Association::CELL_SET == field.GetAssociation();
2048 //bool assoc_mesh = vtkm::cont::Field::ASSOC_WHOLE_MESH == field.GetAssociation();
2049 if(!assoc_points && ! assoc_cells)
2050 {
2051 ASCENT_ERROR("Field must be associtated with cells or points\n");
2052 }
2053 std::string conduit_name;
2054
2055 if(assoc_points) conduit_name = "vertex";
2056 else conduit_name = "element";
2057
2058 output[path + "/association"] = conduit_name;
2059 output[path + "/topology"] = topo_name;
2060
2061 vtkm::cont::UnknownArrayHandle dyn_handle = field.GetData();
2062 //
2063 // this can be literally anything. Lets do some exhaustive casting
2064 //
2065 if(dyn_handle.IsType<vtkm::cont::ArrayHandle<vtkm::Float32>>())
2066 {
2067 using HandleType = vtkm::cont::ArrayHandle<vtkm::Float32>;
2068 HandleType handle;
2069 dyn_handle.AsArrayHandle(handle);
2070 if(zero_copy)
2071 {
2072 output[path + "/values"].
2073 set_external(vtkh::GetVTKMPointer(handle), handle.GetNumberOfValues());
2074 }
2075 else
2076 {
2077 output[path + "/values"]. set(vtkh::GetVTKMPointer(handle), handle.GetNumberOfValues());
2078 }
2079 }
2080 else if(dyn_handle.IsType<vtkm::cont::ArrayHandle<vtkm::Float64>>())
2081 {
2082 using HandleType = vtkm::cont::ArrayHandle<vtkm::Float64>;
2083 HandleType handle;
2084 dyn_handle.AsArrayHandle(handle);
2085 if(zero_copy)
2086 {
2087 output[path + "/values"].
2088 set_external(vtkh::GetVTKMPointer(handle), handle.GetNumberOfValues());
2089 }
2090 else
2091 {
2092 output[path + "/values"]. set(vtkh::GetVTKMPointer(handle), handle.GetNumberOfValues());
2093 }
2094 }
2095 else if(dyn_handle.IsType<vtkm::cont::ArrayHandle<vtkm::Int8>>())
2096 {
2097 using HandleType = vtkm::cont::ArrayHandle<vtkm::Int8>;
2098 HandleType handle;
2099 dyn_handle.AsArrayHandle(handle);
2100 if(zero_copy)
2101 {
2102 output[path + "/values"].
2103 set_external(vtkh::GetVTKMPointer(handle), handle.GetNumberOfValues());
2104 }
2105 else
2106 {
2107 output[path + "/values"]. set(vtkh::GetVTKMPointer(handle), handle.GetNumberOfValues());
2108 }
2109 }
2110 else if(dyn_handle.IsType<vtkm::cont::ArrayHandle<vtkm::Int32>>())
2111 {
2112 using HandleType = vtkm::cont::ArrayHandle<vtkm::Int32>;
2113 HandleType handle;
2114 dyn_handle.AsArrayHandle(handle);
2115 if(zero_copy)
2116 {
2117 output[path + "/values"].
2118 set_external(vtkh::GetVTKMPointer(handle), handle.GetNumberOfValues());
2119 }
2120 else
2121 {
2122 output[path + "/values"]. set(vtkh::GetVTKMPointer(handle), handle.GetNumberOfValues());
2123 }
2124 }
2125 else if(dyn_handle.IsType<vtkm::cont::ArrayHandle<vtkm::Int64>>())
2126 {
2127 using HandleType = vtkm::cont::ArrayHandle<vtkm::Int64>;
2128 HandleType handle;
2129 dyn_handle.AsArrayHandle(handle);
2130 ASCENT_ERROR("Conduit int64 and vtkm::Int64 are different. Cannot convert vtkm::Int64\n");
2131 //output[path + "/values"].set(vtkh::GetVTKMPointer(handle), handle.GetNumberOfValues());
2132 }
2133 else if(dyn_handle.IsType<vtkm::cont::ArrayHandle<vtkm::UInt32>>())
2134 {
2135 using HandleType = vtkm::cont::ArrayHandle<vtkm::UInt32>;
2136 HandleType handle;
2137 dyn_handle.AsArrayHandle(handle);
2138 if(zero_copy)
2139 {
2140 output[path + "/values"].
2141 set_external(vtkh::GetVTKMPointer(handle), handle.GetNumberOfValues());
2142 }
2143 else
2144 {
2145 output[path + "/values"]. set(vtkh::GetVTKMPointer(handle), handle.GetNumberOfValues());
2146 }
2147 }
2148 else if(dyn_handle.IsType<vtkm::cont::ArrayHandle<vtkm::UInt8>>())
2149 {
2150 using HandleType = vtkm::cont::ArrayHandle<vtkm::UInt8>;
2151 HandleType handle;
2152 dyn_handle.AsArrayHandle(handle);
2153 if(zero_copy)
2154 {
2155 output[path + "/values"].
2156 set_external(vtkh::GetVTKMPointer(handle), handle.GetNumberOfValues());
2157 }
2158 else
2159 {
2160 output[path + "/values"]. set(vtkh::GetVTKMPointer(handle), handle.GetNumberOfValues());
2161 }
2162 }
2163 else if(dyn_handle.IsType<vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Float32,3>>>())
2164 {
2165 using HandleType = vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Float32,3>>;
2166 HandleType handle;
2167 dyn_handle.AsArrayHandle(handle);
2168 ConvertVecToNode(output, path, handle, zero_copy);
2169 }
2170 else if(dyn_handle.IsType<vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Float64,3>>>())
2171 {
2172 using HandleType = vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Float64,3>>;
2173 HandleType handle;
2174 dyn_handle.AsArrayHandle(handle);
2175 ConvertVecToNode(output, path, handle, zero_copy);
2176 }
2177 else if(dyn_handle.IsType<vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Int32,3>>>())
2178 {
2179 using HandleType = vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Int32,3>>;
2180 HandleType handle;
2181 dyn_handle.AsArrayHandle(handle);
2182 ConvertVecToNode(output, path, handle, zero_copy);
2183 }
2184 else if(dyn_handle.IsType<vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Float32,2>>>())
2185 {
2186 using HandleType = vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Float32,2>>;
2187 HandleType handle;
2188 dyn_handle.AsArrayHandle(handle);
2189 ConvertVecToNode(output, path, handle, zero_copy);
2190 }
2191 else if(dyn_handle.IsType<vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Float64,2>>>())
2192 {
2193 using HandleType = vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Float64,2>>;
2194 HandleType handle;
2195 dyn_handle.AsArrayHandle(handle);
2196 ConvertVecToNode(output, path, handle, zero_copy);
2197 }
2198 else if(dyn_handle.IsType<vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Int32,2>>>())
2199 {
2200 using HandleType = vtkm::cont::ArrayHandle<vtkm::Vec<vtkm::Int32,2>>;
2201 HandleType handle;
2202 dyn_handle.AsArrayHandle(handle);
2203 ConvertVecToNode(output, path, handle, zero_copy);
2204 }
2205 else
2206 {
2207 std::stringstream msg;
2208 msg<<"Field type unsupported for conversion to blueprint.\n";
2209 field.PrintSummary(msg);
2210 msg<<" Skipping.";
2211 ASCENT_INFO(msg.str());
2212 }
2213 }
2214
VTKHCollectionToBlueprintDataSet(VTKHCollection * collection,conduit::Node & node,bool zero_copy)2215 void VTKHDataAdapter::VTKHCollectionToBlueprintDataSet(VTKHCollection *collection,
2216 conduit::Node &node,
2217 bool zero_copy)
2218 {
2219 node.reset();
2220
2221 bool success = true;
2222 // we have to re-merge the domains so all domains with the same
2223 // domain id end up in a single domain
2224 std::map<int, std::map<std::string,vtkm::cont::DataSet>> domain_map;
2225 domain_map = collection->by_domain_id();
2226
2227 try
2228 {
2229 for(auto domain_it : domain_map)
2230 {
2231 const int domain_id = domain_it.first;
2232
2233 conduit::Node &dom = node.append();
2234 dom["state/domain_id"] = (int) domain_id;
2235
2236 for(auto topo_it : domain_it.second)
2237 {
2238 const std::string topo_name = topo_it.first;
2239 vtkm::cont::DataSet &dataset = topo_it.second;
2240 VTKHDataAdapter::VTKmToBlueprintDataSet(&dataset, dom, topo_name, zero_copy);
2241 }
2242 }
2243 }
2244 catch(...)
2245 {
2246 success = false;
2247 }
2248
2249 success = global_agreement(success);
2250 if(!success)
2251 {
2252 ASCENT_ERROR("Failed to convert vtkm data set to blueprint");
2253 }
2254 }
2255
2256 void
VTKHToBlueprintDataSet(vtkh::DataSet * dset,conduit::Node & node,bool zero_copy)2257 VTKHDataAdapter::VTKHToBlueprintDataSet(vtkh::DataSet *dset,
2258 conduit::Node &node,
2259 bool zero_copy)
2260 {
2261 node.reset();
2262 bool success = true;
2263 try
2264 {
2265 const int num_doms = dset->GetNumberOfDomains();
2266 for(int i = 0; i < num_doms; ++i)
2267 {
2268 conduit::Node &dom = node.append();
2269 vtkm::cont::DataSet vtkm_dom;
2270 vtkm::Id domain_id;
2271 int cycle = dset->GetCycle();
2272 dset->GetDomain(i, vtkm_dom, domain_id);
2273 VTKHDataAdapter::VTKmToBlueprintDataSet(&vtkm_dom,dom, "topo", zero_copy);
2274 dom["state/domain_id"] = (int) domain_id;
2275 dom["state/cycle"] = cycle;
2276 }
2277 }
2278 catch(...)
2279 {
2280 success = false;
2281 }
2282
2283 success = global_agreement(success);
2284 if(!success)
2285 {
2286 ASCENT_ERROR("Failed to convert vtkm data set to blueprint");
2287 }
2288 }
2289
2290 void
VTKmToBlueprintDataSet(const vtkm::cont::DataSet * dset,conduit::Node & node,const std::string topo_name,bool zero_copy)2291 VTKHDataAdapter::VTKmToBlueprintDataSet(const vtkm::cont::DataSet *dset,
2292 conduit::Node &node,
2293 const std::string topo_name,
2294 bool zero_copy)
2295 {
2296 //
2297 // with vtkm, we have no idea what the type is of anything inside
2298 // dataset, so we have to ask all fields, cell sets anc coordinate systems.
2299 //
2300
2301 bool is_empty = VTKmTopologyToBlueprint(node, *dset, topo_name, zero_copy);
2302
2303 if(!is_empty)
2304 {
2305 const vtkm::Id num_fields = dset->GetNumberOfFields();
2306 for(vtkm::Id i = 0; i < num_fields; ++i)
2307 {
2308 vtkm::cont::Field field = dset->GetField(i);
2309 VTKmFieldToBlueprint(node, field, topo_name, zero_copy);
2310 }
2311 }
2312 }
2313
2314
2315 };
2316 //-----------------------------------------------------------------------------
2317 // -- end ascent:: --
2318 //-----------------------------------------------------------------------------
2319