PlotParticlesInVTKFormat.py

Go to the documentation of this file.

# This file is part of the Peano project. For conditions of distribution and
# use, please see the copyright notice at www.peano-framework.org
from peano4.solversteps.ActionSet import ActionSet
 
 
import jinja2
 
 
class PlotParticlesInVTKFormat(ActionSet):
    """
 
    For constructor's time_stamp_evaluation argument.
 
    """
 
    NoMetaFile = "no-meta-file"
 
    """
 
   For constructor's time_stamp_evaluation argument.
 
  """
    CountTimeSteps = "count-time-steps"
 
    """
 
   By default, I plot the point particles (obviously) and drop their
   coordinates as vector, too. I also plot the search radius which
   is a simple scalar. I often use it to colour-code or as radius
   of the rendered sphere.
 
   Finally, I do render the associativity. This is one vector per
   point that points towards the vertex the object belongs to. To
   visualise it, add a Glyph filter and select association as
   orientation array. Make association also the scale array and
   scale by magnitude.
 
  """
 
    def __init__(
        self,
        filename,
        particle_set,
        time_stamp_evaluation=NoMetaFile,
        additional_includes="",
    ):
        """
        Plot only the grid structure
 
        :Attibutes:
 
        filename: String
           Name of the output file.
 
        particle: ParticleSet
           I take this as particle set and I assume that it yields a C++
           container which I can traverse.
 
        time_stamp_evaluation: String yielding a C++ double
           Pass in prefined constants to use default behaviour.
 
        """
        super(PlotParticlesInVTKFormat, self).__init__(descend_invocation_order=1, parallel=False)
 
        self.d = {}
        self.d["FILENAME"] = filename
        self.d["PARTICLE"] = particle_set.particle_model.name
        self.d["PARTICLES_CONTAINER"] = particle_set.name
        self.d["WRITE_BINARY"] = "false"
        self.d["ATTRIBUTE_WRITER_NULL_INITIALISERS"] = ""
        self.d["ATTRIBUTE_WRITER_CLOSES"] = ""
        self.d["ATTRIBUTE_WRITER_DELETES"] = ""
        self.d["ATTRIBUTE_WRITER_INITIALISERS"] = ""
        self.d["ATTRIBUTE_WRITER_PLOT_CALLS"] = ""
        self.d["TIMESTAMP"] = time_stamp_evaluation
 
        self._attributes_to_plot = []
 
        self.additional_includes = additional_includes
 
    __Template_Constructor = jinja2.Template(
        """
  _writer             = nullptr;
  _positionWriter     = nullptr;
  _searchRadiusWriter       = nullptr;
  _associationWriter  = nullptr;
  {{ATTRIBUTE_WRITER_NULL_INITIALISERS}}
  _treeNumber         = treeNumber;
 
  // An MPI lock (critical section) would be important!
 
  logDebug( "PlotGrid2PlotGridInPeanoBlockFormat1()", "created tree instance for " << treeNumber );
"""
    )
 
    def add_attribute_to_plot(self, attribute, cardinality=1, accessor=None):
        """
 
        Add an attribute to the plotter
 
        You have to ensure that each attribute is added at most once.
 
 
        ## Type casts for non-floating point types
 
        The plotter supports only floating point data. If you have attributes that are
        not of floating point type, then you have to convert them into floating point
        data before they are dumped into VTK. For this, the routine offers the accessor
        attribute. As long as it equals None, we assume that there is a pointer p pointing
        to the particle, and that we can use getXXXX() to access the attribute data. The
        XXXX is the name of the attribute with the first letter written uppercase. This
        is how Peano 4 writes down accessors (using DaStGen - but that's a different story
        how and why we use certain code generator conventions).
 
        You can however inject any other accessor expression (using the pointer p) to
        read out particle data. If your particle data is for example a vector over
        integers, you might want to hand in
 
              tarch::la::convertScalar<double>(p->getXXXX())
 
        which converts the data into a proper format.
 
 
        ## Arguments
 
        attribute: Attribute
          The actual attribute that we have to plot.
 
        cardinality: Integer
          The cardinality of the attribute. This is 1 for a scalar value, 2 for a 2d vector,
          and so forth. We support arbitrary cardinalities, but the cardinality hsa to match
          attribute.
 
        accessor: String or None
          If this term is None, then we use getXXX() constructs as they are generated by
          DaStGen by default. If you hand in something not equal to None, then we use this
          particular string to access a field of a pointer p, i.e. we use
 
                p->accessor
 
          to get the data out of the particle.
 
        """
        assert cardinality > 0, "cardinality has to be greater than 0"
 
        self._attributes_to_plot.append(attribute)
        writer_name = f"_{attribute._name}Writer"
        delim = "\n  "
        if accessor != None:
            accessor_name = accessor
        elif cardinality <= 3:
            accessor_name = (
                "p->get" + attribute._name[0].title() + attribute._name[1:] + "()"
            )
        else:
            accessor_name = (
                "p->get"
                + attribute._name[0].title()
                + attribute._name[1:]
                + "().data()"
            )
 
        self.d[
            "ATTRIBUTE_WRITER_NULL_INITIALISERS"
        ] += f"{delim}{writer_name} = nullptr;"
        self.d["ATTRIBUTE_WRITER_CLOSES"] += f"{delim}{writer_name}->close();"
        self.d["ATTRIBUTE_WRITER_DELETES"] += f"{delim}delete {writer_name};"
        self.d[
            "ATTRIBUTE_WRITER_INITIALISERS"
        ] += f'{delim}{writer_name} = _writer->createPointDataWriter( "{attribute._name}", {cardinality} );'
        if cardinality <= 3:
            self.d[
                "ATTRIBUTE_WRITER_PLOT_CALLS"
            ] += f"{delim}  {writer_name}->plot(particleNumber,{accessor_name});"
        else:
            self.d[
                "ATTRIBUTE_WRITER_PLOT_CALLS"
            ] += f"{delim}  {writer_name}->plot(particleNumber,{accessor_name},{cardinality});"
 
    def get_constructor_body(self):
        return self.__Template_Constructor.render(**self.d)
 
    def get_destructor_body(self):
        return ""
 
    def get_body_of_getGridControlEvents(self):
        return "  return std::vector< peano4::grid::GridControlEvent >();\n"
 
    def get_action_set_name(self):
        return __name__.replace(".py", "").replace(".", "_")
 
    def user_should_modify_template(self):
        return False
 
    __Template_TouchVertexFirstTime = jinja2.Template(
        """
  assertion( _positionWriter!=nullptr );
  for (auto* p: fineGridVertex{{PARTICLES_CONTAINER}}) {
    if (
      p->getParallelState()==globaldata::{{PARTICLE}}::ParallelState::Local
    ) {
      int particleNumber = _writer->plotPoint(p->getX());
      _positionWriter->plot(particleNumber,p->getX());
      _searchRadiusWriter->plot(particleNumber,p->getSearchRadius());
      _associationWriter->plot(particleNumber,marker.x()-p->getX());
      {{ATTRIBUTE_WRITER_PLOT_CALLS}}
    }
  }
"""
    )
 
    __Template_BeginTraversal = jinja2.Template(
        """
  tarch::mpi::Lock lock( _semaphore );
 
  static int counter = -1;
  counter++;
 
  std::ostringstream snapshotFileName;
  snapshotFileName << "{{FILENAME}}-" << counter;
 
  if (tarch::mpi::Rank::getInstance().getNumberOfRanks()>1 ) {
    snapshotFileName << "-rank-" << tarch::mpi::Rank::getInstance().getRank();
  }
 
  {% if TIMESTAMP==\""""
        + NoMetaFile
        + """\" %}
  _writer = new tarch::plotter::pointdata::vtk::VTKWriter(
    false, snapshotFileName.str(), "{{FILENAME}}",
    tarch::plotter::PVDTimeSeriesWriter::IndexFileMode::NoIndexFile,
    0.0
  );
  {% elif TIMESTAMP==\""""
        + CountTimeSteps
        + """\" %}
  _writer = new tarch::plotter::pointdata::vtk::VTKWriter(
    false, snapshotFileName.str(), "{{FILENAME}}",
    tarch::plotter::PVDTimeSeriesWriter::IndexFileMode::AppendNewData,
    counter / peano4::parallel::SpacetreeSet::getInstance().getLocalSpacetrees().size()
  );
  {% else %}
  _writer = new tarch::plotter::pointdata::vtk::VTKWriter(
    false, snapshotFileName.str(), "{{FILENAME}}",
    tarch::plotter::PVDTimeSeriesWriter::IndexFileMode::AppendNewData,
    {{TIMESTAMP}}
  );
  {% endif %}
 
  _positionWriter    = _writer->createPointDataWriter( "x", 3 );
  _searchRadiusWriter      = _writer->createPointDataWriter( "search-radius", 1 );
  _associationWriter = _writer->createPointDataWriter( "association", 3 );
  {{ATTRIBUTE_WRITER_INITIALISERS}}
"""
    )
 
    __Template_EndTraversal = jinja2.Template(
        """
  assertion(_positionWriter!=nullptr);
 
  _positionWriter->close();
  _searchRadiusWriter->close();
  _associationWriter->close();
  {{ATTRIBUTE_WRITER_CLOSES}}
  _writer->writeToFile();
 
  delete _positionWriter;
  delete _searchRadiusWriter;
  delete _associationWriter;
  {{ATTRIBUTE_WRITER_DELETES}}
  delete _writer;
 
  _positionWriter    = nullptr;
  _searchRadiusWriter      = nullptr;
  _associationWriter = nullptr;
  {{ATTRIBUTE_WRITER_NULL_INITIALISERS}}
  _writer            = nullptr;
"""
    )
 
    def get_body_of_operation(self, operation_name):
        result = "\n"
        if operation_name == ActionSet.OPERATION_TOUCH_VERTEX_FIRST_TIME:
            result = self.__Template_TouchVertexFirstTime.render(**self.d)
        if operation_name == ActionSet.OPERATION_BEGIN_TRAVERSAL:
            result = self.__Template_BeginTraversal.render(**self.d)
        if operation_name == ActionSet.OPERATION_END_TRAVERSAL:
            result = self.__Template_EndTraversal.render(**self.d)
        return result
 
    def _get_attribute_writer_declarations(self):
        delim = "\n    "
        datatype = "tarch::plotter::pointdata::PointWriter::PointDataWriter*"
        return delim.join(
            f"{datatype}  _{attr._name}Writer;" for attr in self._attributes_to_plot
        )
 
    def get_attributes(self):
        return f"""
    static tarch::mpi::BooleanSemaphore                                              _semaphore;
 
    int                _treeNumber;
 
    tarch::plotter::pointdata::vtk::VTKWriter*                _writer;
    tarch::plotter::pointdata::PointWriter::PointDataWriter*  _positionWriter;
    tarch::plotter::pointdata::PointWriter::PointDataWriter*  _searchRadiusWriter;
    tarch::plotter::pointdata::PointWriter::PointDataWriter*  _associationWriter;
    {self._get_attribute_writer_declarations()}
"""
 
    def get_includes(self):
        return (
            """
#include "tarch/plotter/pointdata/vtk/VTKWriter.h"
#include "tarch/multicore/Lock.h"
#include "tarch/multicore/BooleanSemaphore.h"
#include "tarch/mpi/Lock.h"
#include "tarch/mpi/BooleanSemaphore.h"
#include "peano4/parallel/SpacetreeSet.h"
#include "../vertexdata/"""
            + self.d["PARTICLES_CONTAINER"]
            + """.h"
#include "../globaldata/"""
            + self.d["PARTICLE"]
            + """.h"
"""
            + self.additional_includes
        )
 
    def get_static_initialisations(self, full_qualified_classname):
        return (
            """
tarch::mpi::BooleanSemaphore  """
            + full_qualified_classname
            + """::_semaphore;
"""
        )