Source code for

# This file is part of the pyMOR project (
# Copyright 2013-2018 pyMOR developers and contributors. All rights reserved.
# License: BSD 2-Clause License (

import numpy as np

from pymor.core.config import config
from pymor.grids import referenceelements
from pymor.grids.constructions import flatten_grid

if config.HAVE_PYVTK:
        from evtk.hl import _addDataToFile, _appendDataToFile
        from evtk.vtk import VtkGroup, VtkFile, VtkUnstructuredGrid, VtkTriangle, VtkQuad
    except ImportError:
        from pyevtk.hl import _addDataToFile, _appendDataToFile
        from pyevtk.vtk import VtkGroup, VtkFile, VtkUnstructuredGrid, VtkTriangle, VtkQuad

def _write_vtu_series(grid, coordinates, connectivity, data, filename_base, last_step, is_cell_data):
    steps = last_step + 1 if last_step is not None else len(data)
    fn_tpl = "{}_{:08d}"

    npoints = len(coordinates[0])
    ncells = len(connectivity)

    ref = grid.reference_element
    if ref is ref is referenceelements.triangle:
        points_per_cell = 3
        vtk_el_type = VtkTriangle.tid
    elif ref is referenceelements.square:
        points_per_cell = 4
        vtk_el_type = VtkQuad.tid
        raise NotImplementedError("vtk output only available for grids with triangle or rectangle reference elments")

    connectivity = connectivity.reshape(-1)
    cell_types = np.empty(ncells, dtype='uint8')
    cell_types[:] = vtk_el_type
    offsets = np.arange(start=points_per_cell, stop=ncells*points_per_cell+1, step=points_per_cell, dtype='int32')

    group = VtkGroup(filename_base)
    for i in range(steps):
        fn = fn_tpl.format(filename_base, i)
        vtk_data = data[i, :]
        w = VtkFile(fn, VtkUnstructuredGrid)
        w.openPiece(ncells=ncells, npoints=npoints)

        w.addData("Coordinates", coordinates)
        w.addData("connectivity", connectivity)
        w.addData("offsets", offsets)
        w.addData("types", cell_types)
        if is_cell_data:
            _addDataToFile(w, cellData={"Data": vtk_data}, pointData=None)
            _addDataToFile(w, cellData=None, pointData={"Data": vtk_data})

        if is_cell_data:
            _appendDataToFile(w, cellData={"Data": vtk_data}, pointData=None)
            _appendDataToFile(w, cellData=None, pointData={"Data": vtk_data})
        group.addFile(filepath=fn + '.vtu', sim_time=i)

[docs]def write_vtk(grid, data, filename_base, codim=2, binary_vtk=True, last_step=None): """Output grid-associated data in (legacy) vtk format Parameters ---------- grid A |Grid| with triangular or rectilinear reference element. data |VectorArray| with either cell (ie one datapoint per codim 0 entity) or vertex (ie one datapoint per codim 2 entity) data in each array element. codim the codimension associated with the data filename_base common component for output files in timeseries binary_vtk if false, output files contain human readable inline ascii data, else appended binary last_step if set must be <= len(data) to restrict output of timeseries """ if not config.HAVE_PYVTK: raise ImportError('could not import pyevtk') if grid.dim != 2: raise NotImplementedError if codim not in (0, 2): raise NotImplementedError subentities, coordinates, entity_map = flatten_grid(grid) data = data.to_numpy() if codim == 0 else data.to_numpy()[:, entity_map].copy() x, y, z = coordinates[:, 0].copy(), coordinates[:, 1].copy(), np.zeros(coordinates[:, 1].size) _write_vtu_series(grid, coordinates=(x, y, z), connectivity=subentities, data=data, filename_base=filename_base, last_step=last_step, is_cell_data=(codim == 0))