From a6d564d373aa354e2440ccb2b113ac24f54cd3b7 Mon Sep 17 00:00:00 2001
From: Silvio Rizzi <srizzi@anl.gov>
Date: Tue, 3 Aug 2021 19:23:44 +0000
Subject: [PATCH 01/11] Injecting old lammps files in this branch

---
 CMake/options.cmake                   |   2 +
 miniapps/CMakeLists.txt               |   7 +
 miniapps/lammps/CMakeLists.txt        |  31 ++
 miniapps/lammps/lammpsBridge.cxx      |  55 ++++
 miniapps/lammps/lammpsBridge.h        |  15 +
 miniapps/lammps/lammpsDataAdaptor.cxx | 416 ++++++++++++++++++++++++++
 miniapps/lammps/lammpsDataAdaptor.h   |  78 +++++
 miniapps/lammps/lammpsDriver.cxx      | 218 ++++++++++++++
 8 files changed, 822 insertions(+)
 create mode 100644 miniapps/lammps/CMakeLists.txt
 create mode 100644 miniapps/lammps/lammpsBridge.cxx
 create mode 100644 miniapps/lammps/lammpsBridge.h
 create mode 100644 miniapps/lammps/lammpsDataAdaptor.cxx
 create mode 100644 miniapps/lammps/lammpsDataAdaptor.h
 create mode 100644 miniapps/lammps/lammpsDriver.cxx

diff --git a/CMake/options.cmake b/CMake/options.cmake
index 79cd890b0..6c93d4abd 100644
--- a/CMake/options.cmake
+++ b/CMake/options.cmake
@@ -79,6 +79,7 @@ option(ENABLE_MANDELBROT "Enable Mandelbrot AMR miniapp" ON)
 option(ENABLE_VORTEX "Enable Vortex miniapp (experimental)" OFF)
 option(ENABLE_CONDUITTEST "Enable Conduit miniapp (experimental)" OFF)
 option(ENABLE_KRIPKE "Enable Kripke miniapp (experimental)" OFF)
+option(ENABLE_LAMMPS "Enable LAMMPS driver" OFF)
 option(SENSEI_USE_EXTERNAL_pugixml "Use external pugixml library" OFF)
 
 message(STATUS "ENABLE_SENSEI=${ENABLE_SENSEI}")
@@ -104,6 +105,7 @@ message(STATUS "ENABLE_OPTS=${ENABLE_OPTS}")
 message(STATUS "ENABLE_OSCILLATORS=${ENABLE_OSCILLATORS}")
 message(STATUS "ENABLE_CONDUITTEST=${ENABLE_CONDUITTEST}")
 message(STATUS "ENABLE_KRIPKE=${ENABLE_KRIPKE}")
+message(STATUS "ENABLE_LAMMPS=${ENABLE_LAMMPS}")
 message(STATUS "SENSEI_USE_EXTERNAL_pugixml=${SENSEI_USE_EXTERNAL_pugixml}")
 
 if (ENABLE_ADIOS1 AND ENABLE_ADIOS2)
diff --git a/miniapps/CMakeLists.txt b/miniapps/CMakeLists.txt
index 80e8d14df..ceefee8d7 100644
--- a/miniapps/CMakeLists.txt
+++ b/miniapps/CMakeLists.txt
@@ -33,3 +33,10 @@ else()
   message(STATUS "Disabled: Vortex miniapp.")
 endif()
 
+if(ENABLE_LAMMPS)
+  message(STATUS "Enabled: LAMMPS driver.")
+  add_subdirectory(lammps)
+else()
+  message(STATUS "Disabled: LAMMPS driver.")
+endif()
+
diff --git a/miniapps/lammps/CMakeLists.txt b/miniapps/lammps/CMakeLists.txt
new file mode 100644
index 000000000..10094321d
--- /dev/null
+++ b/miniapps/lammps/CMakeLists.txt
@@ -0,0 +1,31 @@
+if (NOT LAMMPS_DIR)
+  message(FATAL_ERROR "To build LAMMPS driver you must set -DLAMMPS_DIR=<root of lammps dir>")
+else()
+  message(STATUS "LAMMPS_DIR: ${LAMMPS_DIR}")
+endif()
+
+include_directories(${LAMMPS_DIR}/src)
+
+set(sources lammpsDriver.cxx)
+
+set(libs 
+    ${LAMMPS_DIR}/src/liblammps_mpi.a
+)
+
+if(ENABLE_SENSEI)
+  list(APPEND sources
+    lammpsDataAdaptor.h
+    lammpsDataAdaptor.cxx
+    lammpsBridge.h
+    lammpsBridge.cxx
+  )
+  list(APPEND libs sensei)
+endif()
+
+add_executable(lammpsDriver ${sources})
+
+target_link_libraries(lammpsDriver ${libs})
+
+install(TARGETS lammpsDriver
+  RUNTIME DESTINATION bin
+)
diff --git a/miniapps/lammps/lammpsBridge.cxx b/miniapps/lammps/lammpsBridge.cxx
new file mode 100644
index 000000000..15bf094e9
--- /dev/null
+++ b/miniapps/lammps/lammpsBridge.cxx
@@ -0,0 +1,55 @@
+#include "lammpsBridge.h"
+#include "lammpsDataAdaptor.h"
+#include <vtkSmartPointer.h>
+#include <ConfigurableAnalysis.h>
+
+#include <vtkDataObject.h>
+
+namespace lammpsBridge
+{
+  static vtkSmartPointer<senseiLammps::lammpsDataAdaptor>  DataAdaptor;
+  static vtkSmartPointer<sensei::ConfigurableAnalysis>     AnalysisAdaptor;
+
+void Initialize(MPI_Comm world, const std::string& config_file)
+{ 
+  DataAdaptor = vtkSmartPointer<senseiLammps::lammpsDataAdaptor>::New();
+  DataAdaptor->Initialize();
+  DataAdaptor->SetCommunicator(world);
+  DataAdaptor->SetDataTimeStep(-1);
+
+  AnalysisAdaptor = vtkSmartPointer<sensei::ConfigurableAnalysis>::New();
+  AnalysisAdaptor->Initialize(config_file);
+}
+
+void SetData(long ntimestep, int nlocal, int *id, 
+             int nghost, int *type, double **x, 
+             double xsublo, double xsubhi, 
+             double ysublo, double ysubhi, 
+             double zsublo, double zsubhi )
+{
+  DataAdaptor->AddLAMMPSData( ntimestep, nlocal, id, nghost, type, x, \
+                              xsublo, xsubhi, ysublo, ysubhi, zsublo, zsubhi);
+  DataAdaptor->SetDataTimeStep(ntimestep);
+  DataAdaptor->SetDataTime(ntimestep);
+}
+
+void Analyze()
+{
+  AnalysisAdaptor->Execute(DataAdaptor.GetPointer());
+
+  //vtkDataObject *mesh;
+  //DataAdaptor->GetMesh("atoms", false, mesh);
+  //DataAdaptor->AddArray(mesh, "atoms", vtkDataObject::FIELD_ASSOCIATION_POINTS, "type");
+
+  DataAdaptor->ReleaseData();
+}
+
+//-----------------------------------------------------------------------------
+void Finalize()
+{
+  AnalysisAdaptor->Finalize();
+  AnalysisAdaptor = nullptr;
+  DataAdaptor = nullptr;
+}
+
+}	// namespace lammpsBridge
diff --git a/miniapps/lammps/lammpsBridge.h b/miniapps/lammps/lammpsBridge.h
new file mode 100644
index 000000000..fe4987e44
--- /dev/null
+++ b/miniapps/lammps/lammpsBridge.h
@@ -0,0 +1,15 @@
+#pragma once
+#include <mpi.h>
+#include <string>
+
+namespace lammpsBridge
+{
+  void Initialize(MPI_Comm world, const std::string& config_file);
+  void SetData(long ntimestep, int nlocal, int *id, 
+               int    nghost, int *type, double **x, 
+               double xsublo, double xsubhi, 
+               double ysublo, double ysubhi, 
+               double zsublo, double zsubhi );
+  void Analyze();
+  void Finalize();
+}
diff --git a/miniapps/lammps/lammpsDataAdaptor.cxx b/miniapps/lammps/lammpsDataAdaptor.cxx
new file mode 100644
index 000000000..1f8a67a04
--- /dev/null
+++ b/miniapps/lammps/lammpsDataAdaptor.cxx
@@ -0,0 +1,416 @@
+#include "lammpsDataAdaptor.h"
+#include "Error.h"
+#include <vtkObjectFactory.h>
+#include <vtkSmartPointer.h>
+#include <vtkPolyData.h>
+#include <vtkPointData.h>
+#include <vtkDoubleArray.h>
+#include <vtkIntArray.h>
+#include <vtkMultiBlockDataSet.h>
+#include <vtkCellArray.h>
+
+namespace senseiLammps
+{
+
+struct lammpsDataAdaptor::DInternals
+{
+  vtkSmartPointer<vtkMultiBlockDataSet> mesh;
+  vtkSmartPointer<vtkDoubleArray> AtomPositions;
+  vtkSmartPointer<vtkIntArray> AtomTypes;
+  vtkSmartPointer<vtkIntArray> AtomIDs;
+  vtkSmartPointer<vtkCellArray> vertices;
+  double xsublo, ysublo, zsublo, xsubhi, ysubhi, zsubhi;
+  int nlocal, nghost;
+  double **x;
+  int *type;
+  int *id;
+};
+
+//-----------------------------------------------------------------------------
+senseiNewMacro(lammpsDataAdaptor);
+
+//-----------------------------------------------------------------------------
+lammpsDataAdaptor::lammpsDataAdaptor() :
+  Internals(new lammpsDataAdaptor::DInternals())
+{
+}
+
+//-----------------------------------------------------------------------------
+lammpsDataAdaptor::~lammpsDataAdaptor()
+{
+  delete this->Internals;
+}
+
+//-----------------------------------------------------------------------------
+void lammpsDataAdaptor::Initialize()
+{  
+  this->ReleaseData();
+}
+
+//-----------------------------------------------------------------------------
+void lammpsDataAdaptor::AddLAMMPSData( long ntimestep, int nlocal, int *id, 
+                                      int nghost, int *type, double **x, 
+                                      double xsublo, double xsubhi, 
+                                      double ysublo, double ysubhi, 
+                                      double zsublo, double zsubhi)
+{
+  DInternals& internals = (*this->Internals);
+
+  if(!internals.AtomPositions)
+    {
+    internals.AtomPositions = vtkSmartPointer<vtkDoubleArray>::New();
+    }
+
+  if(!internals.AtomTypes)
+    {
+    internals.AtomTypes = vtkSmartPointer<vtkIntArray>::New();
+    }
+
+  if(!internals.AtomIDs)
+    {
+    internals.AtomIDs = vtkSmartPointer<vtkIntArray>::New();
+    }
+
+  if(!internals.vertices)
+    {
+    internals.vertices = vtkSmartPointer<vtkCellArray>::New();
+    }
+
+  // atom coordinates
+  if (internals.AtomPositions)
+    {
+    //long nvals = nlocal + nghost;
+    long nvals = nlocal;
+
+    internals.AtomPositions->SetNumberOfComponents(3);
+    internals.AtomPositions->SetArray(*x, nvals*3, 1); 
+    internals.AtomPositions->SetName("positions");
+
+    internals.x = x;
+    }
+  else 
+    {
+    SENSEI_ERROR("Error. Internal AtomPositions structure not initialized")
+    }
+
+  // atom types
+  if (internals.AtomTypes)
+    {
+    //long nvals = nlocal + nghost;
+    long nvals = nlocal;
+
+    internals.AtomTypes->SetNumberOfComponents(1);
+    internals.AtomTypes->SetArray(type, nvals, 1);
+    internals.AtomTypes->SetName("type");
+
+    internals.type = type;  
+    }
+  else 
+    {
+    SENSEI_ERROR("Error. Internal AtomTypes structure not initialized")
+    }
+
+  // atom IDs
+  if (internals.AtomIDs)
+    {
+    //long nvals = nlocal + nghost;
+    long nvals = nlocal;
+
+    internals.AtomIDs->SetNumberOfComponents(1);
+    internals.AtomIDs->SetArray(id, nvals, 1);
+    internals.AtomIDs->SetName("id");
+
+    internals.id = id;  
+    }
+  else 
+    {
+    SENSEI_ERROR("Error. Internal AtomIDs structure not initialized")
+    }
+
+  // vertices
+  if (internals.vertices)
+    {
+    vtkIdType pid[1] = {0};
+
+    //for( int i=0; i < nlocal+nghost; i++) {
+    for( int i=0; i < nlocal; i++) 
+      {
+      internals.vertices->InsertNextCell (1, pid);
+      pid[0]++;
+      }
+    }
+
+  // number of atoms
+  internals.nlocal = nlocal;
+  internals.nghost = nghost;
+
+  // bounding box
+  internals.xsublo = xsublo;
+  internals.ysublo = ysublo;
+  internals.zsublo = zsublo;
+  internals.xsubhi = xsubhi;
+  internals.ysubhi = ysubhi;
+  internals.zsubhi = zsubhi;
+
+  // timestep
+  this->SetDataTimeStep(ntimestep);
+
+}
+
+//-----------------------------------------------------------------------------
+void lammpsDataAdaptor::GetBounds ( double &xsublo, double &xsubhi, 
+                                   double &ysublo, double &ysubhi, 
+                                   double &zsublo, double &zsubhi)
+{
+  DInternals& internals = (*this->Internals);
+
+  xsublo = internals.xsublo;
+  ysublo = internals.ysublo;
+  zsublo = internals.zsublo;
+  xsubhi = internals.xsubhi;
+  ysubhi = internals.ysubhi;
+  zsubhi = internals.zsubhi;
+}      
+
+void lammpsDataAdaptor::GetN ( int &nlocal, int &nghost )
+{
+  DInternals& internals = (*this->Internals);
+  
+  nlocal = internals.nlocal;
+  nghost = internals.nghost;
+}
+
+void lammpsDataAdaptor::GetPointers ( double **&x, int *&type)
+{
+  DInternals& internals = (*this->Internals);
+  
+  x = internals.x;
+  type = internals.type;
+}
+
+//-----------------------------------------------------------------------------
+void lammpsDataAdaptor::GetAtoms ( vtkDoubleArray *&atoms)
+{
+  DInternals& internals = (*this->Internals);
+
+  if (internals.AtomPositions)
+    atoms = internals.AtomPositions;
+  else
+    SENSEI_ERROR("Trying to get atom position array before setting it")
+}
+
+void lammpsDataAdaptor::GetTypes ( vtkIntArray *&types)
+{
+  DInternals& internals = (*this->Internals);
+
+  if (internals.AtomTypes)
+    types = internals.AtomTypes;
+  else
+    SENSEI_ERROR("Trying to get atom type array before setting it")
+}
+
+void lammpsDataAdaptor::GetIDs ( vtkIntArray *&ids)
+{
+  DInternals& internals = (*this->Internals);
+
+  if (internals.AtomIDs)
+    ids = internals.AtomIDs;
+  else
+    SENSEI_ERROR("Trying to get atom ID array before setting it")
+}
+
+//-----------------------------------------------------------------------------
+int lammpsDataAdaptor::GetNumberOfMeshes(unsigned int &numMeshes)
+{
+  numMeshes = 1;
+  return 0;
+}
+
+//-----------------------------------------------------------------------------
+int lammpsDataAdaptor::GetMesh(const std::string &meshName, bool structureOnly,
+    vtkDataObject *&mesh)
+{
+  if (meshName != "atoms") 
+    {  
+    SENSEI_ERROR("No mesh \"" << meshName << "\"")
+    return -1;
+    }
+
+  DInternals& internals = (*this->Internals);
+
+  if (!internals.mesh)
+    {
+    vtkSmartPointer<vtkPolyData> pd = vtkSmartPointer<vtkPolyData>::New();
+
+    if(!structureOnly)
+      {
+      vtkSmartPointer<vtkPoints> pts = vtkSmartPointer<vtkPoints>::New();
+      pts->SetData(internals.AtomPositions);
+      pd->SetPoints(pts);
+      }
+
+    pd->SetVerts( internals.vertices );
+
+    int rank, size; 
+    MPI_Comm comm;
+  	
+    comm = GetCommunicator();
+    MPI_Comm_rank(comm, &rank);
+    MPI_Comm_size(comm, &size); 
+
+    internals.mesh = vtkSmartPointer<vtkMultiBlockDataSet>::New();
+    internals.mesh->SetNumberOfBlocks(size);
+    internals.mesh->SetBlock(rank, pd);
+    }
+
+  mesh = internals.mesh;
+
+  return 0;
+}
+
+//-----------------------------------------------------------------------------
+int lammpsDataAdaptor::AddArray(vtkDataObject* mesh, const std::string &meshName,
+    int association, const std::string &arrayName)
+{
+  if (meshName != "atoms")
+    {
+    SENSEI_ERROR("No mesh \"" << meshName << "\"")
+    return -1;
+    }
+
+  if (association != vtkDataObject::FIELD_ASSOCIATION_POINTS)
+    {
+    SENSEI_ERROR("No cell data on mesh")
+    return -1;
+    }
+
+  if (arrayName == "type")
+    {  	
+    DInternals& internals = (*this->Internals);
+    vtkMultiBlockDataSet* md = vtkMultiBlockDataSet::SafeDownCast(mesh);
+    vtkSmartPointer<vtkPolyData> pd = vtkSmartPointer<vtkPolyData>::New();
+  
+    int rank;
+    MPI_Comm comm;
+
+    comm = GetCommunicator();
+    MPI_Comm_rank(comm, &rank);
+
+    pd = vtkPolyData::SafeDownCast(md->GetBlock(rank));
+    pd->GetPointData()->AddArray(internals.AtomTypes);
+    }
+
+  if (arrayName == "id")
+    {	
+    DInternals& internals = (*this->Internals);
+    vtkMultiBlockDataSet* md = vtkMultiBlockDataSet::SafeDownCast(mesh);
+    vtkSmartPointer<vtkPolyData> pd = vtkSmartPointer<vtkPolyData>::New();
+  
+    int rank;
+    MPI_Comm comm;
+
+    comm = GetCommunicator();
+    MPI_Comm_rank(comm, &rank);
+
+    pd = vtkPolyData::SafeDownCast(md->GetBlock(rank));
+    pd->GetPointData()->AddArray(internals.AtomIDs);
+    }
+
+  return 0;  
+}
+
+// not implemented
+//----------------------------------------------------------------------------
+int lammpsDataAdaptor::AddGhostCellsArray(vtkDataObject *mesh, const std::string &meshName)
+{
+  (void) mesh;
+  (void) meshName;
+  return 0;
+}
+
+//-----------------------------------------------------------------------------
+int lammpsDataAdaptor::GetMeshMetadata(unsigned int id, sensei::MeshMetadataPtr &metadata)
+{
+
+  if (id > 0)
+    {
+    SENSEI_ERROR("invalid mesh id " << id)
+    return -1;
+    }
+
+  int rank, nRanks;
+  MPI_Comm comm;
+
+  comm = GetCommunicator();
+  MPI_Comm_rank(comm, &rank);
+  MPI_Comm_size(comm, &nRanks);
+
+
+  metadata->MeshName = "atoms";
+  metadata->MeshType = VTK_MULTIBLOCK_DATA_SET;
+  metadata->BlockType = VTK_POLY_DATA;
+  metadata->CoordinateType = VTK_DOUBLE;
+  metadata->NumBlocks = nRanks;
+  metadata->NumBlocksLocal = {1};
+  metadata->NumGhostCells = this->Internals->nghost;
+  metadata->NumArrays = 2;
+  metadata->ArrayName = {"type", "id"};
+  metadata->ArrayCentering = {vtkDataObject::POINT, vtkDataObject::POINT};
+  metadata->ArrayComponents = {1, 1};
+  metadata->ArrayType = {VTK_INT, VTK_INT};
+  metadata->StaticMesh = 1;
+
+  if (metadata->Flags.BlockExtentsSet())
+    {
+    SENSEI_WARNING("lammps data adaptor. Flags.BlockExtentsSet()")
+    }
+
+  if (metadata->Flags.BlockBoundsSet())
+    {
+    SENSEI_WARNING("lammps data adaptor. Flags.BlockBoundsSet()")
+    }
+
+  if (metadata->Flags.BlockSizeSet())
+    {
+    SENSEI_WARNING("lammps data adaptor. Flags.BlockSizeSet()")
+    }
+
+  if (metadata->Flags.BlockDecompSet())
+    {
+    metadata->BlockOwner.push_back(rank);
+    metadata->BlockIds.push_back(rank);
+    }
+
+  if (metadata->Flags.BlockArrayRangeSet())
+    {
+    SENSEI_WARNING("lammps data adaptor. Flags.BlockArrayRangeSet()")
+    }
+
+  return 0;
+}
+
+
+
+//-----------------------------------------------------------------------------
+int lammpsDataAdaptor::ReleaseData()
+{
+  DInternals& internals = (*this->Internals);
+
+  internals.mesh = NULL;
+  internals.AtomPositions = NULL;
+  internals.AtomTypes = NULL;
+  internals.AtomIDs = NULL;
+  internals.nlocal = 0;
+  internals.nghost = 0;
+  internals.xsublo = 0;
+  internals.ysublo = 0;
+  internals.zsublo = 0;
+  internals.xsubhi = 0;
+  internals.ysubhi = 0;
+  internals.zsubhi = 0;
+
+  return 0;
+}
+
+} //senseiLammps
+
diff --git a/miniapps/lammps/lammpsDataAdaptor.h b/miniapps/lammps/lammpsDataAdaptor.h
new file mode 100644
index 000000000..79929e0e8
--- /dev/null
+++ b/miniapps/lammps/lammpsDataAdaptor.h
@@ -0,0 +1,78 @@
+#pragma once
+
+#include <DataAdaptor.h>
+#include <vtkDoubleArray.h>
+#include <vtkIntArray.h>
+
+namespace senseiLammps
+{
+
+class lammpsDataAdaptor : public sensei::DataAdaptor
+{
+public:
+  static lammpsDataAdaptor* New();
+  senseiTypeMacro(lammpsDataAdaptor, sensei::DataAdaptor);
+
+  /// @brief Initialize the data adaptor.
+  ///
+  /// This initializes the data adaptor. This must be called once per simulation run.
+  void Initialize();
+
+  void AddLAMMPSData( long ntimestep, int nlocal, int *id, 
+                      int nghost, int *type, double **x, 
+                      double xsublo, double xsubhi, 
+                      double ysublo, double ysubhi, 
+                      double zsublo, double zsubhi);
+
+  void GetBounds ( double &xsublo, double &xsubhi, 
+                   double &ysublo, double &ysubhi, 
+                   double &zsublo, double &zsubhi);
+
+  void GetN ( int &nlocal, int &nghost );
+
+  void GetPointers ( double **&x, int *&type);
+
+  void GetAtoms ( vtkDoubleArray *&atoms );
+
+  void GetTypes ( vtkIntArray *&types );
+
+  void GetIDs ( vtkIntArray *&ids );
+
+// SENSEI API
+  int GetNumberOfMeshes(unsigned int &numMeshes) override;
+
+  int GetMeshMetadata(unsigned int id, sensei::MeshMetadataPtr &md) override;
+
+  //int GetMeshName(unsigned int id, std::string &meshName) override;
+
+  int GetMesh(const std::string &meshName, bool structureOnly,
+    vtkDataObject *&mesh) override;
+
+  int AddArray(vtkDataObject* mesh, const std::string &meshName,
+    int association, const std::string &arrayName) override;
+
+  //int GetNumberOfArrays(const std::string &meshName, int association,
+  //  unsigned int &numberOfArrays) override;
+
+  //int GetArrayName(const std::string &meshName, int association,
+  //  unsigned int index, std::string &arrayName) override;
+
+  //int GetMeshHasGhostCells(const std::string &meshName, int &nLayers) override;
+
+  int AddGhostCellsArray(vtkDataObject* mesh, const std::string &meshName) override;
+
+  int ReleaseData() override;
+
+protected:
+  lammpsDataAdaptor();
+  ~lammpsDataAdaptor();
+
+private:
+  lammpsDataAdaptor(const lammpsDataAdaptor&); // not implemented.
+  void operator=(const lammpsDataAdaptor&); // not implemented.
+
+  struct DInternals;
+  DInternals* Internals;
+};
+
+}
diff --git a/miniapps/lammps/lammpsDriver.cxx b/miniapps/lammps/lammpsDriver.cxx
new file mode 100644
index 000000000..dc4134b4a
--- /dev/null
+++ b/miniapps/lammps/lammpsDriver.cxx
@@ -0,0 +1,218 @@
+/* ----------------------------------------------------------------------
+   LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
+   www.cs.sandia.gov/~sjplimp/lammps.html
+   Steve Plimpton, sjplimp@sandia.gov, Sandia National Laboratories
+
+   Copyright (2003) Sandia Corporation.  Under the terms of Contract
+   DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
+   certain rights in this software.  This software is distributed under 
+   the GNU General Public License.
+
+   See the README file in the top-level LAMMPS directory.
+------------------------------------------------------------------------- */
+
+#include <string>
+#include <fstream>
+#include <vector>
+#include <iostream>
+#include <mpi.h>
+
+// LAMMPS include files
+#include "library.h"         
+#include "domain.h"         
+#include "lammps.h"
+#include "modify.h"
+#include "fix.h"
+#include "fix_external.h"
+
+// SENSEI bridge 
+#include "lammpsBridge.h"   
+
+using namespace LAMMPS_NS;
+
+struct Info {
+  int me;
+  LAMMPS *lmp;
+} globalInfo;
+
+void insituCallback(void *ptr, bigint ntimestep, int nlocal, 
+                    int *id, double **x, double **f) 
+{
+  (void) x;
+  (void) f;
+ 
+  Info *info = (Info *) ptr;
+
+  double xsublo = info->lmp->domain->sublo[0];
+  double xsubhi = info->lmp->domain->subhi[0];
+  double ysublo = info->lmp->domain->sublo[1];
+  double ysubhi = info->lmp->domain->subhi[1];
+  double zsublo = info->lmp->domain->sublo[2];
+  double zsubhi = info->lmp->domain->subhi[2];
+
+  const char *typ = "type";
+  int *type = (int *)lammps_extract_atom(info->lmp,(char *)typ);
+
+  const char *ng = "nghost";
+  int *nghost = (int*)lammps_extract_global(info->lmp,(char *)ng);
+
+  // The arrays we get through extract atom include the ghost atoms as well
+  const char *xx = "x";
+  double **all_pos = (double**)lammps_extract_atom(info->lmp,(char *)xx);
+
+  lammpsBridge::SetData(ntimestep, nlocal, id, *nghost, type, all_pos, 
+                        xsublo, xsubhi, ysublo, ysubhi, zsublo, zsubhi);
+
+  if (0 == globalInfo.me)
+    std::cout << "###### SENSEI instrumentation: bridge analyze() ######" << std::endl;    
+  lammpsBridge::Analyze();
+
+  if (0 == globalInfo.me)
+    std::cout << "###### SENSEI instrumentation: after bridge analyze() ######" << std::endl;    
+ 
+
+
+}
+
+const static std::string USAGE =
+"Usage: ./driver <lammps input> [options] [-lmp <lammps args>]\n"
+"Options:\n"
+"  -h                  Print this help text\n"
+"  -sensei <file.xml>  Pass this to set the SENSEI xml config file\n"
+"  -n <int>            Specify the number of steps to simulate for. Default is 10000\n"
+"  -lmp <lammps args>  Pass the list of arguments <args> to LAMMPS as if they were\n"
+"                      command line args to LAMMPS. This must be the last argument, all\n"
+"                      following arguments will be passed to lammps.\n"
+"  -log                Generate full time and memory usage log.\n"
+"  -shortlog           Generate a summary time and memory usage log.";
+
+int main(int argc, char **argv) 
+{
+  std::vector<std::string> args(argv, argv + argc);
+  if (args.size() < 2) 
+    {
+    std::cout << USAGE << "\n";
+    return 1;
+    }
+  
+  std::string lammps_input = args[1];
+	
+  // Additional args the user may be passing to lammps
+  std::vector<char*> lammps_args(1, argv[0]);
+  bool log = false;
+  bool shortlog = false;
+  std::string sensei_xml;
+  size_t sim_steps = 10000;
+  for (size_t i = 2; i < args.size(); ++i) 
+    {
+    if (args[i] == "-sensei") 
+      sensei_xml = args[++i];
+    else if (args[i] == "-n" ) 
+      sim_steps = std::stoull(args[++i]);
+    else if (args[i] == "-h") 
+      {
+      std::cout << USAGE << "\n";
+      return 0;
+      } 
+    else if (args[i] == "-log" ) 
+      log = true;
+    else if (args[i] == "-shortlog" ) 
+      shortlog = true;
+    else if (args[i] == "-lmp") 
+      {
+      ++i;
+      for (; i < args.size(); ++i) 
+	lammps_args.push_back(&args[i][0]);
+      break;
+      }
+    }
+
+  MPI_Init(&argc, &argv);
+  MPI_Comm sim_comm = MPI_COMM_WORLD;
+  MPI_Comm_rank(sim_comm, &globalInfo.me);
+
+  // Initialize SENSEI bridge 
+  if (0 == globalInfo.me) 
+    std::cout << "###### SENSEI instrumentation: initialize bridge ######" << std::endl;
+  lammpsBridge::Initialize(sim_comm, sensei_xml );
+
+  LAMMPS *lammps;
+  lammps_open(lammps_args.size(), lammps_args.data(), sim_comm, (void**)&lammps);
+  globalInfo.lmp = lammps;
+
+  // run the input script thru LAMMPS one line at a time until end-of-file
+  // driver proc 0 reads a line, Bcasts it to all procs
+  // (could just send it to proc 0 of comm_lammps and let it Bcast)
+  // all LAMMPS procs call lammps_command() on the line
+
+  if (0 == globalInfo.me) 
+    {
+    std::cout << "Loading lammps input: '" << lammps_input << "'\n";
+    std::ifstream input(lammps_input.c_str());
+    for (std::string line; std::getline(input, line);) 
+      {
+      int len = line.size();
+      
+      // Skip empty lines
+      if (len == 0) 
+        continue;
+      
+      MPI_Bcast(&len, 1, MPI_INT, 0, sim_comm);
+      MPI_Bcast(&line[0], len, MPI_CHAR, 0, sim_comm);
+                lammps_command(lammps, &line[0]);
+      }
+      
+    // Bcast out we're done with the file
+    int len = 0;
+    MPI_Bcast(&len, 1, MPI_INT, 0, sim_comm);
+    } 
+  else 
+    {
+    while (true) 
+      {
+      int len = 0;
+      MPI_Bcast(&len, 1, MPI_INT, 0, sim_comm);
+      if (len == 0) 
+        break;
+      else 
+        {
+	std::vector<char> line(len + 1, '\0');
+	MPI_Bcast(line.data(), len, MPI_CHAR, 0, sim_comm);
+		  lammps_command(lammps, line.data());
+	}
+      }
+    }
+
+  // Setup the fix external callback
+  int ifix = lammps->modify->find_fix_by_style("external");
+	
+  // If there's no external fix, abort
+  if (ifix == -1) 
+    {
+    if (0 == globalInfo.me) 
+      std::cout << "You need to add a fix external line to the input file. Abort" << std::endl;
+    return -1;
+    }
+  FixExternal *fix = (FixExternal*)lammps->modify->fix[ifix];
+  fix->set_callback(insituCallback, &globalInfo);
+
+  // run for a number of steps
+  for (size_t i = 0; i < sim_steps; ++i) 
+    {
+    const char * string = "run 1 pre no post no";
+    lammps_command(lammps, (char *)string);
+    }
+
+  // all LAMMPS timesteps computed
+  lammps_close(lammps);
+
+  // Finalize SENSEI bridge 
+  if (0 == globalInfo.me) 
+    std::cout << "###### SENSEI instrumentation: finalize bridge ######" << std::endl;
+  lammpsBridge::Finalize();
+
+  // close down MPI
+  MPI_Finalize();
+
+  return 0;
+}

From c9b8be955b5e76380ffca2ca7495a0016592fbec Mon Sep 17 00:00:00 2001
From: Silvio Rizzi <srizzi@anl.gov>
Date: Tue, 3 Aug 2021 22:14:55 +0000
Subject: [PATCH 02/11] added metadata changes from victor

---
 miniapps/lammps/lammpsDataAdaptor.cxx | 128 ++++++++++++++++++++------
 miniapps/lammps/lammpsDataAdaptor.h   |   6 ++
 2 files changed, 105 insertions(+), 29 deletions(-)

diff --git a/miniapps/lammps/lammpsDataAdaptor.cxx b/miniapps/lammps/lammpsDataAdaptor.cxx
index 1f8a67a04..1c2bd88f4 100644
--- a/miniapps/lammps/lammpsDataAdaptor.cxx
+++ b/miniapps/lammps/lammpsDataAdaptor.cxx
@@ -9,6 +9,41 @@
 #include <vtkMultiBlockDataSet.h>
 #include <vtkCellArray.h>
 
+static
+std::array<int,2> getArrayRange(unsigned long nSize, int *x) {
+  int xmin = std::numeric_limits<int>::max(); 
+  int xmax = std::numeric_limits<int>::lowest();
+  for(int i=0; i<nSize; ++i) {
+    xmin = std::min(xmin, x[i]);
+    xmax = std::max(xmax, x[i]);
+  }
+
+  return {xmin, xmax};
+}
+
+static
+std::array<double,2> getArrayRange(unsigned long nSize,double *x) {
+  double xmin = std::numeric_limits<double>::max(); 
+  double xmax = std::numeric_limits<double>::lowest();
+  for(int i=0; i<nSize; ++i) {
+    xmin = std::min(xmin, x[i]);
+    xmax = std::max(xmax, x[i]);
+  }
+
+  return {xmin, xmax};
+}
+
+static 
+void getBounds(const sdiy::DiscreteBounds &db, double *ext) 
+{
+  ext[0] = db.min[0];
+  ext[1] = db.max[0];
+  ext[2] = db.min[1];
+  ext[3] = db.max[1];
+  ext[4] = db.min[2];
+  ext[5] = db.max[2];
+}
+
 namespace senseiLammps
 {
 
@@ -19,7 +54,10 @@ struct lammpsDataAdaptor::DInternals
   vtkSmartPointer<vtkIntArray> AtomTypes;
   vtkSmartPointer<vtkIntArray> AtomIDs;
   vtkSmartPointer<vtkCellArray> vertices;
-  double xsublo, ysublo, zsublo, xsubhi, ysubhi, zsubhi;
+  sdiy::DiscreteBounds DomainBounds;
+  sdiy::DiscreteBounds BlockBounds;
+  sdiy::DiscreteBounds typeRange;
+  sdiy::DiscreteBounds idRange;
   int nlocal, nghost;
   double **x;
   int *type;
@@ -144,33 +182,42 @@ void lammpsDataAdaptor::AddLAMMPSData( long ntimestep, int nlocal, int *id,
   internals.nlocal = nlocal;
   internals.nghost = nghost;
 
+  std::array<double,2> x_range = getArrayRange(nlocal, x[0]);
+	std::array<double,2> y_range = getArrayRange(nlocal, x[1]);
+	std::array<double,2> z_range = getArrayRange(nlocal, x[2]);
+
   // bounding box
-  internals.xsublo = xsublo;
-  internals.ysublo = ysublo;
-  internals.zsublo = zsublo;
-  internals.xsubhi = xsubhi;
-  internals.ysubhi = ysubhi;
-  internals.zsubhi = zsubhi;
+  this->SetDomainBounds(xsublo, xsubhi, ysublo, ysubhi, zsublo, zsubhi);
+  this->SetBlockBounds(
+      x_range[0], x_range[1],
+      y_range[0], y_range[1],
+      z_range[0], z_range[1]);
+
+  /// XXX Set type and id range
+  this->Internals->typeRange.min[0] = std::numeric_limits<int>::max();
+  this->Internals->typeRange.max[0] = std::numeric_limits<int>::min();
+  this->Internals->idRange.min[0] = std::numeric_limits<int>::max();
+  this->Internals->idRange.max[0] = std::numeric_limits<int>::min();
 
   // timestep
   this->SetDataTimeStep(ntimestep);
 
 }
 
-//-----------------------------------------------------------------------------
-void lammpsDataAdaptor::GetBounds ( double &xsublo, double &xsubhi, 
-                                   double &ysublo, double &ysubhi, 
-                                   double &zsublo, double &zsubhi)
-{
-  DInternals& internals = (*this->Internals);
+void lammpsDataAdaptor::SetBlockBounds(double *x, int nelem) {
+  this-Internals->  
+}
 
-  xsublo = internals.xsublo;
-  ysublo = internals.ysublo;
-  zsublo = internals.zsublo;
-  xsubhi = internals.xsubhi;
-  ysubhi = internals.ysubhi;
-  zsubhi = internals.zsubhi;
-}      
+void lammpsDataAdaptor::SetDomainBounds(double xmin, double xmax,
+    double ymin, double ymax, double zmin, double zmax) {
+  this->Internals->DomainBounds.min[0] = xmin;
+  this->Internals->DomainBounds.min[1] = ymin;
+  this->Internals->DomainBounds.min[2] = zmin;
+
+  this->Internals->DomainBounds.max[0] = xmax;
+  this->Internals->DomainBounds.max[1] = ymax;
+  this->Internals->DomainBounds.max[2] = zmax;
+}
 
 void lammpsDataAdaptor::GetN ( int &nlocal, int &nghost )
 {
@@ -346,33 +393,52 @@ int lammpsDataAdaptor::GetMeshMetadata(unsigned int id, sensei::MeshMetadataPtr
   MPI_Comm_size(comm, &nRanks);
 
 
+  int nBlocks = 1; // One block per rank
   metadata->MeshName = "atoms";
   metadata->MeshType = VTK_MULTIBLOCK_DATA_SET;
   metadata->BlockType = VTK_POLY_DATA;
   metadata->CoordinateType = VTK_DOUBLE;
   metadata->NumBlocks = nRanks;
-  metadata->NumBlocksLocal = {1};
+  metadata->NumBlocksLocal = {nBlocks};
   metadata->NumGhostCells = this->Internals->nghost;
   metadata->NumArrays = 2;
   metadata->ArrayName = {"type", "id"};
   metadata->ArrayCentering = {vtkDataObject::POINT, vtkDataObject::POINT};
   metadata->ArrayComponents = {1, 1};
   metadata->ArrayType = {VTK_INT, VTK_INT};
-  metadata->StaticMesh = 1;
+  metadata->StaticMesh = 0;
 
   if (metadata->Flags.BlockExtentsSet())
     {
     SENSEI_WARNING("lammps data adaptor. Flags.BlockExtentsSet()")
+    // There should be no extent for a PolyData, but ADIOS2 needs this
+    std::array<int,6> ext = { 0, 0, 0, 0, 0, 0};
+    metadata->Extent = std::move(ext);
+
+    metadata->BlockExtents.reserve(nBlocks);
+    metadata->BlockExtents.emplace_back(std::move(ext));
     }
 
   if (metadata->Flags.BlockBoundsSet())
     {
     SENSEI_WARNING("lammps data adaptor. Flags.BlockBoundsSet()")
+    std::array<double, 6> bounds;
+    getBounds(this->Internals->DomainBounds, bounds.data());
+    metadata->Bounds = std::move(bounds);
+    
+    metadata->BlockBounds.reserve(nBlocks);
+    
+    getBounds(this->Internals->BlockBounds, bounds.data());
+    metadata->BlockBounds.emplace_back(std::move(bounds));
     }
 
   if (metadata->Flags.BlockSizeSet())
     {
+    int nCells = nlocal;
     SENSEI_WARNING("lammps data adaptor. Flags.BlockSizeSet()")
+    metadata->BlockNumCells.push_back(nCells);
+    metadata->BlockNumPoints.push_back(nCells);
+    metadata->BlockCellArraySize.push_back(2 * nCells); // XXX- VTK_POINTS
     }
 
   if (metadata->Flags.BlockDecompSet())
@@ -384,7 +450,16 @@ int lammpsDataAdaptor::GetMeshMetadata(unsigned int id, sensei::MeshMetadataPtr
   if (metadata->Flags.BlockArrayRangeSet())
     {
     SENSEI_WARNING("lammps data adaptor. Flags.BlockArrayRangeSet()")
-    }
+    
+    std::array<int,2> typeBlockRange = getArrayRange(nvals, this->Internals->type);
+    std::array<int,2> idBlockRange = getArrayRange(nvals, this->Internals->id);
+    metadata->BlockArrayRange.push_back({typeBlockRange, idBlockRange});
+
+		std::array<int,2> typeRange = { this->Internals->typeRange.min[0], this->Internals->typeRange.max[0] };
+    std::array<int,2> idRange = { this->Internals->idRange.min[0], this->Internals->idRange.max[0] };
+    metadata->ArrayRange.push_back(typeRange); 
+    metadata->ArrayRange.push_back(idRange); 
+}
 
   return 0;
 }
@@ -402,15 +477,10 @@ int lammpsDataAdaptor::ReleaseData()
   internals.AtomIDs = NULL;
   internals.nlocal = 0;
   internals.nghost = 0;
-  internals.xsublo = 0;
-  internals.ysublo = 0;
-  internals.zsublo = 0;
-  internals.xsubhi = 0;
-  internals.ysubhi = 0;
-  internals.zsubhi = 0;
 
   return 0;
 }
 
 } //senseiLammps
 
+
diff --git a/miniapps/lammps/lammpsDataAdaptor.h b/miniapps/lammps/lammpsDataAdaptor.h
index 79929e0e8..43fe77f84 100644
--- a/miniapps/lammps/lammpsDataAdaptor.h
+++ b/miniapps/lammps/lammpsDataAdaptor.h
@@ -38,6 +38,9 @@ class lammpsDataAdaptor : public sensei::DataAdaptor
 
   void GetIDs ( vtkIntArray *&ids );
 
+  void lammpsDataAdaptor::SetDomainBounds(double xmin, double xmax,
+    double ymin, double ymax, double zmin, double zmax);
+
 // SENSEI API
   int GetNumberOfMeshes(unsigned int &numMeshes) override;
 
@@ -76,3 +79,6 @@ class lammpsDataAdaptor : public sensei::DataAdaptor
 };
 
 }
+
+
+

From 072f22feaa1e31f8d6270236d731b5477f57fc58 Mon Sep 17 00:00:00 2001
From: Silvio Rizzi <srizzi@anl.gov>
Date: Tue, 3 Aug 2021 22:46:14 +0000
Subject: [PATCH 03/11] lammps fix for open_lammps by adding compiler flag
 -DLAMMPS_LIB_MPI. Documented by others here
 https://github.com/lammps/lammps/issues/2494

---
 miniapps/lammps/CMakeLists.txt   | 2 ++
 miniapps/lammps/lammpsDriver.cxx | 3 +--
 2 files changed, 3 insertions(+), 2 deletions(-)

diff --git a/miniapps/lammps/CMakeLists.txt b/miniapps/lammps/CMakeLists.txt
index 10094321d..918ef11f9 100644
--- a/miniapps/lammps/CMakeLists.txt
+++ b/miniapps/lammps/CMakeLists.txt
@@ -6,6 +6,8 @@ endif()
 
 include_directories(${LAMMPS_DIR}/src)
 
+set(CMAKE_CXX_FLAGS  "${CMAKE_CXX_FLAGS} -DLAMMPS_LIB_MPI")
+
 set(sources lammpsDriver.cxx)
 
 set(libs 
diff --git a/miniapps/lammps/lammpsDriver.cxx b/miniapps/lammps/lammpsDriver.cxx
index dc4134b4a..90dc8860e 100644
--- a/miniapps/lammps/lammpsDriver.cxx
+++ b/miniapps/lammps/lammpsDriver.cxx
@@ -136,8 +136,7 @@ int main(int argc, char **argv)
     std::cout << "###### SENSEI instrumentation: initialize bridge ######" << std::endl;
   lammpsBridge::Initialize(sim_comm, sensei_xml );
 
-  LAMMPS *lammps;
-  lammps_open(lammps_args.size(), lammps_args.data(), sim_comm, (void**)&lammps);
+  LAMMPS *lammps = (LAMMPS *)lammps_open(lammps_args.size(), lammps_args.data(), sim_comm, NULL);
   globalInfo.lmp = lammps;
 
   // run the input script thru LAMMPS one line at a time until end-of-file

From c2a0c8b79cc8bb346a0b236cfe7a7830285969b4 Mon Sep 17 00:00:00 2001
From: Silvio Rizzi <srizzi@anl.gov>
Date: Tue, 3 Aug 2021 23:25:40 +0000
Subject: [PATCH 04/11] compile errors fixed in lammpsDataAdaptor so I can
 start debugging

---
 miniapps/lammps/lammpsDataAdaptor.cxx | 46 +++++++++++++++++----------
 miniapps/lammps/lammpsDataAdaptor.h   |  5 ++-
 2 files changed, 34 insertions(+), 17 deletions(-)

diff --git a/miniapps/lammps/lammpsDataAdaptor.cxx b/miniapps/lammps/lammpsDataAdaptor.cxx
index 1c2bd88f4..2032e73d5 100644
--- a/miniapps/lammps/lammpsDataAdaptor.cxx
+++ b/miniapps/lammps/lammpsDataAdaptor.cxx
@@ -8,12 +8,13 @@
 #include <vtkIntArray.h>
 #include <vtkMultiBlockDataSet.h>
 #include <vtkCellArray.h>
+#include <sdiy/master.hpp>
 
 static
 std::array<int,2> getArrayRange(unsigned long nSize, int *x) {
   int xmin = std::numeric_limits<int>::max(); 
   int xmax = std::numeric_limits<int>::lowest();
-  for(int i=0; i<nSize; ++i) {
+  for(unsigned int i=0; i<nSize; ++i) {
     xmin = std::min(xmin, x[i]);
     xmax = std::max(xmax, x[i]);
   }
@@ -25,7 +26,7 @@ static
 std::array<double,2> getArrayRange(unsigned long nSize,double *x) {
   double xmin = std::numeric_limits<double>::max(); 
   double xmax = std::numeric_limits<double>::lowest();
-  for(int i=0; i<nSize; ++i) {
+  for(unsigned int i=0; i<nSize; ++i) {
     xmin = std::min(xmin, x[i]);
     xmax = std::max(xmax, x[i]);
   }
@@ -183,15 +184,15 @@ void lammpsDataAdaptor::AddLAMMPSData( long ntimestep, int nlocal, int *id,
   internals.nghost = nghost;
 
   std::array<double,2> x_range = getArrayRange(nlocal, x[0]);
-	std::array<double,2> y_range = getArrayRange(nlocal, x[1]);
-	std::array<double,2> z_range = getArrayRange(nlocal, x[2]);
+  std::array<double,2> y_range = getArrayRange(nlocal, x[1]);
+  std::array<double,2> z_range = getArrayRange(nlocal, x[2]);
 
   // bounding box
   this->SetDomainBounds(xsublo, xsubhi, ysublo, ysubhi, zsublo, zsubhi);
   this->SetBlockBounds(
-      x_range[0], x_range[1],
-      y_range[0], y_range[1],
-      z_range[0], z_range[1]);
+    x_range[0], x_range[1],
+    y_range[0], y_range[1],
+    z_range[0], z_range[1]);
 
   /// XXX Set type and id range
   this->Internals->typeRange.min[0] = std::numeric_limits<int>::max();
@@ -204,12 +205,21 @@ void lammpsDataAdaptor::AddLAMMPSData( long ntimestep, int nlocal, int *id,
 
 }
 
-void lammpsDataAdaptor::SetBlockBounds(double *x, int nelem) {
-  this-Internals->  
+void lammpsDataAdaptor::SetBlockBounds(double xmin, double xmax,
+    double ymin, double ymax, double zmin, double zmax)
+{
+  this->Internals->BlockBounds.min[0] = xmin;
+  this->Internals->BlockBounds.min[1] = ymin;
+  this->Internals->BlockBounds.min[2] = zmin;
+
+  this->Internals->BlockBounds.max[0] = xmax;
+  this->Internals->BlockBounds.max[1] = ymax;
+  this->Internals->BlockBounds.max[2] = zmax;
 }
 
 void lammpsDataAdaptor::SetDomainBounds(double xmin, double xmax,
-    double ymin, double ymax, double zmin, double zmax) {
+    double ymin, double ymax, double zmin, double zmax) 
+{
   this->Internals->DomainBounds.min[0] = xmin;
   this->Internals->DomainBounds.min[1] = ymin;
   this->Internals->DomainBounds.min[2] = zmin;
@@ -434,7 +444,8 @@ int lammpsDataAdaptor::GetMeshMetadata(unsigned int id, sensei::MeshMetadataPtr
 
   if (metadata->Flags.BlockSizeSet())
     {
-    int nCells = nlocal;
+    DInternals& internals = (*this->Internals);
+    int nCells = internals.nlocal;
     SENSEI_WARNING("lammps data adaptor. Flags.BlockSizeSet()")
     metadata->BlockNumCells.push_back(nCells);
     metadata->BlockNumPoints.push_back(nCells);
@@ -450,15 +461,18 @@ int lammpsDataAdaptor::GetMeshMetadata(unsigned int id, sensei::MeshMetadataPtr
   if (metadata->Flags.BlockArrayRangeSet())
     {
     SENSEI_WARNING("lammps data adaptor. Flags.BlockArrayRangeSet()")
-    
+ 
+    unsigned long nvals;   
     std::array<int,2> typeBlockRange = getArrayRange(nvals, this->Internals->type);
     std::array<int,2> idBlockRange = getArrayRange(nvals, this->Internals->id);
-    metadata->BlockArrayRange.push_back({typeBlockRange, idBlockRange});
+    // fixme
+    //metadata->BlockArrayRange.push_back({typeBlockRange, idBlockRange});
 
-		std::array<int,2> typeRange = { this->Internals->typeRange.min[0], this->Internals->typeRange.max[0] };
+    std::array<int,2> typeRange = { this->Internals->typeRange.min[0], this->Internals->typeRange.max[0] };
     std::array<int,2> idRange = { this->Internals->idRange.min[0], this->Internals->idRange.max[0] };
-    metadata->ArrayRange.push_back(typeRange); 
-    metadata->ArrayRange.push_back(idRange); 
+    // fixme
+    //metadata->ArrayRange.push_back(typeRange); 
+    //metadata->ArrayRange.push_back(idRange); 
 }
 
   return 0;
diff --git a/miniapps/lammps/lammpsDataAdaptor.h b/miniapps/lammps/lammpsDataAdaptor.h
index 43fe77f84..dd214104b 100644
--- a/miniapps/lammps/lammpsDataAdaptor.h
+++ b/miniapps/lammps/lammpsDataAdaptor.h
@@ -38,7 +38,10 @@ class lammpsDataAdaptor : public sensei::DataAdaptor
 
   void GetIDs ( vtkIntArray *&ids );
 
-  void lammpsDataAdaptor::SetDomainBounds(double xmin, double xmax,
+  void SetDomainBounds(double xmin, double xmax,
+    double ymin, double ymax, double zmin, double zmax);
+
+  void SetBlockBounds(double xmin, double xmax,
     double ymin, double ymax, double zmin, double zmax);
 
 // SENSEI API

From 1371a56d669f1aa4679f89eb5cbbd16559ad6eec Mon Sep 17 00:00:00 2001
From: Silvio Rizzi <srizzi@anl.gov>
Date: Thu, 5 Aug 2021 03:45:40 +0000
Subject: [PATCH 05/11] lammps bug fixed. Now I can save VTK files from lammps
 with the configurable analysis adaptor

---
 miniapps/lammps/lammpsDataAdaptor.cxx | 64 +++++++++++++--------------
 1 file changed, 32 insertions(+), 32 deletions(-)

diff --git a/miniapps/lammps/lammpsDataAdaptor.cxx b/miniapps/lammps/lammpsDataAdaptor.cxx
index 2032e73d5..2dd3ec119 100644
--- a/miniapps/lammps/lammpsDataAdaptor.cxx
+++ b/miniapps/lammps/lammpsDataAdaptor.cxx
@@ -50,7 +50,7 @@ namespace senseiLammps
 
 struct lammpsDataAdaptor::DInternals
 {
-  vtkSmartPointer<vtkMultiBlockDataSet> mesh;
+  //vtkSmartPointer<vtkMultiBlockDataSet> mesh;
   vtkSmartPointer<vtkDoubleArray> AtomPositions;
   vtkSmartPointer<vtkIntArray> AtomTypes;
   vtkSmartPointer<vtkIntArray> AtomIDs;
@@ -287,6 +287,8 @@ int lammpsDataAdaptor::GetNumberOfMeshes(unsigned int &numMeshes)
 int lammpsDataAdaptor::GetMesh(const std::string &meshName, bool structureOnly,
     vtkDataObject *&mesh)
 {
+  mesh = nullptr;
+
   if (meshName != "atoms") 
     {  
     SENSEI_ERROR("No mesh \"" << meshName << "\"")
@@ -295,32 +297,29 @@ int lammpsDataAdaptor::GetMesh(const std::string &meshName, bool structureOnly,
 
   DInternals& internals = (*this->Internals);
 
-  if (!internals.mesh)
-    {
-    vtkSmartPointer<vtkPolyData> pd = vtkSmartPointer<vtkPolyData>::New();
+  vtkMultiBlockDataSet *mb = vtkMultiBlockDataSet::New();
+  vtkSmartPointer<vtkPolyData> pd = vtkSmartPointer<vtkPolyData>::New();
 
-    if(!structureOnly)
-      {
-      vtkSmartPointer<vtkPoints> pts = vtkSmartPointer<vtkPoints>::New();
-      pts->SetData(internals.AtomPositions);
-      pd->SetPoints(pts);
-      }
+  if(!structureOnly)
+    {
+    vtkSmartPointer<vtkPoints> pts = vtkSmartPointer<vtkPoints>::New();
+    pts->SetData(internals.AtomPositions);
+    pd->SetPoints(pts);
+    }
 
-    pd->SetVerts( internals.vertices );
+  pd->SetVerts( internals.vertices );
 
-    int rank, size; 
-    MPI_Comm comm;
+  int rank, size; 
+  MPI_Comm comm;
   	
-    comm = GetCommunicator();
-    MPI_Comm_rank(comm, &rank);
-    MPI_Comm_size(comm, &size); 
+  comm = GetCommunicator();
+  MPI_Comm_rank(comm, &rank);
+  MPI_Comm_size(comm, &size); 
 
-    internals.mesh = vtkSmartPointer<vtkMultiBlockDataSet>::New();
-    internals.mesh->SetNumberOfBlocks(size);
-    internals.mesh->SetBlock(rank, pd);
-    }
+  mb->SetNumberOfBlocks(size);
+  mb->SetBlock(rank, pd);
 
-  mesh = internals.mesh;
+  mesh = mb;
 
   return 0;
 }
@@ -329,6 +328,7 @@ int lammpsDataAdaptor::GetMesh(const std::string &meshName, bool structureOnly,
 int lammpsDataAdaptor::AddArray(vtkDataObject* mesh, const std::string &meshName,
     int association, const std::string &arrayName)
 {
+
   if (meshName != "atoms")
     {
     SENSEI_ERROR("No mesh \"" << meshName << "\"")
@@ -343,8 +343,8 @@ int lammpsDataAdaptor::AddArray(vtkDataObject* mesh, const std::string &meshName
 
   if (arrayName == "type")
     {  	
-    DInternals& internals = (*this->Internals);
-    vtkMultiBlockDataSet* md = vtkMultiBlockDataSet::SafeDownCast(mesh);
+    //DInternals& internals = (*this->Internals);
+    vtkMultiBlockDataSet* md = dynamic_cast<vtkMultiBlockDataSet*>(mesh);
     vtkSmartPointer<vtkPolyData> pd = vtkSmartPointer<vtkPolyData>::New();
   
     int rank;
@@ -354,7 +354,7 @@ int lammpsDataAdaptor::AddArray(vtkDataObject* mesh, const std::string &meshName
     MPI_Comm_rank(comm, &rank);
 
     pd = vtkPolyData::SafeDownCast(md->GetBlock(rank));
-    pd->GetPointData()->AddArray(internals.AtomTypes);
+    pd->GetPointData()->AddArray(this->Internals->AtomTypes);
     }
 
   if (arrayName == "id")
@@ -483,14 +483,14 @@ int lammpsDataAdaptor::GetMeshMetadata(unsigned int id, sensei::MeshMetadataPtr
 //-----------------------------------------------------------------------------
 int lammpsDataAdaptor::ReleaseData()
 {
-  DInternals& internals = (*this->Internals);
-
-  internals.mesh = NULL;
-  internals.AtomPositions = NULL;
-  internals.AtomTypes = NULL;
-  internals.AtomIDs = NULL;
-  internals.nlocal = 0;
-  internals.nghost = 0;
+  //DInternals& internals = (*this->Internals);
+
+  //internals.mesh = NULL;
+  //internals.AtomPositions = NULL;
+  //internals.AtomTypes = NULL;
+  //internals.AtomIDs = NULL;
+  //internals.nlocal = 0;
+  //internals.nghost = 0;
 
   return 0;
 }

From 6a07414b83288a8113321a90de4f880775011299 Mon Sep 17 00:00:00 2001
From: Silvio Rizzi <srizzi@anl.gov>
Date: Mon, 9 Aug 2021 03:00:08 +0000
Subject: [PATCH 06/11] added directory with config files for easier debugging
 on other systems

---
 miniapps/lammps/configs/adiostransport.xml    | 17 ++++++++
 miniapps/lammps/configs/adioswrite.xml        | 23 +++++++++++
 miniapps/lammps/configs/adioswriteosc.xml     | 23 +++++++++++
 miniapps/lammps/configs/cmakeconfig_adios.sh  | 11 +++++
 miniapps/lammps/configs/cmakeconfig_sensei.sh | 17 ++++++++
 miniapps/lammps/configs/cmakeconfig_vtk8.sh   | 17 ++++++++
 miniapps/lammps/configs/config.xml            | 34 ++++++++++++++++
 miniapps/lammps/configs/endpoint.sh           |  6 +++
 miniapps/lammps/configs/in.rhodo              | 28 +++++++++++++
 miniapps/lammps/configs/run.sh                |  7 ++++
 miniapps/lammps/configs/sample.osc            |  6 +++
 miniapps/lammps/configs/vtk.xml               |  9 +++++
 miniapps/lammps/configs/vtkosc.xml            |  9 +++++
 miniapps/lammps/lammpsDataAdaptor.cxx         | 40 +++++++++++++++----
 14 files changed, 240 insertions(+), 7 deletions(-)
 create mode 100644 miniapps/lammps/configs/adiostransport.xml
 create mode 100644 miniapps/lammps/configs/adioswrite.xml
 create mode 100644 miniapps/lammps/configs/adioswriteosc.xml
 create mode 100644 miniapps/lammps/configs/cmakeconfig_adios.sh
 create mode 100644 miniapps/lammps/configs/cmakeconfig_sensei.sh
 create mode 100644 miniapps/lammps/configs/cmakeconfig_vtk8.sh
 create mode 100644 miniapps/lammps/configs/config.xml
 create mode 100755 miniapps/lammps/configs/endpoint.sh
 create mode 100644 miniapps/lammps/configs/in.rhodo
 create mode 100755 miniapps/lammps/configs/run.sh
 create mode 100644 miniapps/lammps/configs/sample.osc
 create mode 100644 miniapps/lammps/configs/vtk.xml
 create mode 100644 miniapps/lammps/configs/vtkosc.xml

diff --git a/miniapps/lammps/configs/adiostransport.xml b/miniapps/lammps/configs/adiostransport.xml
new file mode 100644
index 000000000..06f7360c9
--- /dev/null
+++ b/miniapps/lammps/configs/adiostransport.xml
@@ -0,0 +1,17 @@
+<sensei>
+  <transport type="adios2" filename="test_%05d.bp" debug_mode="1" engine="bp4">
+    <engine_parameters>
+        NumAggregators = 0
+        InitialBufferSize = 1Mb
+        BufferGrowthFactor = 2
+        MaxBufferSize = 1024Mb
+        StatsLevel = 0
+        Profile = Off
+      </engine_parameters>
+
+      <!-- subset by mesh and array -->
+      <mesh name="atoms">
+        <point_arrays> type,id </point_arrays>
+      </mesh>
+    </transport>
+</sensei>
diff --git a/miniapps/lammps/configs/adioswrite.xml b/miniapps/lammps/configs/adioswrite.xml
new file mode 100644
index 000000000..5254f425e
--- /dev/null
+++ b/miniapps/lammps/configs/adioswrite.xml
@@ -0,0 +1,23 @@
+<sensei>
+
+  <!-- configure ADIOS2 write -->
+  <analysis type="adios2" filename="test_%05d.bp" engine="BP4"
+    debug_mode="1" enabled="1" steps_per_file="1" >
+
+    <!-- ADIOS2 engine parameters -->
+    <engine_parameters>
+      NumAggregators = 0
+      InitialBufferSize = 1Mb
+      BufferGrowthFactor = 2
+      MaxBufferSize = 1024Mb
+      StatsLevel = 0
+      Profile = Off
+    </engine_parameters>
+
+    <!-- subset by mesh and array -->
+    <mesh name="atoms">
+      <point_arrays>type,id </point_arrays>
+    </mesh>
+
+  </analysis>
+</sensei>
diff --git a/miniapps/lammps/configs/adioswriteosc.xml b/miniapps/lammps/configs/adioswriteosc.xml
new file mode 100644
index 000000000..93f89e05d
--- /dev/null
+++ b/miniapps/lammps/configs/adioswriteosc.xml
@@ -0,0 +1,23 @@
+<sensei>
+
+  <!-- configure ADIOS2 write -->
+  <analysis type="adios2" filename="test_%05d.bp" engine="BP4"
+    debug_mode="1" enabled="1" steps_per_file="1" >
+
+    <!-- ADIOS2 engine parameters -->
+    <engine_parameters>
+      NumAggregators = 0
+      InitialBufferSize = 1Mb
+      BufferGrowthFactor = 2
+      MaxBufferSize = 1024Mb
+      StatsLevel = 0
+      Profile = Off
+    </engine_parameters>
+
+    <!-- subset by mesh and array -->
+    <mesh name="mesh">
+      <cell_arrays>data</cell_arrays>
+    </mesh>
+
+  </analysis>
+</sensei>
diff --git a/miniapps/lammps/configs/cmakeconfig_adios.sh b/miniapps/lammps/configs/cmakeconfig_adios.sh
new file mode 100644
index 000000000..8bb7cb50f
--- /dev/null
+++ b/miniapps/lammps/configs/cmakeconfig_adios.sh
@@ -0,0 +1,11 @@
+#!/bin/bash
+
+export CC=`which gcc`
+export CXX=`which g++`
+
+
+cmake -S ./src/ADIOS2 -B ./build/adios \
+ -DCMAKE_INSTALL_PREFIX=`pwd`/install/adios \
+ -DADIOS2_USE_Fortran=OFF \
+ -DADIOS2_BUILD_EXAMPLES=OFF \
+ -DCMAKE_BUILD_TYPE=Debug
diff --git a/miniapps/lammps/configs/cmakeconfig_sensei.sh b/miniapps/lammps/configs/cmakeconfig_sensei.sh
new file mode 100644
index 000000000..fea20db34
--- /dev/null
+++ b/miniapps/lammps/configs/cmakeconfig_sensei.sh
@@ -0,0 +1,17 @@
+
+export CC=`which gcc`
+export CXX=`which g++`
+export LAMMPS_DIR=`pwd`/src/lammps/
+
+
+
+cmake -S ./src/SENSEI -B ./build/sensei \
+  -DCMAKE_INSTALL_PREFIX=`pwd`/install/sensei \
+  -DVTK_DIR=`pwd`/install/vtk/lib64/cmake/vtk-8.2 \
+  -DENABLE_VTK_IO=ON \
+  -DLAMMPS_DIR=$LAMMPS_DIR \
+  -DENABLE_LAMMPS=ON \
+  -DENABLE_MANDELBROT=OFF \
+  -DENABLE_OSCILLATORS=OFF \
+  -DENABLE_ADIOS2=ON \
+  -DADIOS2_DIR=`pwd`/install/adios/lib64/cmake/adios2
diff --git a/miniapps/lammps/configs/cmakeconfig_vtk8.sh b/miniapps/lammps/configs/cmakeconfig_vtk8.sh
new file mode 100644
index 000000000..13f76fdec
--- /dev/null
+++ b/miniapps/lammps/configs/cmakeconfig_vtk8.sh
@@ -0,0 +1,17 @@
+
+export CC=`which gcc`
+export CXX=`which g++`
+
+cmake -S src/vtk -B build/vtk \
+ -DVTK_Group_Imaging=OFF \
+ -DVTK_Group_MPI=OFF \
+ -DVTK_Group_Qt=OFF \
+ -DVTK_Group_Rendering=OFF \
+ -DVTK_Group_StandAlone=ON \
+ -DVTK_Group_Tk=OFF \
+ -DVTK_Group_Views=OFF \
+ -DVTK_Group_Web=OFF \
+ -DVTK_RENDERING_BACKEND=None \
+ -DCMAKE_BUILD_TYPE=Debug \
+ -DBUILD_TESTING=OFF \
+ -DCMAKE_INSTALL_PREFIX=`pwd`/install/vtk
diff --git a/miniapps/lammps/configs/config.xml b/miniapps/lammps/configs/config.xml
new file mode 100644
index 000000000..10c06d046
--- /dev/null
+++ b/miniapps/lammps/configs/config.xml
@@ -0,0 +1,34 @@
+<sensei>
+  <analysis type="PosthocIO"
+    output_dir="./vtkfiles" file_name="output" mode="paraview"
+    enabled="0">
+    <mesh name="atoms" structure_only="0" ghost_cells="1" ghost_nodes="0">
+      <point_arrays>type,id</point_arrays>
+    </mesh>
+  </analysis>
+
+   <!-- configure ADIOS2 write -->
+   <analysis type="adios2" filename="test_%05d.bp”
+      engine="BP4" enabled="1”
+      steps_per_file="1" >
+
+      <!-- ADIOS2 engine parameters -->
+      <engine_parameters>
+          NumAggregators = 0
+          InitialBufferSize = 1Mb
+          BufferGrowthFactor = 2
+          MaxBufferSize = 1024Mb
+          StatsLevel = 0
+          Profile = Off
+      </engine_parameters>
+
+      <!-- subset by mesh and array -->
+      <mesh name="atoms">
+         <point_arrays>type,id </point_arrays>
+      </mesh>
+
+  </analysis>
+
+
+
+</sensei>
diff --git a/miniapps/lammps/configs/endpoint.sh b/miniapps/lammps/configs/endpoint.sh
new file mode 100755
index 000000000..7238f9ca2
--- /dev/null
+++ b/miniapps/lammps/configs/endpoint.sh
@@ -0,0 +1,6 @@
+
+export LD_LIBRARY_PATH=/lus/grand/projects/visualization/srizzi/BUILDS/sensei-lammps/install/vtk/lib64:/lus/grand/projects/visualization/srizzi/BUILDS/sensei-lammps/install/adios/lib64:$LD_LIBRARY_PATH
+
+#mpiexec -n 4 ./lammpsDriver in.rhodo -n 5 -sensei config.xml
+#mpiexec -n 4 ./lammpsDriver in.rhodo -n 5 -sensei adioswrite.xml
+./SENSEIEndPoint -a vtkosc.xml -t adiostransport.xml
diff --git a/miniapps/lammps/configs/in.rhodo b/miniapps/lammps/configs/in.rhodo
new file mode 100644
index 000000000..d4ec87789
--- /dev/null
+++ b/miniapps/lammps/configs/in.rhodo
@@ -0,0 +1,28 @@
+# Rhodopsin model
+
+units           real  
+neigh_modify    delay 5 every 1   
+
+atom_style      full  
+bond_style      harmonic 
+angle_style     charmm 
+dihedral_style  charmm 
+improper_style  harmonic 
+pair_style      lj/charmm/coul/long 8.0 10.0 
+pair_modify     mix arithmetic 
+kspace_style    pppm 1e-4 
+
+read_data       data.rhodo
+
+fix             1 all shake 0.0001 5 0 m 1.0 a 232
+fix             2 all npt temp 300.0 300.0 100.0 z 0.0 0.0 1000.0 mtk no pchain 0 tchain 1
+fix		3 all external pf/callback 1 10000
+
+special_bonds   charmm
+ 
+thermo          50
+thermo_style    multi 
+timestep        2.0
+
+#run		100
+
diff --git a/miniapps/lammps/configs/run.sh b/miniapps/lammps/configs/run.sh
new file mode 100755
index 000000000..84ea73dbd
--- /dev/null
+++ b/miniapps/lammps/configs/run.sh
@@ -0,0 +1,7 @@
+
+export LD_LIBRARY_PATH=/lus/grand/projects/visualization/srizzi/BUILDS/sensei-lammps/install/vtk/lib64:/lus/grand/projects/visualization/srizzi/BUILDS/sensei-lammps/install/adios/lib64:$LD_LIBRARY_PATH
+
+#mpiexec -n 4 ./lammpsDriver in.rhodo -n 5 -sensei config.xml
+#mpiexec -n 4 ./lammpsDriver in.rhodo -n 5 -sensei adioswrite.xml
+#./lammpsDriver in.rhodo -n 5 -sensei adioswrite.xml
+./oscillator --t-end 0.05 -f adioswriteosc.xml sample.osc
diff --git a/miniapps/lammps/configs/sample.osc b/miniapps/lammps/configs/sample.osc
new file mode 100644
index 000000000..42650fc26
--- /dev/null
+++ b/miniapps/lammps/configs/sample.osc
@@ -0,0 +1,6 @@
+# type      center      r       omega0      zeta
+damped      32 32 32    10.     3.14        .3
+damped      16 32 16    10.     9.5         .1
+damped      48 32 48    5.      3.14        .1
+decaying    16 32 48    15      3.14
+periodic    48 32 16    15      3.14
diff --git a/miniapps/lammps/configs/vtk.xml b/miniapps/lammps/configs/vtk.xml
new file mode 100644
index 000000000..42940f92a
--- /dev/null
+++ b/miniapps/lammps/configs/vtk.xml
@@ -0,0 +1,9 @@
+<sensei>
+  <analysis type="PosthocIO"
+    output_dir="./vtkfiles" file_name="output" mode="paraview"
+    enabled="1">
+    <mesh name="atoms" structure_only="0" ghost_cells="1" ghost_nodes="0">
+      <point_arrays>type,id</point_arrays>
+    </mesh>
+  </analysis>
+</sensei>
diff --git a/miniapps/lammps/configs/vtkosc.xml b/miniapps/lammps/configs/vtkosc.xml
new file mode 100644
index 000000000..ea96e8999
--- /dev/null
+++ b/miniapps/lammps/configs/vtkosc.xml
@@ -0,0 +1,9 @@
+<sensei>
+  <analysis type="PosthocIO"
+    output_dir="./vtkfiles" file_name="output" mode="paraview"
+    enabled="1">
+    <mesh name="mesh" structure_only="0" ghost_cells="1" ghost_nodes="0">
+      <cell_arrays>data</cell_arrays>
+    </mesh>
+  </analysis>
+</sensei>
diff --git a/miniapps/lammps/lammpsDataAdaptor.cxx b/miniapps/lammps/lammpsDataAdaptor.cxx
index 2dd3ec119..777091638 100644
--- a/miniapps/lammps/lammpsDataAdaptor.cxx
+++ b/miniapps/lammps/lammpsDataAdaptor.cxx
@@ -462,18 +462,44 @@ int lammpsDataAdaptor::GetMeshMetadata(unsigned int id, sensei::MeshMetadataPtr
     {
     SENSEI_WARNING("lammps data adaptor. Flags.BlockArrayRangeSet()")
  
-    unsigned long nvals;   
-    std::array<int,2> typeBlockRange = getArrayRange(nvals, this->Internals->type);
-    std::array<int,2> idBlockRange = getArrayRange(nvals, this->Internals->id);
+    //unsigned long nvals;   
+    //std::array<int,2> typeBlockRange = getArrayRange(nvals, this->Internals->type);
+    //std::array<int,2> idBlockRange = getArrayRange(nvals, this->Internals->id);
     // fixme
     //metadata->BlockArrayRange.push_back({typeBlockRange, idBlockRange});
 
-    std::array<int,2> typeRange = { this->Internals->typeRange.min[0], this->Internals->typeRange.max[0] };
-    std::array<int,2> idRange = { this->Internals->idRange.min[0], this->Internals->idRange.max[0] };
+    //std::array<int,2> typeRange = { this->Internals->typeRange.min[0], this->Internals->typeRange.max[0] };
+    //std::array<int,2> idRange = { this->Internals->idRange.min[0], this->Internals->idRange.max[0] };
     // fixme
     //metadata->ArrayRange.push_back(typeRange); 
-    //metadata->ArrayRange.push_back(idRange); 
-}
+    //metadata->ArrayRange.push_back(idRange);
+    
+    //float gmin = std::numeric_limits<float>::lowest();
+    //float gmax = std::numeric_limits<float>::max();
+    
+
+    //std::vector<std::array<double,2>> blkRange{{bmin,bmax}};
+    
+    //metadata->BlockArrayRange.push_back(blkRange);
+    //metadata->BlockArrayRange.push_back(blkRange);
+    //metadata->ArrayRange.push_back({gmin, gmax});
+    //metadata->ArrayRange.push_back({gmin, gmax});
+
+
+    float bmin = std::numeric_limits<float>::lowest();
+    float bmax = std::numeric_limits<float>::max();
+
+    std::vector<std::array<double,2>> typeBlockRange{{bmin,bmax}};
+    std::vector<std::array<double,2>> idBlockRange ={{bmin,bmax}};
+
+    metadata->BlockArrayRange.push_back(typeBlockRange);
+    metadata->BlockArrayRange.push_back(idBlockRange);
+
+    metadata->ArrayRange.push_back({bmin,bmax}); 
+ 
+
+
+    }
 
   return 0;
 }

From bc3adedf690dcc6fdd5e38f87012f80a7ab12d35 Mon Sep 17 00:00:00 2001
From: Silvio Rizzi <srizzi@anl.gov>
Date: Mon, 9 Aug 2021 21:06:51 +0000
Subject: [PATCH 07/11] debugging and warning clean up. Added /fixme tags where
 some attention may be needed

---
 miniapps/lammps/lammpsBridge.cxx      |  8 +--
 miniapps/lammps/lammpsDataAdaptor.cxx | 99 +++++++++------------------
 miniapps/lammps/lammpsDataAdaptor.h   | 15 +---
 miniapps/lammps/lammpsDriver.cxx      | 35 ++++------
 4 files changed, 49 insertions(+), 108 deletions(-)

diff --git a/miniapps/lammps/lammpsBridge.cxx b/miniapps/lammps/lammpsBridge.cxx
index 15bf094e9..18c5857e5 100644
--- a/miniapps/lammps/lammpsBridge.cxx
+++ b/miniapps/lammps/lammpsBridge.cxx
@@ -2,7 +2,6 @@
 #include "lammpsDataAdaptor.h"
 #include <vtkSmartPointer.h>
 #include <ConfigurableAnalysis.h>
-
 #include <vtkDataObject.h>
 
 namespace lammpsBridge
@@ -36,11 +35,6 @@ void SetData(long ntimestep, int nlocal, int *id,
 void Analyze()
 {
   AnalysisAdaptor->Execute(DataAdaptor.GetPointer());
-
-  //vtkDataObject *mesh;
-  //DataAdaptor->GetMesh("atoms", false, mesh);
-  //DataAdaptor->AddArray(mesh, "atoms", vtkDataObject::FIELD_ASSOCIATION_POINTS, "type");
-
   DataAdaptor->ReleaseData();
 }
 
@@ -52,4 +46,4 @@ void Finalize()
   DataAdaptor = nullptr;
 }
 
-}	// namespace lammpsBridge
+} // namespace lammpsBridge
diff --git a/miniapps/lammps/lammpsDataAdaptor.cxx b/miniapps/lammps/lammpsDataAdaptor.cxx
index 777091638..3b42679b4 100644
--- a/miniapps/lammps/lammpsDataAdaptor.cxx
+++ b/miniapps/lammps/lammpsDataAdaptor.cxx
@@ -10,26 +10,33 @@
 #include <vtkCellArray.h>
 #include <sdiy/master.hpp>
 
+// fixme - unused
+/***
 static
-std::array<int,2> getArrayRange(unsigned long nSize, int *x) {
+std::array<int,2> getArrayRange(unsigned long nSize, int *x) 
+{
   int xmin = std::numeric_limits<int>::max(); 
   int xmax = std::numeric_limits<int>::lowest();
-  for(unsigned int i=0; i<nSize; ++i) {
+  for(unsigned int i=0; i<nSize; ++i) 
+    {
     xmin = std::min(xmin, x[i]);
     xmax = std::max(xmax, x[i]);
-  }
+    }
 
   return {xmin, xmax};
 }
+***/
 
 static
-std::array<double,2> getArrayRange(unsigned long nSize,double *x) {
+std::array<double,2> getArrayRange(unsigned long nSize,double *x) 
+{
   double xmin = std::numeric_limits<double>::max(); 
   double xmax = std::numeric_limits<double>::lowest();
-  for(unsigned int i=0; i<nSize; ++i) {
+  for(unsigned int i=0; i<nSize; ++i) 
+    {
     xmin = std::min(xmin, x[i]);
     xmax = std::max(xmax, x[i]);
-  }
+    }
 
   return {xmin, xmax};
 }
@@ -50,7 +57,6 @@ namespace senseiLammps
 
 struct lammpsDataAdaptor::DInternals
 {
-  //vtkSmartPointer<vtkMultiBlockDataSet> mesh;
   vtkSmartPointer<vtkDoubleArray> AtomPositions;
   vtkSmartPointer<vtkIntArray> AtomTypes;
   vtkSmartPointer<vtkIntArray> AtomIDs;
@@ -118,6 +124,7 @@ void lammpsDataAdaptor::AddLAMMPSData( long ntimestep, int nlocal, int *id,
   // atom coordinates
   if (internals.AtomPositions)
     {
+    // fixme - ghost arrays
     //long nvals = nlocal + nghost;
     long nvals = nlocal;
 
@@ -135,6 +142,7 @@ void lammpsDataAdaptor::AddLAMMPSData( long ntimestep, int nlocal, int *id,
   // atom types
   if (internals.AtomTypes)
     {
+    // fixme - ghost arrays
     //long nvals = nlocal + nghost;
     long nvals = nlocal;
 
@@ -152,6 +160,7 @@ void lammpsDataAdaptor::AddLAMMPSData( long ntimestep, int nlocal, int *id,
   // atom IDs
   if (internals.AtomIDs)
     {
+    // fixme - ghost arryas
     //long nvals = nlocal + nghost;
     long nvals = nlocal;
 
@@ -171,6 +180,7 @@ void lammpsDataAdaptor::AddLAMMPSData( long ntimestep, int nlocal, int *id,
     {
     vtkIdType pid[1] = {0};
 
+    // fixme - ghost arrays
     //for( int i=0; i < nlocal+nghost; i++) {
     for( int i=0; i < nlocal; i++) 
       {
@@ -194,7 +204,7 @@ void lammpsDataAdaptor::AddLAMMPSData( long ntimestep, int nlocal, int *id,
     y_range[0], y_range[1],
     z_range[0], z_range[1]);
 
-  /// XXX Set type and id range
+  // Set type and id range
   this->Internals->typeRange.min[0] = std::numeric_limits<int>::max();
   this->Internals->typeRange.max[0] = std::numeric_limits<int>::min();
   this->Internals->idRange.min[0] = std::numeric_limits<int>::max();
@@ -202,7 +212,6 @@ void lammpsDataAdaptor::AddLAMMPSData( long ntimestep, int nlocal, int *id,
 
   // timestep
   this->SetDataTimeStep(ntimestep);
-
 }
 
 void lammpsDataAdaptor::SetBlockBounds(double xmin, double xmax,
@@ -343,7 +352,6 @@ int lammpsDataAdaptor::AddArray(vtkDataObject* mesh, const std::string &meshName
 
   if (arrayName == "type")
     {  	
-    //DInternals& internals = (*this->Internals);
     vtkMultiBlockDataSet* md = dynamic_cast<vtkMultiBlockDataSet*>(mesh);
     vtkSmartPointer<vtkPolyData> pd = vtkSmartPointer<vtkPolyData>::New();
   
@@ -376,6 +384,7 @@ int lammpsDataAdaptor::AddArray(vtkDataObject* mesh, const std::string &meshName
   return 0;  
 }
 
+// fixme
 // not implemented
 //----------------------------------------------------------------------------
 int lammpsDataAdaptor::AddGhostCellsArray(vtkDataObject *mesh, const std::string &meshName)
@@ -403,14 +412,13 @@ int lammpsDataAdaptor::GetMeshMetadata(unsigned int id, sensei::MeshMetadataPtr
   MPI_Comm_size(comm, &nRanks);
 
 
-  int nBlocks = 1; // One block per rank
   metadata->MeshName = "atoms";
   metadata->MeshType = VTK_MULTIBLOCK_DATA_SET;
   metadata->BlockType = VTK_POLY_DATA;
   metadata->CoordinateType = VTK_DOUBLE;
   metadata->NumBlocks = nRanks;
-  metadata->NumBlocksLocal = {nBlocks};
-  metadata->NumGhostCells = this->Internals->nghost;
+  metadata->NumBlocksLocal = {1};     // One block per rank
+  metadata->NumGhostCells = 0;
   metadata->NumArrays = 2;
   metadata->ArrayName = {"type", "id"};
   metadata->ArrayCentering = {vtkDataObject::POINT, vtkDataObject::POINT};
@@ -420,23 +428,25 @@ int lammpsDataAdaptor::GetMeshMetadata(unsigned int id, sensei::MeshMetadataPtr
 
   if (metadata->Flags.BlockExtentsSet())
     {
-    SENSEI_WARNING("lammps data adaptor. Flags.BlockExtentsSet()")
+    //SENSEI_WARNING("lammps data adaptor. Flags.BlockExtentsSet()")
+    
+    // fixme
     // There should be no extent for a PolyData, but ADIOS2 needs this
     std::array<int,6> ext = { 0, 0, 0, 0, 0, 0};
     metadata->Extent = std::move(ext);
 
-    metadata->BlockExtents.reserve(nBlocks);
+    metadata->BlockExtents.reserve(1);	// One block per rank
     metadata->BlockExtents.emplace_back(std::move(ext));
     }
 
   if (metadata->Flags.BlockBoundsSet())
     {
-    SENSEI_WARNING("lammps data adaptor. Flags.BlockBoundsSet()")
+    //SENSEI_WARNING("lammps data adaptor. Flags.BlockBoundsSet()")
     std::array<double, 6> bounds;
     getBounds(this->Internals->DomainBounds, bounds.data());
     metadata->Bounds = std::move(bounds);
     
-    metadata->BlockBounds.reserve(nBlocks);
+    metadata->BlockBounds.reserve(1);  // One block per rank
     
     getBounds(this->Internals->BlockBounds, bounds.data());
     metadata->BlockBounds.emplace_back(std::move(bounds));
@@ -449,7 +459,7 @@ int lammpsDataAdaptor::GetMeshMetadata(unsigned int id, sensei::MeshMetadataPtr
     SENSEI_WARNING("lammps data adaptor. Flags.BlockSizeSet()")
     metadata->BlockNumCells.push_back(nCells);
     metadata->BlockNumPoints.push_back(nCells);
-    metadata->BlockCellArraySize.push_back(2 * nCells); // XXX- VTK_POINTS
+    metadata->BlockCellArraySize.push_back(2 * nCells); // fixme - VTK_POINTS
     }
 
   if (metadata->Flags.BlockDecompSet())
@@ -460,64 +470,23 @@ int lammpsDataAdaptor::GetMeshMetadata(unsigned int id, sensei::MeshMetadataPtr
 
   if (metadata->Flags.BlockArrayRangeSet())
     {
-    SENSEI_WARNING("lammps data adaptor. Flags.BlockArrayRangeSet()")
+    //SENSEI_WARNING("lammps data adaptor. Flags.BlockArrayRangeSet()")
  
-    //unsigned long nvals;   
-    //std::array<int,2> typeBlockRange = getArrayRange(nvals, this->Internals->type);
-    //std::array<int,2> idBlockRange = getArrayRange(nvals, this->Internals->id);
-    // fixme
-    //metadata->BlockArrayRange.push_back({typeBlockRange, idBlockRange});
-
-    //std::array<int,2> typeRange = { this->Internals->typeRange.min[0], this->Internals->typeRange.max[0] };
-    //std::array<int,2> idRange = { this->Internals->idRange.min[0], this->Internals->idRange.max[0] };
-    // fixme
-    //metadata->ArrayRange.push_back(typeRange); 
-    //metadata->ArrayRange.push_back(idRange);
-    
-    //float gmin = std::numeric_limits<float>::lowest();
-    //float gmax = std::numeric_limits<float>::max();
-    
-
-    //std::vector<std::array<double,2>> blkRange{{bmin,bmax}};
-    
-    //metadata->BlockArrayRange.push_back(blkRange);
-    //metadata->BlockArrayRange.push_back(blkRange);
-    //metadata->ArrayRange.push_back({gmin, gmax});
-    //metadata->ArrayRange.push_back({gmin, gmax});
-
-
-    float bmin = std::numeric_limits<float>::lowest();
-    float bmax = std::numeric_limits<float>::max();
-
-    std::vector<std::array<double,2>> typeBlockRange{{bmin,bmax}};
-    std::vector<std::array<double,2>> idBlockRange ={{bmin,bmax}};
-
-    metadata->BlockArrayRange.push_back(typeBlockRange);
-    metadata->BlockArrayRange.push_back(idBlockRange);
+    // fixme - forced values
+    double bmin = (double) std::numeric_limits<int>::lowest();
+    double bmax = (double) std::numeric_limits<int>::max();
 
+    std::vector<std::array<double,2>> typeBlockRange ={{bmin,bmax} , {bmin,bmax}};
+    metadata->BlockArrayRange.push_back({typeBlockRange});
     metadata->ArrayRange.push_back({bmin,bmax}); 
- 
-
-
     }
 
   return 0;
 }
 
-
-
 //-----------------------------------------------------------------------------
 int lammpsDataAdaptor::ReleaseData()
 {
-  //DInternals& internals = (*this->Internals);
-
-  //internals.mesh = NULL;
-  //internals.AtomPositions = NULL;
-  //internals.AtomTypes = NULL;
-  //internals.AtomIDs = NULL;
-  //internals.nlocal = 0;
-  //internals.nghost = 0;
-
   return 0;
 }
 
diff --git a/miniapps/lammps/lammpsDataAdaptor.h b/miniapps/lammps/lammpsDataAdaptor.h
index dd214104b..11662ddff 100644
--- a/miniapps/lammps/lammpsDataAdaptor.h
+++ b/miniapps/lammps/lammpsDataAdaptor.h
@@ -13,9 +13,6 @@ class lammpsDataAdaptor : public sensei::DataAdaptor
   static lammpsDataAdaptor* New();
   senseiTypeMacro(lammpsDataAdaptor, sensei::DataAdaptor);
 
-  /// @brief Initialize the data adaptor.
-  ///
-  /// This initializes the data adaptor. This must be called once per simulation run.
   void Initialize();
 
   void AddLAMMPSData( long ntimestep, int nlocal, int *id, 
@@ -44,27 +41,17 @@ class lammpsDataAdaptor : public sensei::DataAdaptor
   void SetBlockBounds(double xmin, double xmax,
     double ymin, double ymax, double zmin, double zmax);
 
-// SENSEI API
+  // SENSEI API
   int GetNumberOfMeshes(unsigned int &numMeshes) override;
 
   int GetMeshMetadata(unsigned int id, sensei::MeshMetadataPtr &md) override;
 
-  //int GetMeshName(unsigned int id, std::string &meshName) override;
-
   int GetMesh(const std::string &meshName, bool structureOnly,
     vtkDataObject *&mesh) override;
 
   int AddArray(vtkDataObject* mesh, const std::string &meshName,
     int association, const std::string &arrayName) override;
 
-  //int GetNumberOfArrays(const std::string &meshName, int association,
-  //  unsigned int &numberOfArrays) override;
-
-  //int GetArrayName(const std::string &meshName, int association,
-  //  unsigned int index, std::string &arrayName) override;
-
-  //int GetMeshHasGhostCells(const std::string &meshName, int &nLayers) override;
-
   int AddGhostCellsArray(vtkDataObject* mesh, const std::string &meshName) override;
 
   int ReleaseData() override;
diff --git a/miniapps/lammps/lammpsDriver.cxx b/miniapps/lammps/lammpsDriver.cxx
index 90dc8860e..f78f7e1c6 100644
--- a/miniapps/lammps/lammpsDriver.cxx
+++ b/miniapps/lammps/lammpsDriver.cxx
@@ -1,4 +1,14 @@
 /* ----------------------------------------------------------------------
+   This file contains code from Will Usher ( https://www.willusher.io ) 
+   used originally in this publication:
+   Usher, Will, Silvio Rizzi, Ingo Wald, Jefferson Amstutz, Joseph Insley, 
+   Venkatram Vishwanath, Nicola Ferrier, Michael E. Papka, and Valerio Pascucci. 
+   "libIS: a lightweight library for flexible in transit visualization." 
+   In Proceedings of the Workshop on In Situ Infrastructures for Enabling 
+   Extreme-Scale Analysis and Visualization, pp. 33-38. 2018.
+
+   This code also borrowed from an example file in the LAMMPS source code.   
+
    LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
    www.cs.sandia.gov/~sjplimp/lammps.html
    Steve Plimpton, sjplimp@sandia.gov, Sandia National Laboratories
@@ -30,7 +40,8 @@
 
 using namespace LAMMPS_NS;
 
-struct Info {
+struct Info 
+{
   int me;
   LAMMPS *lmp;
 } globalInfo;
@@ -63,15 +74,7 @@ void insituCallback(void *ptr, bigint ntimestep, int nlocal,
   lammpsBridge::SetData(ntimestep, nlocal, id, *nghost, type, all_pos, 
                         xsublo, xsubhi, ysublo, ysubhi, zsublo, zsubhi);
 
-  if (0 == globalInfo.me)
-    std::cout << "###### SENSEI instrumentation: bridge analyze() ######" << std::endl;    
   lammpsBridge::Analyze();
-
-  if (0 == globalInfo.me)
-    std::cout << "###### SENSEI instrumentation: after bridge analyze() ######" << std::endl;    
- 
-
-
 }
 
 const static std::string USAGE =
@@ -82,9 +85,7 @@ const static std::string USAGE =
 "  -n <int>            Specify the number of steps to simulate for. Default is 10000\n"
 "  -lmp <lammps args>  Pass the list of arguments <args> to LAMMPS as if they were\n"
 "                      command line args to LAMMPS. This must be the last argument, all\n"
-"                      following arguments will be passed to lammps.\n"
-"  -log                Generate full time and memory usage log.\n"
-"  -shortlog           Generate a summary time and memory usage log.";
+"                      following arguments will be passed to lammps.\n";
 
 int main(int argc, char **argv) 
 {
@@ -99,8 +100,6 @@ int main(int argc, char **argv)
 	
   // Additional args the user may be passing to lammps
   std::vector<char*> lammps_args(1, argv[0]);
-  bool log = false;
-  bool shortlog = false;
   std::string sensei_xml;
   size_t sim_steps = 10000;
   for (size_t i = 2; i < args.size(); ++i) 
@@ -114,10 +113,6 @@ int main(int argc, char **argv)
       std::cout << USAGE << "\n";
       return 0;
       } 
-    else if (args[i] == "-log" ) 
-      log = true;
-    else if (args[i] == "-shortlog" ) 
-      shortlog = true;
     else if (args[i] == "-lmp") 
       {
       ++i;
@@ -132,8 +127,6 @@ int main(int argc, char **argv)
   MPI_Comm_rank(sim_comm, &globalInfo.me);
 
   // Initialize SENSEI bridge 
-  if (0 == globalInfo.me) 
-    std::cout << "###### SENSEI instrumentation: initialize bridge ######" << std::endl;
   lammpsBridge::Initialize(sim_comm, sensei_xml );
 
   LAMMPS *lammps = (LAMMPS *)lammps_open(lammps_args.size(), lammps_args.data(), sim_comm, NULL);
@@ -206,8 +199,6 @@ int main(int argc, char **argv)
   lammps_close(lammps);
 
   // Finalize SENSEI bridge 
-  if (0 == globalInfo.me) 
-    std::cout << "###### SENSEI instrumentation: finalize bridge ######" << std::endl;
   lammpsBridge::Finalize();
 
   // close down MPI

From 4a83a1e72efc8a585282ff22e3f741c64c155b1d Mon Sep 17 00:00:00 2001
From: Silvio Rizzi <srizzi@anl.gov>
Date: Tue, 10 Aug 2021 03:11:10 +0000
Subject: [PATCH 08/11] Removed a SENSEI warning message in
 lammpsDataAdaptor.cxx

---
 miniapps/lammps/lammpsDataAdaptor.cxx | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/miniapps/lammps/lammpsDataAdaptor.cxx b/miniapps/lammps/lammpsDataAdaptor.cxx
index 3b42679b4..4ddf5fd69 100644
--- a/miniapps/lammps/lammpsDataAdaptor.cxx
+++ b/miniapps/lammps/lammpsDataAdaptor.cxx
@@ -456,7 +456,7 @@ int lammpsDataAdaptor::GetMeshMetadata(unsigned int id, sensei::MeshMetadataPtr
     {
     DInternals& internals = (*this->Internals);
     int nCells = internals.nlocal;
-    SENSEI_WARNING("lammps data adaptor. Flags.BlockSizeSet()")
+    //SENSEI_WARNING("lammps data adaptor. Flags.BlockSizeSet()")
     metadata->BlockNumCells.push_back(nCells);
     metadata->BlockNumPoints.push_back(nCells);
     metadata->BlockCellArraySize.push_back(2 * nCells); // fixme - VTK_POINTS

From 45aca31ac83c0ee4647495e2bfb20f6b00c4a84f Mon Sep 17 00:00:00 2001
From: Isaac Nealey <inealey@ucsd.edu>
Date: Thu, 10 Feb 2022 14:33:49 -0800
Subject: [PATCH 09/11] added profiler initialize

---
 miniapps/lammps/lammpsDriver.cxx | 3 +++
 1 file changed, 3 insertions(+)

diff --git a/miniapps/lammps/lammpsDriver.cxx b/miniapps/lammps/lammpsDriver.cxx
index f78f7e1c6..0f526c874 100644
--- a/miniapps/lammps/lammpsDriver.cxx
+++ b/miniapps/lammps/lammpsDriver.cxx
@@ -126,6 +126,9 @@ int main(int argc, char **argv)
   MPI_Comm sim_comm = MPI_COMM_WORLD;
   MPI_Comm_rank(sim_comm, &globalInfo.me);
 
+  // Initialize Profiler
+  sensei::Profiler::Initialize();
+
   // Initialize SENSEI bridge 
   lammpsBridge::Initialize(sim_comm, sensei_xml );
 

From dfeabae8d013251d403f76f6517930776c9f2036 Mon Sep 17 00:00:00 2001
From: Silvio Rizzi <31414702+srizzi88@users.noreply.github.com>
Date: Sun, 20 Feb 2022 20:07:55 -0600
Subject: [PATCH 10/11] Revert "added profiler initialize"

---
 miniapps/lammps/lammpsDriver.cxx | 3 ---
 1 file changed, 3 deletions(-)

diff --git a/miniapps/lammps/lammpsDriver.cxx b/miniapps/lammps/lammpsDriver.cxx
index 0f526c874..f78f7e1c6 100644
--- a/miniapps/lammps/lammpsDriver.cxx
+++ b/miniapps/lammps/lammpsDriver.cxx
@@ -126,9 +126,6 @@ int main(int argc, char **argv)
   MPI_Comm sim_comm = MPI_COMM_WORLD;
   MPI_Comm_rank(sim_comm, &globalInfo.me);
 
-  // Initialize Profiler
-  sensei::Profiler::Initialize();
-
   // Initialize SENSEI bridge 
   lammpsBridge::Initialize(sim_comm, sensei_xml );
 

From 976033f95e19dccc6405aeec99fbe4af11c3177f Mon Sep 17 00:00:00 2001
From: srizzi <srizzi@anl.gov>
Date: Mon, 28 Feb 2022 22:21:04 -0600
Subject: [PATCH 11/11] added profiler in lammps miniapp and sensei endpoint

---
 endpoints/SENSEIEndPoint.cpp     | 1 +
 miniapps/lammps/lammpsDriver.cxx | 7 +++++++
 2 files changed, 8 insertions(+)

diff --git a/endpoints/SENSEIEndPoint.cpp b/endpoints/SENSEIEndPoint.cpp
index f95e1fdb2..2ff7207fe 100644
--- a/endpoints/SENSEIEndPoint.cpp
+++ b/endpoints/SENSEIEndPoint.cpp
@@ -17,6 +17,7 @@ using AnalysisAdaptorPtr = vtkSmartPointer<sensei::ConfigurableAnalysis>;
 
 int main(int argc, char **argv)
 {
+  sensei::Profiler::Initialize();
   sensei::MPIManager mpiMan(argc, argv);
   int rank = mpiMan.GetCommRank();
 
diff --git a/miniapps/lammps/lammpsDriver.cxx b/miniapps/lammps/lammpsDriver.cxx
index f78f7e1c6..48b139188 100644
--- a/miniapps/lammps/lammpsDriver.cxx
+++ b/miniapps/lammps/lammpsDriver.cxx
@@ -37,6 +37,10 @@
 
 // SENSEI bridge 
 #include "lammpsBridge.h"   
+// SENSEI profiler
+#include "Profiler.h"
+
+
 
 using namespace LAMMPS_NS;
 
@@ -126,6 +130,9 @@ int main(int argc, char **argv)
   MPI_Comm sim_comm = MPI_COMM_WORLD;
   MPI_Comm_rank(sim_comm, &globalInfo.me);
 
+  // Initialize SENSEI profiler
+  sensei::Profiler::Initialize();
+
   // Initialize SENSEI bridge 
   lammpsBridge::Initialize(sim_comm, sensei_xml );