The TS-GitLab will have to shut down towards the end of the year — please think about migrating your projects to GitLab.com or GitHub.
(This is still a very early message, meant to keep you informed. There will be more urgent ones in the future.)

Commit d8eebfa9 authored by Santiago Ospina De Los Ríos's avatar Santiago Ospina De Los Ríos
Browse files

Merge branch 'master' into 137-add-a-finite-volume-local-operator-for-the-richards-equation

parents bb773998 b0290153
......@@ -10,7 +10,7 @@ variables:
BASE_IMAGE: dorie/dune-env
# Use semantic versioning (not the version of DUNE) and bump according to
# to whether changes are backwards-compatible or not.
IMAGE_VERSION: "1.0"
IMAGE_VERSION: "1.1"
DUNE_ENV_IMAGE: ${BASE_IMAGE}:img-v${IMAGE_VERSION}
CMAKE_FLAGS:
......@@ -184,6 +184,13 @@ test:unit-tests:
- $CI_PROJECT_DIR/build-cmake/dune/dorie/test
expire_in: 1 day
test:python-tests:
<<: *test
dependencies: []
script:
- $DUNECONTROL --only=dorie configure
- $DUNECONTROL --only=dorie make test_python
# --- Deploy jobs ---
deploy:dockerhub-devel: &deploy
stage: deploy
......
......@@ -39,6 +39,7 @@
[muparser](http://beltoforion.de/article.php?a=muparser) !131
* Coupling between transient models for water flow and solute transport !96
* Initial conditions generated from H5 input data !130
* Generic Python VTK file reader !143
### Changed
* `Simulation` is renamed `RichardsSimulation` and moved to
......
......@@ -78,6 +78,7 @@ by CI tests.
| [yaml-cpp](https://github.com/jbeder/yaml-cpp) | >= 5.2.0 |
| [spdlog](https://github.com/gabime/spdlog) | 1.1.0 | Included as Git Submodule
| [muparser](http://beltoforion.de/article.php?a=muparser) | master |
| [VTK](https://vtk.org/) | >= 7.1.1 | For the Python module only
| [dune-common](https://gitlab.dune-project.org/core/dune-common) | releases/2.6
| [dune-geometry](https://gitlab.dune-project.org/core/dune-geometry) | releases/2.6
| [dune-grid](https://gitlab.dune-project.org/core/dune-grid) | releases/2.6
......@@ -95,8 +96,8 @@ by CI tests.
| ---------| -------------- | -------- |
| [dune-testtools](https://gitlab.dune-project.org/quality/dune-testtools) | releases/2.6 | Handles system tests
| [doxygen](http://www.stack.nl/~dimitri/doxygen/) | 1.8.13 | Builds documentation
| [METIS](http://glaros.dtc.umn.edu/gkhome/views/metis) | 5 | For parallel runs
| [ParMETIS](http://glaros.dtc.umn.edu/gkhome/views/metis) | 4 | For parallel runs
| [METIS](http://glaros.dtc.umn.edu/gkhome/views/metis) | 5 | For parallel runs
| [ParMETIS](http://glaros.dtc.umn.edu/gkhome/views/metis) | 4 | For parallel runs
#### Step-by-step Instructions
......@@ -119,7 +120,7 @@ If you installed [Anaconda](https://conda.io/docs/user-guide/install/download.ht
git libatlas-base-dev libfftw3-dev libfftw3-mpi-dev \
libfreetype6-dev libhdf5-mpi-dev libmuparser-dev \
libopenmpi-dev libpng-dev libsuperlu-dev libyaml-cpp-dev \
libxft-dev python3-dev python3-pip
libxft-dev python3-dev python3-pip python3-vtk7
**macOS:**
......@@ -302,12 +303,14 @@ or former versions of itself. This ensures that DORiE is running correctly and
producing the expected results. We distinguish _unit tests_ for testing certain
features of the code, and _system tests_ for verifying the results of the
final application. As system tests require executing the DUNE solvers,
it is recommended to build them in a `Release` environment.
it is recommended to build them in a `Release` environment. Additionaly, there
is a set of tests for the Python module.
| Test category | Build tests | Execute tests | Recommended build type |
| ------------- | ----------- | ------------- | ---------------------- |
| Unit tests | `make build_unit_tests` | `make unit_tests` | `Debug`
| System tests | `make build_system_tests` | `make system_tests` | `Release`
| Python tests | N/A | `make test_python` | any
The `make` commands are to be executed from within the `build-cmake` directory.
......
......@@ -49,6 +49,17 @@ Package reference
:members:
:undoc-members:
Generic VTK Reader
------------------
This is a flexible VTK XML file reader which can evaluate all data sets written
by DORiE.
.. autoclass:: dorie.utilities.vtktools.vtkreader.VTKReader
:members:
.. automethod:: __init__
dorie.utilities.check_path
--------------------------
......
......@@ -6,30 +6,35 @@ ARG PROCNUM=1
ARG CC=gcc
ARG CXX=g++
RUN apt-get clean && apt-get update && apt-get install -y \
clang \
cmake \
doxygen \
gcc \
g++ \
gfortran \
git \
libatlas-base-dev \
libfftw3-dev \
libfftw3-mpi-dev \
libfreetype6-dev \
libhdf5-mpi-dev \
libmetis-dev \
libmuparser-dev \
libopenmpi-dev \
libpng-dev \
libparmetis-dev \
libsuperlu-dev \
libxft-dev \
libyaml-cpp-dev \
locales \
python3-dev \
python3-pip \
# disable any prompts while installing packages
ENV DEBIAN_FRONTEND=noninteractive
RUN apt-get clean \
&& apt-get update \
&& apt-get install -y \
clang \
cmake \
doxygen \
gcc \
g++ \
gfortran \
git \
libatlas-base-dev \
libfftw3-dev \
libfftw3-mpi-dev \
libfreetype6-dev \
libhdf5-mpi-dev \
libmetis-dev \
libmuparser-dev \
libopenmpi-dev \
libpng-dev \
libparmetis-dev \
libsuperlu-dev \
libxft-dev \
libyaml-cpp-dev \
locales \
python3-dev \
python3-pip \
python3-vtk7 \
&& apt-get clean
RUN rm -rf /var/lib/apt/lists/* \
......
......@@ -8,3 +8,12 @@ dune_python_install_package(PATH dorie)
# set cache variable
set(DORIE_PYTHON_MODULES ${python_paths} CACHE PATH "Paths to the python modules")
# add Python tests
dune_python_add_test(NAME python_tests
COMMAND ${CMAKE_BINARY_DIR}/run-in-dune-env
python -m pytest -v
${CMAKE_CURRENT_SOURCE_DIR}/dorie/dorie/test/
WORKING_DIRECTORY
${CMAKE_CURRENT_SOURCE_DIR}/dorie/dorie/test/
)
import pytest
import numpy as np
from dorie.utilities.vtktools.vtkreader import VTKReader
VTK_CELL = "vtkreader-cell.vtu"
VTK_VERTEX = "vtkreader-vertex.vtu"
DEFAULT_ARRAY = "head"
ARRAYS = ["head", "flux", "K_0", "Theta", "theta_w"]
# Fixtures ---------------------------------------------------------------
@pytest.fixture
def vtk_cell ():
"""Create a VTKReader for the cell VTK test file"""
return VTKReader(VTK_CELL, DEFAULT_ARRAY)
@pytest.fixture
def vtk_vertex ():
"""Create a VTKReader for the vertex VTK test file"""
return VTKReader(VTK_VERTEX, DEFAULT_ARRAY)
# Tests ------------------------------------------------------------------
def test_init ():
"""Test if initialization works"""
VTKReader(VTK_CELL, DEFAULT_ARRAY)
VTKReader(VTK_VERTEX, DEFAULT_ARRAY)
def test_array_set (vtk_cell, vtk_vertex):
"""Test if all included arrays are found"""
for array in ARRAYS:
vtk_cell.set_data_array(array)
vtk_vertex.set_data_array(array)
with pytest.raises(RuntimeError):
vtk_cell.set_data_array("other")
def test_properties (vtk_cell):
"""Test the VTKReader properties"""
assert np.allclose(vtk_cell.bounds, [[0, 1],
[0, 1],
[0, 0]])
assert vtk_cell.number_of_components == 1
vtk_cell.set_data_array("flux")
assert vtk_cell.number_of_components == 3
def test_evaluate (vtk_cell, vtk_vertex):
"""Test the evaluation function"""
center = np.mean(vtk_cell.bounds, axis=1)
# should return scalar
with pytest.raises(TypeError):
len(vtk_cell.evaluate(center))
# check return of vector
vtk_vertex.set_data_array("flux")
assert len(vtk_vertex.evaluate(center)) == 3
# check values are equal
for array in ["head", "K_0", "flux"]:
vtk_cell.set_data_array(array)
vtk_vertex.set_data_array(array)
assert np.allclose(vtk_cell.evaluate(center),
vtk_vertex.evaluate(center))
# check evaluation errors
# position has wrong shape
with pytest.raises(ValueError):
vtk_vertex.evaluate([0.0, 0.0])
# position out of bounds
pos = vtk_vertex.bounds[:, 1]
pos = pos + [1.0, 0.0, 1.0]
with pytest.raises(ValueError):
vtk_vertex.evaluate(pos)
<?xml version="1.0"?>
<VTKFile type="UnstructuredGrid" version="0.1" byte_order="LittleEndian">
<UnstructuredGrid>
<Piece NumberOfCells="1" NumberOfPoints="4">
<CellData Scalars="head" Vectors="flux">
<DataArray type="Float32" Name="head" NumberOfComponents="1" format="binary">
BAAAAA==AAAAvw==
</DataArray>
<DataArray type="Float32" Name="flux" NumberOfComponents="3" format="binary">
DAAAAA==AAAAgAAAAIAAAAAA
</DataArray>
<DataArray type="Float32" Name="K_0" NumberOfComponents="1" format="binary">
BAAAAA==pIy4N4==
</DataArray>
<DataArray type="Float32" Name="theta_w" NumberOfComponents="1" format="binary">
BAAAAA==Jh4iPh==
</DataArray>
<DataArray type="Float32" Name="Theta" NumberOfComponents="1" format="binary">
BAAAAA==jKPqPq==
</DataArray>
</CellData>
<Points>
<DataArray type="Float32" Name="Coordinates" NumberOfComponents="3" format="binary">
MAAAAA==AAAAAAAAAAAAAAAAAACAPwAAAAAAAAAAAAAAAAAAgD8AAAAAAACAPwAAgD8AAAAA
</DataArray>
</Points>
<Cells>
<DataArray type="Int32" Name="connectivity" NumberOfComponents="1" format="binary">
EAAAAA==AAAAAAEAAAADAAAAAgAAAA==
</DataArray>
<DataArray type="Int32" Name="offsets" NumberOfComponents="1" format="binary">
BAAAAA==BAAAAA==
</DataArray>
<DataArray type="UInt8" Name="types" NumberOfComponents="1" format="binary">
AQAAAA==CQ==
</DataArray>
</Cells>
</Piece>
</UnstructuredGrid>
</VTKFile>
<?xml version="1.0"?>
<VTKFile type="UnstructuredGrid" version="0.1" byte_order="LittleEndian">
<UnstructuredGrid>
<Piece NumberOfCells="1" NumberOfPoints="4">
<PointData Scalars="head" Vectors="flux">
<DataArray type="Float32" Name="head" NumberOfComponents="1" format="binary">
EAAAAA==AAAAAAAAAAAAAIC/AACAvw==
</DataArray>
<DataArray type="Float32" Name="flux" NumberOfComponents="3" format="binary">
MAAAAA==AAAAgAAAAIAAAAAAAAAAgAAAAIAAAAAAAAAAgAAAAIAAAAAAAAAAgAAAAIAAAAAA
</DataArray>
<DataArray type="Float32" Name="K_0" NumberOfComponents="1" format="binary">
EAAAAA==pIy4N6SMuDekjLg3pIy4N4==
</DataArray>
<DataArray type="Float32" Name="theta_w" NumberOfComponents="1" format="binary">
EAAAAA==UriePlK4nj6N0Uo9jdFKPd==
</DataArray>
<DataArray type="Float32" Name="Theta" NumberOfComponents="1" format="binary">
EAAAAA==AACAPwAAgD9Xv449V7+OPb==
</DataArray>
</PointData>
<Points>
<DataArray type="Float32" Name="Coordinates" NumberOfComponents="3" format="binary">
MAAAAA==AAAAAAAAAAAAAAAAAACAPwAAAAAAAAAAAAAAAAAAgD8AAAAAAACAPwAAgD8AAAAA
</DataArray>
</Points>
<Cells>
<DataArray type="Int32" Name="connectivity" NumberOfComponents="1" format="binary">
EAAAAA==AAAAAAEAAAADAAAAAgAAAA==
</DataArray>
<DataArray type="Int32" Name="offsets" NumberOfComponents="1" format="binary">
BAAAAA==BAAAAA==
</DataArray>
<DataArray type="UInt8" Name="types" NumberOfComponents="1" format="binary">
AQAAAA==CQ==
</DataArray>
</Cells>
</Piece>
</UnstructuredGrid>
</VTKFile>
from vtk import vtkXMLGenericDataObjectReader, \
vtkCellTreeLocator, \
vtkDataSet, \
vtkGenericCell, \
reference
import numpy as np
class VTKReader:
"""A generic VTK reader for any type of VTK XML files.
It is intended to evaluate cell or vertex (point) data arrays on any point
of the grid.
Attributes:
_reader: Generic VTK XML data file reader.
_dataset: The complete data set loaded by the reader.
_locator: Fast locator of grid cells based on global coordinates.
_cell: Generic grid cell implementation.
_cell_id: ID of the cell on the grid.
_eval: Evaluation method based on the selected data array type.
_data_array: The actual data to evaluate.
"""
def __init__(self, file_name, array_name, data_type=None):
"""Open the target VTK file, build a cell locator and open the data array.
See :func:`~set_data_array` for details on the arguments regarding the
data array.
Args:
file_name (str): Path to the VTK file.
array_name (str): Name of the data array to evaluate.
data_type (str): Type of the data array.
"""
self._eval = None
self._data_array = None
# set up the reader
self._reader = vtkXMLGenericDataObjectReader()
self._reader.SetFileName(file_name)
self._reader.Update()
self._dataset = vtkDataSet.SafeDownCast(self._reader.GetOutput())
# build the locator
self._locator = vtkCellTreeLocator()
self._locator.SetDataSet(self._dataset)
self._locator.BuildLocator()
# cache for cell and ID
self._cell = vtkGenericCell()
self._cell_id = -1 # invalid cell
# set the evaluation data array
self.set_data_array(array_name, data_type)
def set_data_array (self, array_name, data_type=None):
"""Switch the data array (dataset) to evaluate.
Args:
array_name (str): The name of the data array to evaluate.
data_type (str, optional): The type of the data array (either
``cell`` or ``vertex``). If this is ``None``, vertex arrays
take precedence.
Raises:
RuntimeError: Data array was not found.
"""
self._data_array = None
def retrieve_array (data, name):
"""Return the array if it exists in the data."""
if data.HasArray(name) == 1:
return data.GetArray(name)
else:
return None
if data_type == "cell" or data_type is None:
# retrieve cell data array...
self._data_array = retrieve_array(self._dataset.GetCellData(),
array_name)
self._eval = self.__eval_cell__
if self._data_array is None \
and (data_type == "vertex" or data_type is None):
# ...or vertex data array
self._data_array = retrieve_array(self._dataset.GetPointData(),
array_name)
self._eval = self.__eval_vertex__
if self._data_array is None:
if data_type is None:
raise RuntimeError("Data array '{}' not found in VTK file"
.format(array_name))
else:
raise RuntimeError("Data array '{}' of type '{}' not found in "
"VTK file. Supported types are: 'cell', "
"'vertex'".format(array_name, data_type))
self._num_comp = self._data_array.GetNumberOfComponents()
@property
def bounds (self):
"""Return the outer bounds of the grid.
Returns:
np.array: The rectangular bounds in any spatial direction.
Shape: ``(3, 2)``.
"""
bounds = self._dataset.GetBounds()
return np.reshape(bounds, (3, 2))
@property
def number_of_components (self):
"""Number of components of the evaluated quantity."""
return self._num_comp
def evaluate (self, position):
"""Evaluate the selected data array at a certain position.
Args:
position (list): 3D coordinate probing the data array.
Returns:
np.array: The value of the data array at the given position.
The output is squeezed: Scalar data arrays return scalar values,
vector arrays return numpy arrays with shape ``(3)``.
Raises:
ValueError: ``position`` argument not a 3D list or array.
ValueError: Position outside grid.
RuntimeError: Numerical error in cell locator algorithm.
"""
if len(position) != 3:
raise ValueError("Position must be a 3D lest/array")
# unused stubs
zeros = [0, 0, 0]
zero_ref = reference(0)
# interpolation weights
weights = np.zeros(10)
# check if we cached the correct cell
on_cell = self._cell.EvaluatePosition(position,
zeros,
zero_ref,
zeros,
zero_ref,
weights)
if on_cell == 0 or self._cell_id < 0:
# locate correct cell
self._cell_id = self._locator.FindCell(position,
0,
self._cell,
zeros,
weights)
if self._cell_id < 0:
raise ValueError("Position outside domain: {}".format(position))
elif on_cell < 0:
raise RuntimeError("Error while evaluating if position {} is on "
"cell {}".format(position, self._cell_id))
# evaluate the data
data = self._eval(weights)
return np.squeeze(data)
def __eval_cell__ (self, weights):
"""Evaluation function for a cell-type data array."""
return np.asarray(self._data_array.GetTuple(self._cell_id))
def __eval_vertex__ (self, weights):
"""Evaluation function for a vertex-type data array."""
data = np.zeros(self._num_comp)
for point in range(self._cell.GetNumberOfPoints()):
point_id = self._cell.GetPointId(point)
tup = np.asarray(self._data_array.GetTuple(point_id))
data += weights[point] * tup
return data
......@@ -9,7 +9,8 @@ setup(name='dorie',
author='Lukas Riedel',
author_email='dorieteam@iup.uni-heidelberg.de',
url='https://ts-gitlab.iup.uni-heidelberg.de/dorie/dorie',
packages=find_packages(),
packages=find_packages(exclude=["dorie.test"]),
test_suite='py.test',
install_requires= ['cycler',
'configparser',
'h5py',
......@@ -18,11 +19,13 @@ setup(name='dorie',
'Pillow',
'pygmsh',
'pyparsing',
'pytest',
'python-dateutil',
'pytz',
'pyyaml',
'scipy',
'six',
'vtk',
'wheel'
],
scripts=[
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment