From 54b9f7037966ba9cb257bd1f66665c7b5af461c3 Mon Sep 17 00:00:00 2001 From: Nicole Marsaglia Date: Thu, 30 Nov 2023 00:31:55 +0000 Subject: [PATCH] make all fields optional --- .../ascent_runtime_vtkh_filters.cpp | 59 ++++++++----------- src/libs/vtkh/filters/WarpXStreamline.cpp | 54 +++++++++-------- src/libs/vtkh/filters/WarpXStreamline.hpp | 8 +++ src/tests/vtkh/t_vtk-h_warpx_streamline.cpp | 11 +--- 4 files changed, 63 insertions(+), 69 deletions(-) diff --git a/src/libs/ascent/runtimes/flow_filters/ascent_runtime_vtkh_filters.cpp b/src/libs/ascent/runtimes/flow_filters/ascent_runtime_vtkh_filters.cpp index e3bd0c2b5..02a143807 100644 --- a/src/libs/ascent/runtimes/flow_filters/ascent_runtime_vtkh_filters.cpp +++ b/src/libs/ascent/runtimes/flow_filters/ascent_runtime_vtkh_filters.cpp @@ -3965,6 +3965,10 @@ VTKHWarpXStreamline::verify_params(const conduit::Node ¶ms, std::vector valid_paths; valid_paths.push_back("b_field"); valid_paths.push_back("e_field"); + valid_paths.push_back("charge_field"); + valid_paths.push_back("mass_field"); + valid_paths.push_back("momentum_field"); + valid_paths.push_back("weighting_field"); valid_paths.push_back("num_steps"); valid_paths.push_back("step_size"); @@ -4000,19 +4004,23 @@ VTKHWarpXStreamline::execute() std::string b_field = "B"; std::string e_field = "E"; + std::string charge_field = "Charge"; + std::string mass_field = "Mass"; + std::string momentum_field = "Momentum"; + std::string weighting_field = "Weighting"; if(params().has_path("b_field")) b_field = params()["b_field"].as_string(); if(params().has_path("e_field")) e_field = params()["e_field"].as_string(); + if(params().has_path("charge_field")) + charge_field = params()["charge_field"].as_string(); + if(params().has_path("mass_field")) + mass_field = params()["mass_field"].as_string(); + if(params().has_path("momentum_field")) + momentum_field = params()["momentum_field"].as_string(); + if(params().has_path("weighting_field")) + weighting_field = params()["weighting_field"].as_string(); - if(!collection->has_field(e_field)) - { - bool throw_error = false; - detail::field_error(e_field, this->name(), collection, throw_error); - // this creates a data object with an invalid soource - set_output(new DataObject()); - return; - } if(!collection->has_field(b_field)) { bool throw_error = false; @@ -4021,40 +4029,15 @@ VTKHWarpXStreamline::execute() set_output(new DataObject()); return; } - - if(!collection->has_field("Momentum")) - { - bool throw_error = false; - detail::field_error("Momentum", this->name(), collection, throw_error); - // this creates a data object with an invalid soource - set_output(new DataObject()); - return; - } - if(!collection->has_field("Mass")) - { - bool throw_error = false; - detail::field_error("Mass", this->name(), collection, throw_error); - // this creates a data object with an invalid soource - set_output(new DataObject()); - return; - } - if(!collection->has_field("Charge")) - { - bool throw_error = false; - detail::field_error("Charge", this->name(), collection, throw_error); - // this creates a data object with an invalid soource - set_output(new DataObject()); - return; - } - if(!collection->has_field("Weighting")) + if(!collection->has_field(e_field)) { bool throw_error = false; - detail::field_error("Weighting", this->name(), collection, throw_error); + detail::field_error(e_field, this->name(), collection, throw_error); // this creates a data object with an invalid soource set_output(new DataObject()); return; } - + std::string topo_name = collection->field_topology(b_field); vtkh::DataSet &data = collection->dataset_by_topology(topo_name); @@ -4069,6 +4052,10 @@ VTKHWarpXStreamline::execute() sl.SetNumberOfSteps(numSteps); sl.SetBField(b_field); sl.SetEField(e_field); + sl.SetChargeField(charge_field); + sl.SetMassField(mass_field); + sl.SetMomentumField(momentum_field); + sl.SetWeightingField(weighting_field); sl.SetInput(&data); sl.Update(); output = sl.GetOutput(); diff --git a/src/libs/vtkh/filters/WarpXStreamline.cpp b/src/libs/vtkh/filters/WarpXStreamline.cpp index cf081d5c6..acfe7033c 100644 --- a/src/libs/vtkh/filters/WarpXStreamline.cpp +++ b/src/libs/vtkh/filters/WarpXStreamline.cpp @@ -50,7 +50,12 @@ void GenerateChargedParticles(const vtkm::cont::ArrayHandle& pos, } //end detail WarpXStreamline::WarpXStreamline() -: m_e_field_name("E"), m_b_field_name("B") +: m_e_field_name("E"), + m_b_field_name("B"), + m_charge_field_name("Charge"), + m_mass_field_name("Mass"), + m_momentum_field_name("Momentum"), + m_weighting_field_name("Weighting") { } @@ -63,8 +68,13 @@ WarpXStreamline::~WarpXStreamline() void WarpXStreamline::PreExecute() { Filter::PreExecute(); - Filter::CheckForRequiredField(m_e_field_name); Filter::CheckForRequiredField(m_b_field_name); + Filter::CheckForRequiredField(m_e_field_name); + Filter::CheckForRequiredField(m_charge_field_name); + Filter::CheckForRequiredField(m_mass_field_name); + Filter::CheckForRequiredField(m_momentum_field_name); + Filter::CheckForRequiredField(m_weighting_field_name); + } void WarpXStreamline::PostExecute() @@ -100,7 +110,7 @@ void WarpXStreamline::DoExecute() vtkm::cont::ArrayHandle seeds; //Create charged particles for all domains with the particle spec fields. //TODO: user specified momentum,mass,charge,weighting? - if (this->m_input->FieldExists("Momentum")) + if (this->m_input->FieldExists(m_momentum_field_name)) { const int num_domains = this->m_input->GetNumberOfDomains(); int id_offset = 0; @@ -109,15 +119,15 @@ void WarpXStreamline::DoExecute() vtkm::Id domain_id; vtkm::cont::DataSet dom; this->m_input->GetDomain(i, dom, domain_id); - if(dom.HasField("Momentum")) + if(dom.HasField(m_momentum_field_name)) { vtkm::cont::ArrayHandle pos, mom; vtkm::cont::ArrayHandle mass, charge, w; dom.GetCoordinateSystem().GetData().AsArrayHandle(pos); - dom.GetField("Momentum").GetData().AsArrayHandle(mom); - dom.GetField("Mass").GetData().AsArrayHandle(mass); - dom.GetField("Charge").GetData().AsArrayHandle(charge); - dom.GetField("Weighting").GetData().AsArrayHandle(w); + dom.GetField(m_momentum_field_name).GetData().AsArrayHandle(mom); + dom.GetField(m_mass_field_name).GetData().AsArrayHandle(mass); + dom.GetField(m_charge_field_name).GetData().AsArrayHandle(charge); + dom.GetField(m_weighting_field_name).GetData().AsArrayHandle(w); detail::GenerateChargedParticles(pos,mom,mass,charge,w,seeds, id_offset); //Actual: local unique ids //Question: do we global unique ids? @@ -125,25 +135,20 @@ void WarpXStreamline::DoExecute() } if(dom.HasField(m_b_field_name)) { - using vectorField_d = vtkm::cont::ArrayHandle>; - using vectorField_f = vtkm::cont::ArrayHandle>; - std::cerr << "HERE 1" << std::endl; auto field = dom.GetField(m_b_field_name).GetData(); - std::cerr << "HERE 2" << std::endl; - if(field.IsType()) - std::cerr << "Vector field is DOUBLE" << std::endl; - if(field.IsType()) - std::cerr << "Vector field is FLOAT" << std::endl; - if(field.IsType() && !field.IsType()) - { - inputs.AppendPartition(dom); - } + inputs.AppendPartition(dom); } } } + else + { + throw Error("Domain is missing one or more neccessary fields to create a charged particle: \ + Charge, Mass, Momentum, and/or Weighting."); + } bool validField = (inputs.GetNumberOfPartitions() > 0); - +//Don't really need this check +//since we got rid of the other check #ifdef VTKH_PARALLEL int localNum = static_cast(inputs.GetNumberOfPartitions()); int globalNum = 0; @@ -158,7 +163,7 @@ void WarpXStreamline::DoExecute() if (!validField) { - throw Error("Vector field type does not match > or >"); + throw Error("Vector field type does not match a supportable type."); } //Everything is valid. Call the VTKm filter. @@ -171,10 +176,11 @@ void WarpXStreamline::DoExecute() warpxStreamlineFilter.SetSeeds(seeds); warpxStreamlineFilter.SetNumberOfSteps(m_num_steps); auto out = warpxStreamlineFilter.Execute(inputs); - + //out.PrintSummary(std::cerr); + int num_domains = m_output->GetNumberOfDomains(); for (vtkm::Id i = 0; i < out.GetNumberOfPartitions(); i++) { - this->m_output->AddDomain(out.GetPartition(i), i); + this->m_output->AddDomain(out.GetPartition(i), num_domains + i); } #endif } diff --git a/src/libs/vtkh/filters/WarpXStreamline.hpp b/src/libs/vtkh/filters/WarpXStreamline.hpp index 444670436..eab180ef3 100644 --- a/src/libs/vtkh/filters/WarpXStreamline.hpp +++ b/src/libs/vtkh/filters/WarpXStreamline.hpp @@ -19,6 +19,10 @@ class VTKH_API WarpXStreamline : public Filter std::string GetName() const override { return "vtkh::WarpXStreamline";} void SetBField(const std::string &field_name) { m_b_field_name = field_name; } void SetEField(const std::string &field_name) { m_e_field_name = field_name; } + void SetChargeField(const std::string &field_name) { m_charge_field_name = field_name; } + void SetMassField(const std::string &field_name) { m_mass_field_name = field_name; } + void SetMomentumField(const std::string &field_name) { m_momentum_field_name = field_name; } + void SetWeightingField(const std::string &field_name) { m_weighting_field_name = field_name; } void SetStepSize(const double &step_size) { m_step_size = step_size; } void SetNumberOfSteps(int numSteps) { m_num_steps = numSteps; } @@ -29,6 +33,10 @@ class VTKH_API WarpXStreamline : public Filter std::string m_b_field_name; std::string m_e_field_name; + std::string m_charge_field_name; + std::string m_mass_field_name; + std::string m_momentum_field_name; + std::string m_weighting_field_name; double m_step_size; int m_num_steps; }; diff --git a/src/tests/vtkh/t_vtk-h_warpx_streamline.cpp b/src/tests/vtkh/t_vtk-h_warpx_streamline.cpp index 3d5c0d9fe..27a9a5cf0 100644 --- a/src/tests/vtkh/t_vtk-h_warpx_streamline.cpp +++ b/src/tests/vtkh/t_vtk-h_warpx_streamline.cpp @@ -165,23 +165,16 @@ TEST(vtkh_particle_advection, vtkh_serial_particle_advection) vtkm::cont::CoordinateSystem coords = fieldsData.GetCoordinateSystem(); auto w_bounds = coords.GetBounds(); - std::cout << "Bounds : " << w_bounds << std::endl; using Structured3DType = vtkm::cont::CellSetStructured<3>; - std::cerr << "HERE 1" << std::endl; Structured3DType castedCells; - std::cerr << "HERE 2" << std::endl; cells.AsCellSet(castedCells); - std::cerr << "HERE 3" << std::endl; auto dims = castedCells.GetSchedulingRange(vtkm::TopologyElementTagPoint()); - std::cerr << "HERE 4" << std::endl; vtkm::Vec3f spacing = { static_cast(w_bounds.X.Length()) / (dims[0] - 1), static_cast(w_bounds.Y.Length()) / (dims[1] - 1), static_cast(w_bounds.Z.Length()) / (dims[2] - 1) }; - std::cerr << "HERE 5" << std::endl; constexpr static vtkm::FloatDefault SPEED_OF_LIGHT = static_cast(2.99792458e8); spacing = spacing * spacing; - std::cerr << "HERE 6" << std::endl; vtkm::FloatDefault length = static_cast( 1.0 / (SPEED_OF_LIGHT * vtkm::Sqrt(1. / spacing[0] + 1. / spacing[1] + 1. / spacing[2]))); @@ -190,9 +183,9 @@ TEST(vtkh_particle_advection, vtkh_serial_particle_advection) vtkh::DataSet *outWSL=NULL; - warpx_data_set.PrintSummary(std::cerr); + //warpx_data_set.PrintSummary(std::cerr); outWSL = RunWFilter(warpx_data_set, maxAdvSteps, length); - outWSL->PrintSummary(std::cout); + //outWSL->PrintSummary(std::cerr); checkValidity(outWSL, maxAdvSteps+1, true); writeDataSet(outWSL, "warpx_streamline", rank);