/*=========================================================================

  Program:   Visualization Toolkit
  Module:    vtkLSDynaReader.h

  Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
  All rights reserved.
  See Copyright.txt or http://www.kitware.com/Copyright.htm for details.

     This software is distributed WITHOUT ANY WARRANTY; without even
     the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
     PURPOSE.  See the above copyright notice for more information.

=========================================================================*/
/*----------------------------------------------------------------------------
 Copyright (c) Sandia Corporation
 See Copyright.txt or http://www.paraview.org/HTML/Copyright.html for details.
----------------------------------------------------------------------------*/

/**
 * @class   vtkLSDynaReader
 * @brief   Read LS-Dyna databases (d3plot)
 *
 * This filter reads LS-Dyna databases.
 *
 * The Set/GetFileName() routines are actually wrappers around the
 * Set/GetDatabaseDirectory() members; the actual filename you choose is
 * irrelevant -- only the directory name is used.  This is done in order to
 * accommodate ParaView.
 *
 * Note that this reader produces 7 output meshes.
 * These meshes are required as several attributes are defined on subsets
 * of the mesh.  Below is a list of meshes in the order they are output and
 * an explanation of which attributes are unique to each mesh:
 * - solid (3D) elements: number of integration points are different than 2D
 * - thick shell elements: number of integration points are different than
 * planar 2D
 * - shell (2D) elements: number of integration points are different than 3D
 * - rigid surfaces: can't have deflection, only velocity, accel, etc.
 * - road surfaces: have only a "segment ID" (serves as material ID) and a
 * velocity.
 * - beam elements: have Frenet (TNB) frame and cross-section attributes
 * (shape and size)
 * - spherical particle hydrodynamics (SPH) elements: have a radius of
 * influence, internal energy, etc.
 * Because each mesh has its own cell attributes, the vtkLSDynaReader has a
 * rather large API.  Instead of a single set of routines to query and set
 * cell array names and status, one exists for each possible output mesh.
 * Also, GetNumberOfCells() will return the sum of all the cells in all 7
 * meshes.  If you want the number of cells in a specific mesh, there are
 * separate routines for each mesh type.
 *
 *
 * @attention
 * LSDyna files contain 3 different types of sections: control, data, and
 * state.  Control sections contain constants that describe the type of
 * simulation data in a file or group of files.  Data sections contain
 * simulation information that is invariant across individual time steps
 * (but can vary when a mesh adaptation occurs).  This information includes
 * material, connectivity, and undeformed geometry.  Finally, state data is
 * information that varies with each time step.  Unless a mesh adaptation
 * occurs, there will be a single control and data section, and they will
 * be located at the start of the database (the first file).
 *
 * @attention
 * In their infinite wisdom, LSDyna developers decided to split simulation
 * data into multiple files, each no larger than some predetermined limit.
 * Each file can contain one section, a partial section (if it would not
 * fit into a single file), or multiple sections. Files are padded with
 * zeros so that their lengths will be multiples of 512*512.  The size of
 * each section is determined by constants in the control and data
 * sections, which means that these must be parsed carefully in order to
 * correctly locate desired information.  Unfortunately, the constants are
 * not terribly well-documented and in some cases the documentation is in
 * error.
 *
 * @par "Open Issues":
 * The LS-Dyna file format document leaves a good bit open to
 * interpretation.  In addition to the "documentation vs. files in the
 * wild" issues there are also implementation problems.
 *
 * @par "Open Issues":
 * - Where exactly may breaks to a new file occur in the pre-state
 * information? At each section?
 * - Will state data sections (node/cell data, element deletion, sph data,
 * rigid body motion) be moved to the beginning of a new file if their data
 * will be too large for a given file, or are all the sections
 * counted together as a single state (makes more sense for keeping time
 * word at start of every file).
 *  The questions above arise because the docs (p. 3) state "There are 3
 * sections in this database." but then call many smaller pieces of data
 * "sections". Should they be subsections? The docs are quiet about whether
 * the second section (of 3) is ever split across multiple files and, if
 * so, whether it is done at (sub)section boundaries when possible or just
 * wherever it needs to occur.
 * - How many components does Eddy Viscosity have? It's shown as 7 bits in
 * NCFDV1 which makes no sense at all.
 * - Why is NARBS larger than 10+NUMNP+NEL8+NEL2+NEL4+NELT (which is the
 * value specified by the documentation)?
 * Obviously, NARBS is definitive, but what are the extra numbers at the end?
 * - Is there a difference between rigid body elements NUMRBE and rigid road
 * surfaces? It appears that the nodes and connectivity of the road surface
 * are given separately (p.13) while on p.7 the Material
 *   Type Data subsection says that shells in a rigid body will just have a
 * certain material ID but be interspersed among deformable shell elements.
 * - Word 37 of the control section serves two possible purposes... it can
 * mean NMSPH or EDLOPT.
 *   I assume that different versions of the code use that word differently.
 * How do we know the difference?
 * - It's unclear how much state isn't stored when a shell element is marked
 * as rigid. Specifically, is element deletion data stored for rigid shells?
 * Page 21 of the spec is mute on this.
 * - The loop to read cell User IDs won't work if Rigid Body and Shell
 * elements are interleaved (which I now believe they are).
 *
 * @par "Open Issues":
 * On the VTK side of things:
 * - Berk has nudged me towards multiblock outputs but hasn't committed to
 * exactly how things can be made efficient for a parallel version of the
 * reader.
 * - This reader will eventually need to respond to a second output port for
 * "small spatial, large temporal" queries.
 * - The reader doesn't handle crack files (d3crck)
 * - The reader doesn't handle interface force files (no default name)
 * - The reader doesn't handle time history (abbreviated output) files (d3thdt)
 * - The reader doesn't handle dynamic relaxation files (d3drfl)
 * - The reader doesn't handle reduced parts (state for a subset of parts) files (d3part)
 * - The reader doesn't handle mode shape files (d3eigv)
 * - The reader doesn't handle equilibrium iteration files (d3iter)
 * - The reader doesn't handle extra time data files (d3xtf)
 * - The reader doesn't handle printer files (d3hsp)
 * - The reader doesn't handle modal neutral files (d3mnf)
 * - The reader doesn't handle packed connectivity.
 * - The reader doesn't handle adapted element parent lists (but the 2002 specification says LSDyna
 * doesn't implement it).
 * - All the sample datasets have MATTYP = 0. Need something to test MATTYP = 1.
 * - I have no test datasets with rigid body and/or road surfaces, so the
 * implementation is half-baked.
 * - It's unclear how some of the data should be presented. Although blindly
 * tacking the numbers into a large chuck of cell data is better than nothing,
 * some attributes (e.g., forces & moments) lend themselves to more elaborate
 * presentation. Also, shell and thick shell elements have stresses that
 * belong to a particular side of an element or have a finite thickness that
 * could be rendered.
 *   Finally, beam elements have cross sections that could be rendered.
 * Some of these operations require numerical processing of the results and
 * so we shouldn't eliminate the ability to get at the raw simulation data.
 * Perhaps a filter could be applied to "fancify" the geometry.
 *
 */

#ifndef vtkLSDynaReader_h
#define vtkLSDynaReader_h

#include "vtkIOLSDynaModule.h" // For export macro
#include "vtkLegacy.h"         // For VTK_LEGACY_REMOVE
#include "vtkMultiBlockDataSetAlgorithm.h"
#include <string> // for method signature

class LSDynaMetaData;
class vtkLSDynaPartCollection;
class vtkPoints;
class vtkDataArray;
class vtkUnsignedCharArray;
class vtkUnstructuredGrid;

class VTKIOLSDYNA_EXPORT vtkLSDynaReader : public vtkMultiBlockDataSetAlgorithm
{
public:
  vtkTypeMacro(vtkLSDynaReader, vtkMultiBlockDataSetAlgorithm);
  void PrintSelf(ostream& os, vtkIndent indent) override;
  static vtkLSDynaReader* New();

  /**
   * Print out more complete information about the dataset
   * (and less complete information about the VTK hierarchy) than PrintSelf.
   */
  void Dump(ostream& os);

  /**
   * A routine to call Dump() from within a lame debugger that won't
   * properly pass a C++ iostream object like cout.
   */
  void DebugDump();

  /**
   * Determine if the file can be read with this reader.
   */
  virtual int CanReadFile(const char* fname);

  ///@{
  /**
   * Get/Set the directory containing the LS-Dyna database and determine
   * whether it is valid.
   */
  virtual void SetDatabaseDirectory(const std::string&);
  virtual void SetDatabaseDirectory(const char*);
#ifdef VTK_LEGACY_REMOVE
  std::string GetDatabaseDirectory();
#else
  VTK_LEGACY(const char* GetDatabaseDirectory());
#endif
  int IsDatabaseValid();
  ///@}

  ///@{
  /**
   * Get/Set the filename. The Set/GetFileName() routines are actually
   * wrappers around the Set/GetDatabaseDirectory() members; the actual
   * filename you choose is irrelevant -- only the directory name is used.
   * This is done in order to accommodate ParaView.
   */
  virtual void SetFileName(const std::string&);
  virtual void SetFileName(const char*);
#ifdef VTK_LEGACY_REMOVE
  std::string GetFileName();
#else
  VTK_LEGACY(const char* GetFileName());
#endif
  ///@}

  /**
   * The title of the database is a 40 or 80 character text description
   * stored at the front of a d3plot file.  Do not call this function
   * before setting the database directory and calling UpdateInformation().
   */
  char* GetTitle();

  /**
   * Retrieve the dimension of points in the database. This should return 2
   * or 3.  Do not call this function before setting the database directory
   * and calling UpdateInformation().
   */
  int GetDimensionality();

  /**
   * Retrieve the number of points in the database.  Do not call this
   * function before setting the database directory and calling
   * UpdateInformation().
   */
  vtkIdType GetNumberOfNodes();

  /**
   * Retrieve the number of cells of a given type in the database.  Do not
   * call this function before setting the database directory and calling
   * UpdateInformation().

   * Note that GetNumberOfCells() returns the sum of
   * GetNumberOfContinuumCells() and GetNumberOfParticleCells().
   */
  vtkIdType GetNumberOfCells();

  /**
   * Retrieve the number of cells of a given type in the database.  Do not
   * call this function before setting the database directory and calling
   * UpdateInformation().

   * Note that GetNumberOfContinuumCells() returns the sum of
   * GetNumberOfSolidCells(), GetNumberOfThickShellCells(),
   * GetNumberOfShellCells(), GetNumberOfRigidBodyCells(),
   * GetNumberOfRoadSurfaceCells(), and GetNumberOfBeamCells().
   */
  vtkIdType GetNumberOfContinuumCells();

  /**
   * Retrieve the number of cells of a given type in the database.  Do not
   * call this function before setting the database directory and calling
   * UpdateInformation().
   */
  vtkIdType GetNumberOfSolidCells();

  /**
   * Retrieve the number of cells of a given type in the database.  Do not
   * call this function before setting the database directory and calling
   * UpdateInformation().
   */
  vtkIdType GetNumberOfThickShellCells();

  /**
   * Retrieve the number of cells of a given type in the database.  Do not
   * call this function before setting the database directory and calling
   * UpdateInformation().
   */
  vtkIdType GetNumberOfShellCells();

  /**
   * Retrieve the number of cells of a given type in the database.  Do not
   * call this function before setting the database directory and calling
   * UpdateInformation().
   */
  vtkIdType GetNumberOfRigidBodyCells();

  /**
   * Retrieve the number of cells of a given type in the database.  Do not
   * call this function before setting the database directory and calling
   * UpdateInformation().
   */
  vtkIdType GetNumberOfRoadSurfaceCells();

  /**
   * Retrieve the number of cells of a given type in the database.
   * Do not call this function before setting the database directory and calling
   * UpdateInformation().
   */
  vtkIdType GetNumberOfBeamCells();

  /**
   * Retrieve the number of cells of a given type in the database.  Do not
   * call this function before setting the database directory and calling
   * UpdateInformation().
   */
  vtkIdType GetNumberOfParticleCells();

  ///@{
  /**
   * Retrieve information about the time extents of the LS-Dyna database.
   * Do not call these functions before setting the database directory and
   * calling UpdateInformation().
   */
  vtkIdType GetNumberOfTimeSteps();
  virtual void SetTimeStep(vtkIdType);
  vtkIdType GetTimeStep();
  double GetTimeValue(vtkIdType);
  vtkGetVector2Macro(TimeStepRange, int);
  vtkSetVector2Macro(TimeStepRange, int);
  ///@}

  ///@{
  /**
   * These methods allow you to load only selected subsets of the nodal
   * variables defined over the mesh.
   */
  int GetNumberOfPointArrays();
  const char* GetPointArrayName(int);
  virtual void SetPointArrayStatus(int arr, int status);
  virtual void SetPointArrayStatus(const char* arrName, int status);
  int GetPointArrayStatus(int arr);
  int GetPointArrayStatus(const char* arrName);
  int GetNumberOfComponentsInPointArray(int arr);
  int GetNumberOfComponentsInPointArray(const char* arrName);
  ///@}

  ///@{
  /**
   * Routines that allow the status of a cell variable to be adjusted or
   * queried independent of the output mesh.  The \a cellType parameter
   * should be one of: LS_POINT, LS_BEAM, LS_SHELL, LS_THICK_SHELL,
   * LS_SOLID, LS_RIGID_BODY, or LS_ROAD_SURFACE
   */
  int GetNumberOfCellArrays(int cellType);
  const char* GetCellArrayName(int cellType, int arr);
  virtual void SetCellArrayStatus(int cellType, int arr, int status);
  virtual void SetCellArrayStatus(int cellType, const char* arrName, int status);
  int GetCellArrayStatus(int cellType, int arr);
  int GetCellArrayStatus(int cellType, const char* arrName);
  int GetNumberOfComponentsInCellArray(int cellType, int arr);
  int GetNumberOfComponentsInCellArray(int cellType, const char* arrName);
  ///@}

  ///@{
  /**
   * These methods allow you to load only selected subsets of the cell
   * variables defined over the mesh.
   */
  int GetNumberOfSolidArrays();
  const char* GetSolidArrayName(int);
  virtual void SetSolidArrayStatus(int arr, int status);
  virtual void SetSolidArrayStatus(const char* arrName, int status);
  int GetSolidArrayStatus(int arr);
  int GetSolidArrayStatus(const char* arrName);
  ///@}

  int GetNumberOfComponentsInSolidArray(int a);
  int GetNumberOfComponentsInSolidArray(const char* arrName);

  ///@{
  /**
   * These methods allow you to load only selected subsets of the cell
   * variables defined over the mesh.
   */
  int GetNumberOfThickShellArrays();
  const char* GetThickShellArrayName(int);
  virtual void SetThickShellArrayStatus(int arr, int status);
  virtual void SetThickShellArrayStatus(const char* arrName, int status);
  int GetThickShellArrayStatus(int arr);
  int GetThickShellArrayStatus(const char* arrName);
  ///@}

  int GetNumberOfComponentsInThickShellArray(int a);
  int GetNumberOfComponentsInThickShellArray(const char* arrName);

  ///@{
  /**
   * These methods allow you to load only selected subsets of the cell
   * variables defined over the mesh.
   */
  int GetNumberOfShellArrays();
  const char* GetShellArrayName(int);
  virtual void SetShellArrayStatus(int arr, int status);
  virtual void SetShellArrayStatus(const char* arrName, int status);
  int GetShellArrayStatus(int arr);
  int GetShellArrayStatus(const char* arrName);
  ///@}

  int GetNumberOfComponentsInShellArray(int a);
  int GetNumberOfComponentsInShellArray(const char* arrName);

  ///@{
  /**
   * These methods allow you to load only selected subsets of the cell
   * variables defined over the mesh.
   */
  int GetNumberOfRigidBodyArrays();
  const char* GetRigidBodyArrayName(int);
  virtual void SetRigidBodyArrayStatus(int arr, int status);
  virtual void SetRigidBodyArrayStatus(const char* arrName, int status);
  int GetRigidBodyArrayStatus(int arr);
  int GetRigidBodyArrayStatus(const char* arrName);
  ///@}

  int GetNumberOfComponentsInRigidBodyArray(int a);
  int GetNumberOfComponentsInRigidBodyArray(const char* arrName);

  ///@{
  /**
   * These methods allow you to load only selected subsets of the cell
   * variables defined over the mesh.
   */
  int GetNumberOfRoadSurfaceArrays();
  const char* GetRoadSurfaceArrayName(int);
  virtual void SetRoadSurfaceArrayStatus(int arr, int status);
  virtual void SetRoadSurfaceArrayStatus(const char* arrName, int status);
  int GetRoadSurfaceArrayStatus(int arr);
  int GetRoadSurfaceArrayStatus(const char* arrName);
  ///@}

  int GetNumberOfComponentsInRoadSurfaceArray(int a);
  int GetNumberOfComponentsInRoadSurfaceArray(const char* arrName);

  ///@{
  /**
   * These methods allow you to load only selected subsets of the cell
   * variables defined over the mesh.
   */
  int GetNumberOfBeamArrays();
  const char* GetBeamArrayName(int);
  virtual void SetBeamArrayStatus(int arr, int status);
  virtual void SetBeamArrayStatus(const char* arrName, int status);
  int GetBeamArrayStatus(int arr);
  int GetBeamArrayStatus(const char* arrName);
  ///@}

  int GetNumberOfComponentsInBeamArray(int a);
  int GetNumberOfComponentsInBeamArray(const char* arrName);

  ///@{
  /**
   * These methods allow you to load only selected subsets of the cell
   * variables defined over the mesh.
   */
  int GetNumberOfParticleArrays();
  const char* GetParticleArrayName(int);
  virtual void SetParticleArrayStatus(int arr, int status);
  virtual void SetParticleArrayStatus(const char* arrName, int status);
  int GetParticleArrayStatus(int arr);
  int GetParticleArrayStatus(const char* arrName);
  ///@}

  int GetNumberOfComponentsInParticleArray(int a);
  int GetNumberOfComponentsInParticleArray(const char* arrName);

  ///@{
  /**
   * Should deflected coordinates be used, or should the mesh remain
   * undeflected?  By default, this is true but its value is ignored if the
   * nodal "Deflected Coordinates" array is not set to be loaded.
   */
  void SetDeformedMesh(vtkTypeBool);
  vtkGetMacro(DeformedMesh, vtkTypeBool);
  vtkBooleanMacro(DeformedMesh, vtkTypeBool);
  ///@}

  ///@{
  /**
   * Should dead cells be removed from the mesh?  Cells are marked dead by
   * setting the corresponding entry in the <b>cell</b> array "Death" to 0.
   * Cells that are not dead have the corresponding entry in the cell array
   * "Death" set to their material ID.  By default, this is true but its
   * value is ignored if the cell "Death" array is not set to be loaded.
   * It is also ignored if the database's element deletion option is set to
   * denote <b>points</b> (not cells) as deleted; in that case, "Death"
   * will appear to be a point array.
   */
  vtkSetMacro(RemoveDeletedCells, vtkTypeBool);
  vtkGetMacro(RemoveDeletedCells, vtkTypeBool);
  vtkBooleanMacro(RemoveDeletedCells, vtkTypeBool);
  ///@}

  ///@{
  /**
   * Instead of removing the cells that are dead, hide them by setting
   * the array as the ghost levels arrays
   */
  vtkSetMacro(DeletedCellsAsGhostArray, vtkTypeBool);
  vtkGetMacro(DeletedCellsAsGhostArray, vtkTypeBool);
  vtkBooleanMacro(DeletedCellsAsGhostArray, vtkTypeBool);
  ///@}

  ///@{
  /**
   * The name of the input deck corresponding to the current database.
   * This is used to determine the part names associated with each material ID.
   * This file may be in two formats: a valid LSDyna input deck or a
   * short XML summary.
   * If the file begins with "<?xml" then the summary format is used.
   * Otherwise, the keyword format is used and a summary file will be
   * created if write permissions exist in the directory containing
   * the keyword file. The newly created summary will have ".k" or ".key"
   * stripped from the end of the keyword filename and ".lsdyna" appended.
   */
  vtkSetStringMacro(InputDeck);
  vtkGetStringMacro(InputDeck);
  ///@}

  ///@{
  /**
   * These methods allow you to load only selected parts of the input.
   * If InputDeck points to a valid keyword file (or summary), then part
   * names will be taken from that file.
   * Otherwise, when arbitrary material numbering is used, parts will be named
   * "PartXXX (MatlYYY)" where XXX is an increasing sequential number and YYY
   * is the respective material ID. If no input deck is specified and arbitrary
   * arbitrary material numbering is not used, parts will be named
   * "PartXXX" where XXX is a sequential material ID.
   */
  int GetNumberOfPartArrays();
  const char* GetPartArrayName(int);
  virtual void SetPartArrayStatus(int arr, int status);
  virtual void SetPartArrayStatus(const char* partName, int status);
  int GetPartArrayStatus(int arr);
  int GetPartArrayStatus(const char* partName);
  ///@}

protected:
  // holds all the parts and all the properties for each part
  vtkLSDynaPartCollection* Parts;

  /**
   * Should deflected coordinates be used, or should the mesh remain
   * undeflected?  By default, this is true.
   */
  vtkTypeBool DeformedMesh;

  ///@{
  /**
   * Should cells marked as deleted be removed from the mesh?
   * By default, this is true.
   */
  vtkTypeBool RemoveDeletedCells;
  vtkTypeBool DeletedCellsAsGhostArray;
  ///@}

  /**
   * The range of time steps available within a database.
   * Only valid after UpdateInformation() is called on the reader.
   */
  int TimeStepRange[2];

  /**
   * The name of a file containing part names and IDs.
   */
  char* InputDeck;

  vtkLSDynaReader();
  ~vtkLSDynaReader() override;

  /**
   * This function populates the reader's private dictionary with
   * information about the database.  It is called once from
   * RequestInformation() and once any time the adaptation level changes.
   * The adaptation level can change any time the current state(time) is
   * modified.  Upon success, -1 is returned. "Soft" failures return 0 and
   * "hard" failures return 1.
   */
  int ReadHeaderInformation(int currentAdaptLevel);

  /**
   * This function scans the list of files in the database and bookmarks the
   * start of each time step's state information.
   * Before this function is called:
   * - The database directory name must be set,
   * - ReadHeaderInformation must have been called for adaptation level 0, and
   * - The "read head" must be placed at the end of the first adaptation's geometry section.
   * Upon success, -1 is returned. "Soft" failures return 0 and "hard" failures return 1.
   */
  int ScanDatabaseTimeSteps();

  int RequestInformation(vtkInformation*, vtkInformationVector**, vtkInformationVector*) override;
  int RequestData(vtkInformation*, vtkInformationVector**, vtkInformationVector*) override;

  ///@{
  /**
   * These functions read various parts of the database.
   * The functions that take a vtkIdType argument must be passed
   * the current timestep.
   * Functions that do not take a timestep must have the read head
   * positioned to the start of their data sections.
   * These functions should only be called from within RequestData() since
   * they require the various output meshes to exist.
   */
  virtual int ReadTopology();
  virtual int ReadNodes();
  virtual int ReadPartSizes();
  virtual int ReadConnectivityAndMaterial();
  virtual int ReadUserIds();
  virtual int ReadState(vtkIdType);
  virtual int ReadNodeStateInfo(vtkIdType);
  virtual int ReadCellStateInfo(vtkIdType);
  virtual int ReadDeletion();
  virtual int ReadSPHState(vtkIdType);
  virtual int ComputeDeflectionAndUpdateGeometry(vtkUnstructuredGrid* grid);
  ///@}

  /**
   * Resets the Part information to the default state
   */
  virtual void ResetPartInfo();

  /**
   * Called from within ReadHeaderInformation() to read part names
   * associated with material IDs.
   */
  virtual int ReadInputDeck();

  /**
   * Called from within ReadHeaderInformation to read part names
   * from the end of the first d3plot file. Used if ReadInputDeck
   * fails.
   */
  virtual int ReadPartTitlesFromRootFile();

  /**
   * Called from within ReadHeaderInformation() to read arbitrary material
   * IDs (if present) or manufacture sequential material IDs (if not
   * present).
   */
  virtual int ReadUserMaterialIds();

  ///@{
  /**
   * ReadInputDeck determines the type of file (keyword or XML summary) and
   * calls one of these two routines to read the file.
   */
  int ReadInputDeckXML(istream& deck);
  int ReadInputDeckKeywords(istream& deck);
  ///@}

  /**
   * ReadInputDeckKeywords calls this function if it was successful in reading
   * part names for materials.
   */
  int WriteInputDeckSummary(const char* fname);

  /**
   * Read an array of deletion data.
   * This is used by ReadDeletion to actually read the data from the file
   * (as opposed to attach it to the proper place in the VTK dataset)
   * depending on the value of "MDLOPT".
   * The array passed to this routine is filled with ones if deleted, zero
   * it is not deleted
   * The number of tuples must be set on the array previous to calling
   * this routine.
   * Note: pos is the position in the size that the death value is store at
   */
  virtual void ReadDeletionArray(vtkUnsignedCharArray* arr, const int& pos, const int& size);

  /**
   * Read all the cell properties of a given part type
   */
  virtual void ReadCellProperties(const int& type, const int& numTuples);

  LSDynaMetaData* P;

  void ResetPartsCache();

private:
  // Helper templated methods to optimize reading. We cast the entire buffer
  // to a given type instead of casting each element to improve performance
  template <typename T>
  void FillDeletionArray(T* buffer, vtkUnsignedCharArray* arr, const vtkIdType& start,
    const vtkIdType& numCells, const int& deathPos, const int& cellSize);

  template <int wordSize, typename T>
  int FillTopology(T* buffer);

  template <typename T, int blockType, vtkIdType numWordsPerCell, vtkIdType cellLength>
  void ReadBlockCellSizes();

  template <typename T>
  int FillPartSizes();

  vtkLSDynaReader(const vtkLSDynaReader&) = delete;
  void operator=(const vtkLSDynaReader&) = delete;
};

inline void vtkLSDynaReader::SetPointArrayStatus(const char* arrName, int status)
{
  for (int a = 0; a < this->GetNumberOfPointArrays(); ++a)
  {
    if (strcmp(arrName, this->GetPointArrayName(a)) == 0)
    {
      this->SetPointArrayStatus(a, status);
      return;
    }
  }
  vtkWarningMacro("Point array \"" << arrName << "\" does not exist");
}

inline int vtkLSDynaReader::GetPointArrayStatus(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfPointArrays(); ++a)
  {
    if (strcmp(arrName, this->GetPointArrayName(a)) == 0)
    {
      return this->GetPointArrayStatus(a);
    }
  }
  // vtkWarningMacro( "Point array \"" << arrName << "\" does not exist" );
  return 0;
}

inline int vtkLSDynaReader::GetNumberOfComponentsInPointArray(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfPointArrays(); ++a)
  {
    if (strcmp(arrName, this->GetPointArrayName(a)) == 0)
    {
      return this->GetNumberOfComponentsInPointArray(a);
    }
  }
  // vtkWarningMacro( "Point array \"" << arrName << "\" does not exist" );
  return 0;
}

inline void vtkLSDynaReader::SetCellArrayStatus(int cellType, const char* arrName, int status)
{
  for (int a = 0; a < this->GetNumberOfCellArrays(cellType); ++a)
  {
    if (strcmp(arrName, this->GetCellArrayName(cellType, a)) == 0)
    {
      this->SetCellArrayStatus(cellType, a, status);
      return;
    }
  }
  vtkWarningMacro("Cell array \"" << arrName << "\" (type " << cellType << ") does not exist");
}

inline int vtkLSDynaReader::GetCellArrayStatus(int cellType, const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfCellArrays(cellType); ++a)
  {
    if (strcmp(arrName, this->GetCellArrayName(cellType, a)) == 0)
    {
      return this->GetCellArrayStatus(cellType, a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline int vtkLSDynaReader::GetNumberOfComponentsInCellArray(int cellType, const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfCellArrays(cellType); ++a)
  {
    if (strcmp(arrName, this->GetCellArrayName(cellType, a)) == 0)
    {
      return this->GetNumberOfComponentsInCellArray(cellType, a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline void vtkLSDynaReader::SetSolidArrayStatus(const char* arrName, int status)
{
  for (int a = 0; a < this->GetNumberOfSolidArrays(); ++a)
  {
    if (strcmp(arrName, this->GetSolidArrayName(a)) == 0)
    {
      this->SetSolidArrayStatus(a, status);
      return;
    }
  }
  vtkWarningMacro("Solid array \"" << arrName << "\" does not exist");
}

inline int vtkLSDynaReader::GetSolidArrayStatus(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfSolidArrays(); ++a)
  {
    if (strcmp(arrName, this->GetSolidArrayName(a)) == 0)
    {
      return this->GetSolidArrayStatus(a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline int vtkLSDynaReader::GetNumberOfComponentsInSolidArray(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfSolidArrays(); ++a)
  {
    if (strcmp(arrName, this->GetSolidArrayName(a)) == 0)
    {
      return this->GetNumberOfComponentsInSolidArray(a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline void vtkLSDynaReader::SetThickShellArrayStatus(const char* arrName, int status)
{
  for (int a = 0; a < this->GetNumberOfThickShellArrays(); ++a)
  {
    if (strcmp(arrName, this->GetThickShellArrayName(a)) == 0)
    {
      this->SetThickShellArrayStatus(a, status);
      return;
    }
  }
  vtkWarningMacro("Thick shell array \"" << arrName << "\" does not exist");
}

inline int vtkLSDynaReader::GetThickShellArrayStatus(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfThickShellArrays(); ++a)
  {
    if (strcmp(arrName, this->GetThickShellArrayName(a)) == 0)
    {
      return this->GetThickShellArrayStatus(a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline int vtkLSDynaReader::GetNumberOfComponentsInThickShellArray(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfThickShellArrays(); ++a)
  {
    if (strcmp(arrName, this->GetThickShellArrayName(a)) == 0)
    {
      return this->GetNumberOfComponentsInThickShellArray(a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline void vtkLSDynaReader::SetShellArrayStatus(const char* arrName, int status)
{
  for (int a = 0; a < this->GetNumberOfShellArrays(); ++a)
  {
    if (strcmp(arrName, this->GetShellArrayName(a)) == 0)
    {
      this->SetShellArrayStatus(a, status);
      return;
    }
  }
  vtkWarningMacro("Shell array \"" << arrName << "\" does not exist");
}

inline int vtkLSDynaReader::GetShellArrayStatus(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfShellArrays(); ++a)
  {
    if (strcmp(arrName, this->GetShellArrayName(a)) == 0)
    {
      return this->GetShellArrayStatus(a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline int vtkLSDynaReader::GetNumberOfComponentsInShellArray(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfShellArrays(); ++a)
  {
    if (strcmp(arrName, this->GetShellArrayName(a)) == 0)
    {
      return this->GetNumberOfComponentsInShellArray(a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline void vtkLSDynaReader::SetBeamArrayStatus(const char* arrName, int status)
{
  for (int a = 0; a < this->GetNumberOfBeamArrays(); ++a)
  {
    if (strcmp(arrName, this->GetBeamArrayName(a)) == 0)
    {
      this->SetBeamArrayStatus(a, status);
      return;
    }
  }
  vtkWarningMacro("Beam array \"" << arrName << "\" does not exist");
}

inline int vtkLSDynaReader::GetBeamArrayStatus(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfBeamArrays(); ++a)
  {
    if (strcmp(arrName, this->GetBeamArrayName(a)) == 0)
    {
      return this->GetBeamArrayStatus(a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline int vtkLSDynaReader::GetNumberOfComponentsInBeamArray(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfBeamArrays(); ++a)
  {
    if (strcmp(arrName, this->GetBeamArrayName(a)) == 0)
    {
      return this->GetNumberOfComponentsInBeamArray(a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline void vtkLSDynaReader::SetParticleArrayStatus(const char* arrName, int status)
{
  for (int a = 0; a < this->GetNumberOfParticleArrays(); ++a)
  {
    if (strcmp(arrName, this->GetParticleArrayName(a)) == 0)
    {
      this->SetParticleArrayStatus(a, status);
      return;
    }
  }
  vtkWarningMacro("Particle array \"" << arrName << "\" does not exist");
}

inline int vtkLSDynaReader::GetParticleArrayStatus(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfParticleArrays(); ++a)
  {
    if (strcmp(arrName, this->GetParticleArrayName(a)) == 0)
    {
      return this->GetParticleArrayStatus(a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline int vtkLSDynaReader::GetNumberOfComponentsInParticleArray(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfParticleArrays(); ++a)
  {
    if (strcmp(arrName, this->GetParticleArrayName(a)) == 0)
    {
      return this->GetNumberOfComponentsInParticleArray(a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline void vtkLSDynaReader::SetRigidBodyArrayStatus(const char* arrName, int status)
{
  for (int a = 0; a < this->GetNumberOfRigidBodyArrays(); ++a)
  {
    if (strcmp(arrName, this->GetRigidBodyArrayName(a)) == 0)
    {
      this->SetRigidBodyArrayStatus(a, status);
      return;
    }
  }
  vtkWarningMacro("Rigid body array \"" << arrName << "\" does not exist");
}

inline int vtkLSDynaReader::GetRigidBodyArrayStatus(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfRigidBodyArrays(); ++a)
  {
    if (strcmp(arrName, this->GetRigidBodyArrayName(a)) == 0)
    {
      return this->GetRigidBodyArrayStatus(a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline int vtkLSDynaReader::GetNumberOfComponentsInRigidBodyArray(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfRigidBodyArrays(); ++a)
  {
    if (strcmp(arrName, this->GetRigidBodyArrayName(a)) == 0)
    {
      return this->GetNumberOfComponentsInRigidBodyArray(a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline void vtkLSDynaReader::SetRoadSurfaceArrayStatus(const char* arrName, int status)
{
  for (int a = 0; a < this->GetNumberOfRoadSurfaceArrays(); ++a)
  {
    if (strcmp(arrName, this->GetRoadSurfaceArrayName(a)) == 0)
    {
      this->SetRoadSurfaceArrayStatus(a, status);
      return;
    }
  }
  vtkWarningMacro("Road surface array \"" << arrName << "\" does not exist");
}

inline int vtkLSDynaReader::GetRoadSurfaceArrayStatus(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfRoadSurfaceArrays(); ++a)
  {
    if (strcmp(arrName, this->GetRoadSurfaceArrayName(a)) == 0)
    {
      return this->GetRoadSurfaceArrayStatus(a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline int vtkLSDynaReader::GetNumberOfComponentsInRoadSurfaceArray(const char* arrName)
{
  for (int a = 0; a < this->GetNumberOfRoadSurfaceArrays(); ++a)
  {
    if (strcmp(arrName, this->GetRoadSurfaceArrayName(a)) == 0)
    {
      return this->GetNumberOfComponentsInRoadSurfaceArray(a);
    }
  }
  // vtkWarningMacro( "Cell array \"" << arrName << "\" does not exist" );
  return 0;
}

inline void vtkLSDynaReader::SetPartArrayStatus(const char* arrName, int status)
{
  for (int a = 0; a < this->GetNumberOfPartArrays(); ++a)
  {
    if (strcmp(arrName, this->GetPartArrayName(a)) == 0)
    {
      this->SetPartArrayStatus(a, status);
      return;
    }
  }
  vtkWarningMacro("Part \"" << arrName << "\" does not exist");
}

inline int vtkLSDynaReader::GetPartArrayStatus(const char* partName)
{
  for (int a = 0; a < this->GetNumberOfPartArrays(); ++a)
  {
    if (strcmp(partName, this->GetPartArrayName(a)) == 0)
    {
      return this->GetPartArrayStatus(a);
    }
  }
  // vtkWarningMacro( "PartArray \"" << partName << "\" does not exist" );
  return 0;
}

#endif // vtkLSDynaReader_h
