Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found

Target

Select target project
No results found
Show changes
Commits on Source (65)
Showing
with 1459 additions and 473 deletions
image: quantumtinkerer/research
image: gitlab.kwant-project.org:5005/qt/adaptive:latest
stages:
- prebuild
- test
build docker:
stage: prebuild
image:
name: gcr.io/kaniko-project/executor:debug
entrypoint: [""]
before_script:
- echo "{\"auths\":{\"$CI_REGISTRY\":{\"username\":\"$CI_REGISTRY_USER\",\"password\":\"$CI_REGISTRY_PASSWORD\"}}}" > /kaniko/.docker/config.json
- export CI_REF=${CI_COMMIT_TAG:-latest}
script:
- /kaniko/executor
--context $CI_PROJECT_DIR
--dockerfile $CI_PROJECT_DIR/Dockerfile
--destination $CI_REGISTRY_IMAGE:$CI_REF
only:
changes:
- Dockerfile
- environment.yml
- test-requirements.txt
test:
stage: test
script:
- pip install -r test-requirements.txt
- conda install -y -q scikit-optimize
- py.test --verbose --cov=adaptive --cov-report term --cov-report html adaptive
artifacts:
paths:
- htmlcov
authors check:
stage: test
script:
- MISSING_AUTHORS=$(git shortlog -s HEAD | sed -e "s/^[0-9\t ]*//"| xargs -i sh -c 'grep -q "{}" AUTHORS.md || echo "{} missing from authors"')
- if [ ! -z "$MISSING_AUTHORS" ]; then { echo $MISSING_AUTHORS; exit 1; }; fi
allow_failure: true
check whitespace style:
stage: test
script: ./check_whitespace
allow_failure: true
# Adaptive Authors
## Authors
Below is a list of the contributors to Adaptive:
+ [Anton Akhmerov](<https://antonakhmerov.org>)
......@@ -6,9 +6,3 @@ Below is a list of the contributors to Adaptive:
+ [Christoph Groth](<http://inac.cea.fr/Pisp/christoph.groth/>)
+ Jorn Hoofwijk
+ [Joseph Weston](<https://joseph.weston.cloud>)
For a full list of contributors run
```
git log --pretty=format:"%an" | sort | uniq
```
# Docker image for testing adaptive
FROM conda/miniconda3:latest
# make our environment sane
ENV LANG C.UTF-8
ENV LC_ALL C.UTF-8
# install git
RUN apt-get update && \
apt-get install -y --no-install-recommends git && \
apt-get clean && \
rm -rf /var/lib/apt/lists/*
COPY environment.yml test-requirements.txt /
RUN conda env update --quiet -n root -f environment.yml
RUN conda clean --yes --all
RUN pip install -r test-requirements.txt
# ![][logo] adaptive
[![PyPI](https://img.shields.io/pypi/v/adaptive.svg)](https://pypi.python.org/pypi/adaptive)
[![Conda](https://anaconda.org/conda-forge/adaptive/badges/installer/conda.svg)](https://anaconda.org/conda-forge/adaptive)
[![Downloads](https://anaconda.org/conda-forge/adaptive/badges/downloads.svg)](https://anaconda.org/conda-forge/adaptive)
[![pipeline status](https://gitlab.kwant-project.org/qt/adaptive/badges/master/pipeline.svg)](https://gitlab.kwant-project.org/qt/adaptive/pipelines)
[![DOI](https://zenodo.org/badge/113714660.svg)](https://zenodo.org/badge/latestdoi/113714660)
[![Binder](https://mybinder.org/badge.svg)](https://mybinder.org/v2/gh/python-adaptive/adaptive/master?filepath=learner.ipynb)
[![Join the chat at https://gitter.im/python-adaptive/adaptive](https://img.shields.io/gitter/room/nwjs/nw.js.svg)](https://gitter.im/python-adaptive/adaptive)
**Tools for adaptive parallel sampling of mathematical functions.**
`adaptive` is an [open-source](LICENSE) Python library designed to make adaptive parallel function evaluation simple.
With `adaptive` you just supply a function with its bounds, and it will be evaluated at the "best" points in parameter space.
With just a few lines of code you can evaluate functions on a computing cluster, live-plot the data as it returns, and fine-tune the adaptive sampling algorithm.
Check out the `adaptive` [example notebook `learner.ipynb`](learner.ipynb) (or run it [live on Binder](https://mybinder.org/v2/gh/python-adaptive/adaptive/master?filepath=learner.ipynb)) to see examples of how to use `adaptive`.
**WARNING: `adaptive` is still in a beta development stage**
## Implemented algorithms
The core concept in `adaptive` is that of a *learner*. A *learner* samples
a function at the best places in its parameter space to get maximum
"information" about the function. As it evaluates the function
at more and more points in the parameter space, it gets a better idea of where
the best places are to sample next.
Of course, what qualifies as the "best places" will depend on your application domain!
`adaptive` makes some reasonable default choices, but the details of the adaptive
sampling are completely customizable.
The following learners are implemented:
* `Learner1D`, for 1D functions `f: ℝ → ℝ^N`,
* `Learner2D`, for 2D functions `f: ℝ^2 → ℝ^N`,
* `LearnerND`, for ND functions `f: ℝ^N → ℝ^M`,
* `AverageLearner`, For stochastic functions where you want to average the result over many evaluations,
* `IntegratorLearner`, for when you want to intergrate a 1D function `f: ℝ → ℝ`,
* `BalancingLearner`, for when you want to run several learners at once, selecting the "best" one each time you get more points.
In addition to the learners, `adaptive` also provides primitives for running
the sampling across several cores and even several machines, with built-in support
for [`concurrent.futures`](https://docs.python.org/3/library/concurrent.futures.html),
[`ipyparallel`](https://ipyparallel.readthedocs.io/en/latest/)
and [`distributed`](https://distributed.readthedocs.io/en/latest/).
## Examples
<img src="https://user-images.githubusercontent.com/6897215/38739170-6ac7c014-3f34-11e8-9e8f-93b3a3a3d61b.gif" width='20%'> </img>
<img src="https://user-images.githubusercontent.com/6897215/35219611-ac8b2122-ff73-11e7-9332-adffab64a8ce.gif" width='40%'> </img>
## Installation
`adaptive` works with Python 3.6 and higher on Linux, Windows, or Mac, and provides optional extensions for working with the Jupyter/IPython Notebook.
The recommended way to install adaptive is using `conda`:
```bash
conda install -c conda-forge adaptive
```
`adaptive` is also available on PyPI:
```bash
pip install adaptive[notebook]
```
The `[notebook]` above will also install the optional dependencies for running `adaptive` inside
a Jupyter notebook.
## Development
Clone the repository and run `setup.py develop` to add a link to the cloned repo into your
Python path:
```
git clone git@github.com:python-adaptive/adaptive.git
cd adaptive
python3 setup.py develop
```
We highly recommend using a Conda environment or a virtualenv to manage the versions of your installed
packages while working on `adaptive`.
In order to not pollute the history with the output of the notebooks, please setup the git filter by executing
```bash
python ipynb_filter.py
```
in the repository.
## Credits
We would like to give credits to the following people:
- Pedro Gonnet for his implementation of [`CQUAD`](https://www.gnu.org/software/gsl/manual/html_node/CQUAD-doubly_002dadaptive-integration.html), "Algorithm 4" as described in "Increasing the Reliability of Adaptive Quadrature Using Explicit Interpolants", P. Gonnet, ACM Transactions on Mathematical Software, 37 (3), art. no. 26, 2010.
- Pauli Virtanen for his `AdaptiveTriSampling` script (no longer available online since SciPy Central went down) which served as inspiration for the [`Learner2D`](adaptive/learner/learner2D.py).
For general discussion, we have a [Gitter chat channel](https://gitter.im/python-adaptive/adaptive). If you find any bugs or have any feature suggestions please file a GitLab [issue](https://gitlab.kwant-project.org/qt/adaptive/issues/new?issue) or submit a [merge request](https://gitlab.kwant-project.org/qt/adaptive/merge_requests).
[logo]: https://gitlab.kwant-project.org/qt/adaptive/uploads/d20444093920a4a0499e165b5061d952/logo.png "adaptive logo"
.. summary-start
|logo| adaptive
===============
|PyPI| |Conda| |Downloads| |Pipeline status| |DOI| |Binder| |Gitter|
|Documentation| |GitHub|
**Tools for adaptive parallel sampling of mathematical functions.**
``adaptive`` is an open-source Python library designed to
make adaptive parallel function evaluation simple. With ``adaptive`` you
just supply a function with its bounds, and it will be evaluated at the
“best” points in parameter space. With just a few lines of code you can
evaluate functions on a computing cluster, live-plot the data as it
returns, and fine-tune the adaptive sampling algorithm.
Run the ``adaptive`` example notebook `live on
Binder <https://mybinder.org/v2/gh/python-adaptive/adaptive/master?filepath=learner.ipynb>`_
to see examples of how to use ``adaptive`` or visit the
`tutorial on Read the Docs <https://adaptive.readthedocs.io/en/latest/tutorial/tutorial.html>`__.
.. summary-end
**WARNING: adaptive is still in a beta development stage**
.. not-in-documentation-start
Implemented algorithms
----------------------
The core concept in ``adaptive`` is that of a *learner*. A *learner*
samples a function at the best places in its parameter space to get
maximum “information” about the function. As it evaluates the function
at more and more points in the parameter space, it gets a better idea of
where the best places are to sample next.
Of course, what qualifies as the “best places” will depend on your
application domain! ``adaptive`` makes some reasonable default choices,
but the details of the adaptive sampling are completely customizable.
The following learners are implemented:
- ``Learner1D``, for 1D functions ``f: ℝ → ℝ^N``,
- ``Learner2D``, for 2D functions ``f: ℝ^2 → ℝ^N``,
- ``LearnerND``, for ND functions ``f: ℝ^N → ℝ^M``,
- ``AverageLearner``, For stochastic functions where you want to
average the result over many evaluations,
- ``IntegratorLearner``, for
when you want to intergrate a 1D function ``f: ℝ → ℝ``,
- ``BalancingLearner``, for when you want to run several learners at once,
selecting the “best” one each time you get more points.
In addition to the learners, ``adaptive`` also provides primitives for
running the sampling across several cores and even several machines,
with built-in support for
`concurrent.futures <https://docs.python.org/3/library/concurrent.futures.html>`_,
`ipyparallel <https://ipyparallel.readthedocs.io/en/latest/>`_ and
`distributed <https://distributed.readthedocs.io/en/latest/>`_.
Examples
--------
.. raw:: html
<img src="https://user-images.githubusercontent.com/6897215/38739170-6ac7c014-3f34-11e8-9e8f-93b3a3a3d61b.gif" width='20%'> </img> <img src="https://user-images.githubusercontent.com/6897215/35219611-ac8b2122-ff73-11e7-9332-adffab64a8ce.gif" width='40%'> </img> <img src="https://user-images.githubusercontent.com/6897215/47256441-d6d53700-d480-11e8-8224-d1cc49dbdcf5.gif" width='20%'> </img>
.. not-in-documentation-end
Installation
------------
``adaptive`` works with Python 3.6 and higher on Linux, Windows, or Mac,
and provides optional extensions for working with the Jupyter/IPython
Notebook.
The recommended way to install adaptive is using ``conda``:
.. code:: bash
conda install -c conda-forge adaptive
``adaptive`` is also available on PyPI:
.. code:: bash
pip install adaptive[notebook]
The ``[notebook]`` above will also install the optional dependencies for
running ``adaptive`` inside a Jupyter notebook.
Development
-----------
Clone the repository and run ``setup.py develop`` to add a link to the
cloned repo into your Python path:
.. code:: bash
git clone git@github.com:python-adaptive/adaptive.git
cd adaptive
python3 setup.py develop
We highly recommend using a Conda environment or a virtualenv to manage
the versions of your installed packages while working on ``adaptive``.
In order to not pollute the history with the output of the notebooks,
please setup the git filter by executing
.. code:: bash
python ipynb_filter.py
in the repository.
Credits
-------
We would like to give credits to the following people:
- Pedro Gonnet for his implementation of `CQUAD <https://www.gnu.org/software/gsl/manual/html_node/CQUAD-doubly_002dadaptive-integration.html>`_,
“Algorithm 4” as described in “Increasing the Reliability of Adaptive
Quadrature Using Explicit Interpolants”, P. Gonnet, ACM Transactions on
Mathematical Software, 37 (3), art. no. 26, 2010.
- Pauli Virtanen for his ``AdaptiveTriSampling`` script (no longer
available online since SciPy Central went down) which served as
inspiration for the `~adaptive.Learner2D`.
.. credits-end
For general discussion, we have a `Gitter chat
channel <https://gitter.im/python-adaptive/adaptive>`_. If you find any
bugs or have any feature suggestions please file a GitLab
`issue <https://gitlab.kwant-project.org/qt/adaptive/issues/new?issue>`_
or submit a `merge
request <https://gitlab.kwant-project.org/qt/adaptive/merge_requests>`_.
.. references-start
.. |logo| image:: https://adaptive.readthedocs.io/en/latest/_static/logo.png
.. |PyPI| image:: https://img.shields.io/pypi/v/adaptive.svg
:target: https://pypi.python.org/pypi/adaptive
.. |Conda| image:: https://img.shields.io/badge/install%20with-conda-green.svg
:target: https://anaconda.org/conda-forge/adaptive
.. |Downloads| image:: https://img.shields.io/conda/dn/conda-forge/adaptive.svg
:target: https://anaconda.org/conda-forge/adaptive
.. |Pipeline status| image:: https://gitlab.kwant-project.org/qt/adaptive/badges/master/pipeline.svg
:target: https://gitlab.kwant-project.org/qt/adaptive/pipelines
.. |DOI| image:: https://img.shields.io/badge/doi-10.5281%2Fzenodo.1182437-blue.svg
:target: https://doi.org/10.5281/zenodo.1182437
.. |Binder| image:: https://mybinder.org/badge.svg
:target: https://mybinder.org/v2/gh/python-adaptive/adaptive/master?filepath=learner.ipynb
.. |Gitter| image:: https://img.shields.io/gitter/room/nwjs/nw.js.svg
:target: https://gitter.im/python-adaptive/adaptive
.. |Documentation| image:: https://readthedocs.org/projects/adaptive/badge/?version=latest
:target: https://adaptive.readthedocs.io/en/latest/?badge=latest
.. |GitHub| image:: https://img.shields.io/github/stars/python-adaptive/adaptive.svg?style=social
:target: https://github.com/python-adaptive/adaptive/stargazers
.. references-end
......@@ -8,18 +8,17 @@ from . import learner
from . import runner
from . import utils
from .learner import (Learner1D, Learner2D, LearnerND, AverageLearner,
BalancingLearner, make_datasaver, DataSaver,
IntegratorLearner)
from .learner import (BaseLearner, Learner1D, Learner2D, LearnerND,
AverageLearner, BalancingLearner, make_datasaver,
DataSaver, IntegratorLearner)
with suppress(ImportError):
# Only available if 'scikit-optimize' is installed
from .learner import SKOptLearner
from .runner import Runner, BlockingRunner
from . import version
from .runner import Runner, AsyncRunner, BlockingRunner
__version__ = version.version
from ._version import __version__
del _version
del notebook_integration # to avoid confusion with `notebook_extension`
del version
# This file is part of 'miniver': https://github.com/jbweston/miniver
#
# This file will be overwritten by setup.py when a source or binary
# distribution is made. The magic value "__use_git__" is interpreted by
# version.py.
# _version.py.
version = "__use_git__"
# These values are only set if the distribution was created with 'git archive'
refnames = "$Format:%D$"
refnames = "$Format:%D$"
git_hash = "$Format:%h$"
# -*- coding: utf-8 -*-
# This file is part of 'miniver': https://github.com/jbweston/miniver
#
from collections import namedtuple
import os
import subprocess
import sys
from distutils.command.build import build as build_orig
from distutils.command.build_py import build_py as build_py_orig
from setuptools.command.sdist import sdist as sdist_orig
Version = namedtuple('Version', ('release', 'dev', 'labels'))
......@@ -19,31 +21,39 @@ STATIC_VERSION_FILE = '_static_version.py'
def get_version(version_file=STATIC_VERSION_FILE):
version_info = {}
with open(os.path.join(package_root, version_file), 'rb') as f:
exec(f.read(), {}, version_info)
version_info = get_static_version_info(version_file)
version = version_info['version']
version_is_from_git = (version == "__use_git__")
if version_is_from_git:
if version == "__use_git__":
version = get_version_from_git()
if not version:
version = get_version_from_git_archive(version_info)
if not version:
version = Version("unknown", None, None)
return semver_format(version)
return pep440_format(version)
else:
return version
def semver_format(version_info):
def get_static_version_info(version_file=STATIC_VERSION_FILE):
version_info = {}
with open(os.path.join(package_root, version_file), 'rb') as f:
exec(f.read(), {}, version_info)
return version_info
def version_is_from_git(version_file=STATIC_VERSION_FILE):
return get_static_version_info(version_file)['version'] == '__use_git__'
def pep440_format(version_info):
release, dev, labels = version_info
version_parts = [release]
if dev:
if release.endswith('-dev'):
if release.endswith('-dev') or release.endswith('.dev'):
version_parts.append(dev)
else:
version_parts.append('-dev{}'.format(dev))
else: # prefer PEP440 over strict adhesion to semver
version_parts.append('.dev{}'.format(dev))
if labels:
version_parts.append('+')
......@@ -64,26 +74,42 @@ def get_version_from_git():
if not os.path.samefile(p.communicate()[0].decode().rstrip('\n'),
distr_root):
# The top-level directory of the current Git repository is not the same
# as the root directory of the Kwant distribution: do not extract the
# as the root directory of the distribution: do not extract the
# version from Git.
return
# git describe --first-parent does not take into account tags from branches
# that were merged-in.
# that were merged-in. The '--long' flag gets us the 'dev' version and
# git hash, '--always' returns the git hash even if there are no tags.
for opts in [['--first-parent'], []]:
try:
p = subprocess.Popen(['git', 'describe', '--long'] + opts,
cwd=distr_root,
stdout=subprocess.PIPE, stderr=subprocess.PIPE)
p = subprocess.Popen(
['git', 'describe', '--long', '--always'] + opts,
cwd=distr_root,
stdout=subprocess.PIPE, stderr=subprocess.PIPE)
except OSError:
return
if p.wait() == 0:
break
else:
return
description = p.communicate()[0].decode().strip('v').rstrip('\n')
release, dev, git = description.rsplit('-', 2)
description = (
p.communicate()[0]
.decode()
.strip('v') # Tags can have a leading 'v', but the version should not
.rstrip('\n')
.rsplit('-', 2) # Split the latest tag, commits since tag, and hash
)
try:
release, dev, git = description
except ValueError: # No tags, only the git hash
# prepend 'g' to match with format returned by 'git describe'
git = 'g{}'.format(*description)
release = 'unknown'
dev = None
labels = []
if dev == "0":
dev = None
......@@ -93,7 +119,7 @@ def get_version_from_git():
try:
p = subprocess.Popen(['git', 'diff', '--quiet'], cwd=distr_root)
except OSError:
labels.append('confused') # This should never happen.
labels.append('confused') # This should never happen.
else:
if p.wait() == 1:
labels.append('dirty')
......@@ -126,14 +152,15 @@ def get_version_from_git_archive(version_info):
release, *_ = sorted(version_tags) # prefer e.g. "2.0" over "2.0rc1"
return Version(release, dev=None, labels=None)
else:
return Version('unknown', dev=None, labels=[f'g{git_hash}'])
return Version('unknown', dev=None, labels=['g{}'.format(git_hash)])
__version__ = get_version()
version = get_version()
# The following section defines a module global 'cmdclass',
# which can be used from setup.py. The 'package_name' and
# 'version' module globals are used (but not modified).
# '__version__' module globals are used (but not modified).
def _write_version(fname):
# This could be a hard link, so try to delete it first. Is there any way
......@@ -144,10 +171,10 @@ def _write_version(fname):
pass
with open(fname, 'w') as f:
f.write("# This file has been created by setup.py.\n"
"version = '{}'\n".format(version))
"version = '{}'\n".format(__version__))
class _build(build_orig):
class _build_py(build_py_orig):
def run(self):
super().run()
_write_version(os.path.join(self.build_lib, package_name,
......@@ -161,4 +188,4 @@ class _sdist(sdist_orig):
STATIC_VERSION_FILE))
cmdclass = dict(sdist=_sdist, build=_build)
cmdclass = dict(sdist=_sdist, build_py=_build_py)
......@@ -17,9 +17,9 @@ class AverageLearner(BaseLearner):
Parameters
----------
atol : float
Desired absolute tolerance
Desired absolute tolerance.
rtol : float
Desired relative tolerance
Desired relative tolerance.
Attributes
----------
......@@ -27,6 +27,8 @@ class AverageLearner(BaseLearner):
Sampled points and values.
pending_points : set
Points that still have to be evaluated.
npoints : int
Number of evaluated points.
"""
def __init__(self, function, atol=None, rtol=None):
......@@ -48,7 +50,7 @@ class AverageLearner(BaseLearner):
@property
def n_requested(self):
return len(self.data) + len(self.pending_points)
return self.npoints + len(self.pending_points)
def ask(self, n, tell_pending=True):
points = list(range(self.n_requested, self.n_requested + n))
......@@ -59,7 +61,7 @@ class AverageLearner(BaseLearner):
- set(self.data)
- set(self.pending_points))[:n]
loss_improvements = [self.loss_improvement(n) / n] * n
loss_improvements = [self._loss_improvement(n) / n] * n
if tell_pending:
for p in points:
self.tell_pending(p)
......@@ -81,10 +83,13 @@ class AverageLearner(BaseLearner):
@property
def mean(self):
"""The average of all values in `data`."""
return self.sum_f / self.npoints
@property
def std(self):
"""The corrected sample standard deviation of the values
in `data`."""
n = self.npoints
if n < 2:
return np.inf
......@@ -106,7 +111,7 @@ class AverageLearner(BaseLearner):
return max(standard_error / self.atol,
standard_error / abs(self.mean) / self.rtol)
def loss_improvement(self, n):
def _loss_improvement(self, n):
loss = self.loss()
if np.isfinite(loss):
return loss - self.loss(n=self.npoints + n)
......@@ -118,6 +123,12 @@ class AverageLearner(BaseLearner):
self.pending_points = set()
def plot(self):
"""Returns a histogram of the evaluated data.
Returns
-------
holoviews.element.Histogram
A histogram of the evaluated data."""
hv = ensure_holoviews()
vals = [v for v in self.data.values() if v is not None]
if not vals:
......@@ -125,3 +136,9 @@ class AverageLearner(BaseLearner):
num_bins = int(max(5, sqrt(self.npoints)))
vals = hv.Points(vals)
return hv.operation.histogram(vals, num_bins=num_bins, dimension=1)
def _get_data(self):
return (self.data, self.npoints, self.sum_f, self.sum_f_sq)
def _set_data(self, data):
self.data, self.npoints, self.sum_f, self.sum_f_sq = data
......@@ -3,6 +3,7 @@ from collections import defaultdict
from contextlib import suppress
from functools import partial
from operator import itemgetter
import os.path
import numpy as np
......@@ -21,37 +22,49 @@ class BalancingLearner(BaseLearner):
Parameters
----------
learners : sequence of BaseLearner
learners : sequence of `~adaptive.BaseLearner`\s
The learners from which to choose. These must all have the same type.
cdims : sequence of dicts, or (keys, iterable of values), optional
Constant dimensions; the parameters that label the learners. Used
in `plot`.
Example inputs that all give identical results:
- sequence of dicts:
>>> cdims = [{'A': True, 'B': 0},
... {'A': True, 'B': 1},
... {'A': False, 'B': 0},
... {'A': False, 'B': 1}]`
- tuple with (keys, iterable of values):
>>> cdims = (['A', 'B'], itertools.product([True, False], [0, 1]))
>>> cdims = (['A', 'B'], [(True, 0), (True, 1),
... (False, 0), (False, 1)])
Attributes
----------
learners : list
The sequence of `~adaptive.BaseLearner`\s.
function : callable
A function that calls the functions of the underlying learners.
Its signature is ``function(learner_index, point)``.
strategy : 'loss_improvements' (default), 'loss', or 'npoints'
The points that the 'BalancingLearner' choses can be either based on:
The points that the `BalancingLearner` choses can be either based on:
the best 'loss_improvements', the smallest total 'loss' of the
child learners, or the number of points per learner, using 'npoints'.
One can dynamically change the strategy while the simulation is
running by changing the 'learner.strategy' attribute.
running by changing the ``learner.strategy`` attribute.
Notes
-----
This learner compares the 'loss' calculated from the "child" learners.
This learner compares the `loss` calculated from the "child" learners.
This requires that the 'loss' from different learners *can be meaningfully
compared*. For the moment we enforce this restriction by requiring that
all learners are the same type but (depending on the internals of the
learner) it may be that the loss cannot be compared *even between learners
of the same type*. In this case the BalancingLearner will behave in an
undefined way.
of the same type*. In this case the `~adaptive.BalancingLearner` will
behave in an undefined way. Change the `strategy` in that case.
"""
def __init__(self, learners, *, cdims=None, strategy='loss_improvements'):
......@@ -75,6 +88,12 @@ class BalancingLearner(BaseLearner):
@property
def strategy(self):
"""Can be either 'loss_improvements' (default), 'loss', or 'npoints'
The points that the `BalancingLearner` choses can be either based on:
the best 'loss_improvements', the smallest total 'loss' of the
child learners, or the number of points per learner, using 'npoints'.
One can dynamically change the strategy while the simulation is
running by changing the ``learner.strategy`` attribute."""
return self._strategy
@strategy.setter
......@@ -116,7 +135,7 @@ class BalancingLearner(BaseLearner):
points = []
loss_improvements = []
for _ in range(n):
losses = self.losses(real=False)
losses = self._losses(real=False)
max_ind = np.argmax(losses)
xs, ls = self.learners[max_ind].ask(1)
points.append((max_ind, xs[0]))
......@@ -159,7 +178,7 @@ class BalancingLearner(BaseLearner):
self._loss.pop(index, None)
self.learners[index].tell_pending(x)
def losses(self, real=True):
def _losses(self, real=True):
losses = []
loss_dict = self._loss if real else self._pending_loss
......@@ -172,7 +191,7 @@ class BalancingLearner(BaseLearner):
@cache_latest
def loss(self, real=True):
losses = self.losses(real)
losses = self._losses(real)
return max(losses)
def plot(self, cdims=None, plotter=None, dynamic=True):
......@@ -183,28 +202,34 @@ class BalancingLearner(BaseLearner):
cdims : sequence of dicts, or (keys, iterable of values), optional
Constant dimensions; the parameters that label the learners.
Example inputs that all give identical results:
- sequence of dicts:
>>> cdims = [{'A': True, 'B': 0},
... {'A': True, 'B': 1},
... {'A': False, 'B': 0},
... {'A': False, 'B': 1}]`
- tuple with (keys, iterable of values):
>>> cdims = (['A', 'B'], itertools.product([True, False], [0, 1]))
>>> cdims = (['A', 'B'], [(True, 0), (True, 1),
... (False, 0), (False, 1)])
plotter : callable, optional
A function that takes the learner as a argument and returns a
holoviews object. By default learner.plot() will be called.
holoviews object. By default ``learner.plot()`` will be called.
dynamic : bool, default True
Return a holoviews.DynamicMap if True, else a holoviews.HoloMap.
The DynamicMap is rendered as the sliders change and can therefore
not be exported to html. The HoloMap does not have this problem.
Return a `holoviews.core.DynamicMap` if True, else a
`holoviews.core.HoloMap`. The `~holoviews.core.DynamicMap` is
rendered as the sliders change and can therefore not be exported
to html. The `~holoviews.core.HoloMap` does not have this problem.
Returns
-------
dm : holoviews.DynamicMap object (default) or holoviews.HoloMap object
A DynamicMap (dynamic=True) or HoloMap (dynamic=False) with
sliders that are defined by 'cdims'.
dm : `holoviews.core.DynamicMap` (default) or `holoviews.core.HoloMap`
A `DynamicMap` ``(dynamic=True)`` or `HoloMap`
``(dynamic=False)`` with sliders that are defined by `cdims`.
"""
hv = ensure_holoviews()
cdims = cdims or self._cdims_default
......@@ -248,13 +273,13 @@ class BalancingLearner(BaseLearner):
def from_product(cls, f, learner_type, learner_kwargs, combos):
"""Create a `BalancingLearner` with learners of all combinations of
named variables’ values. The `cdims` will be set correctly, so calling
`learner.plot` will be a `holoviews.HoloMap` with the correct labels.
`learner.plot` will be a `holoviews.core.HoloMap` with the correct labels.
Parameters
----------
f : callable
Function to learn, must take arguments provided in in `combos`.
learner_type : BaseLearner
learner_type : `BaseLearner`
The learner that should wrap the function. For example `Learner1D`.
learner_kwargs : dict
Keyword argument for the `learner_type`. For example `dict(bounds=[0, 1])`.
......@@ -283,7 +308,7 @@ class BalancingLearner(BaseLearner):
Notes
-----
The order of the child learners inside `learner.learners` is the same
as `adaptive.utils.named_product(**combos)`.
as ``adaptive.utils.named_product(**combos)``.
"""
learners = []
arguments = named_product(**combos)
......@@ -291,3 +316,75 @@ class BalancingLearner(BaseLearner):
learner = learner_type(function=partial(f, **combo), **learner_kwargs)
learners.append(learner)
return cls(learners, cdims=arguments)
def save(self, folder, compress=True):
"""Save the data of the child learners into pickle files
in a directory.
Parameters
----------
folder : str
Directory in which the learners's data will be saved.
compress : bool, default True
Compress the data upon saving using `gzip`. When saving
using compression, one must load it with compression too.
Notes
-----
The child learners need to have a 'fname' attribute in order to use
this method.
Example
-------
>>> def combo_fname(val):
... return '__'.join([f'{k}_{v}.p' for k, v in val.items()])
...
... def f(x, a, b): return a * x**2 + b
...
>>> learners = []
>>> for combo in adaptive.utils.named_product(a=[1, 2], b=[1]):
... l = Learner1D(functools.partial(f, combo=combo))
... l.fname = combo_fname(combo) # 'a_1__b_1.p', 'a_2__b_1.p' etc.
... learners.append(l)
... learner = BalancingLearner(learners)
... # Run the learner
... runner = adaptive.Runner(learner)
... # Then save
... learner.save('data_folder') # use 'load' in the same way
"""
if len(self.learners) != len(set(l.fname for l in self.learners)):
raise RuntimeError("The 'learner.fname's are not all unique.")
for l in self.learners:
l.save(os.path.join(folder, l.fname), compress=compress)
def load(self, folder, compress=True):
"""Load the data of the child learners from pickle files
in a directory.
Parameters
----------
folder : str
Directory from which the learners's data will be loaded.
compress : bool, default True
If the data is compressed when saved, one must load it
with compression too.
Notes
-----
The child learners need to have a 'fname' attribute in order to use
this method.
Example
-------
See the example in the `BalancingLearner.save` doc-string.
"""
for l in self.learners:
l.load(os.path.join(folder, l.fname), compress=compress)
def _get_data(self):
return [l._get_data() for l in learner.learners]
def _set_data(self, data):
for l, _data in zip(self.learners, data):
l._set_data(_data)
# -*- coding: utf-8 -*-
import abc
from contextlib import suppress
from copy import deepcopy
from ..utils import save, load
class BaseLearner(metaclass=abc.ABCMeta):
"""Base class for algorithms for learning a function 'f: X → Y'.
......@@ -11,14 +14,19 @@ class BaseLearner(metaclass=abc.ABCMeta):
function : callable: X → Y
The function to learn.
data : dict: X → Y
'function' evaluated at certain points.
`function` evaluated at certain points.
The values can be 'None', which indicates that the point
will be evaluated, but that we do not have the result yet.
npoints : int, optional
The number of evaluated points that have been added to the learner.
Subclasses do not *have* to implement this attribute.
pending_points : set, optional
Points that have been requested but have not been evaluated yet.
Subclasses do not *have* to implement this attribute.
Subclasses may define a 'plot' method that takes no parameters
Notes
-----
Subclasses may define a ``plot`` method that takes no parameters
and returns a holoviews plot.
"""
......@@ -75,15 +83,94 @@ class BaseLearner(metaclass=abc.ABCMeta):
n : int
The number of points to choose.
tell_pending : bool, default: True
If True, add the chosen points to this
learner's 'data' with 'None' for the 'y'
values. Set this to False if you do not
If True, add the chosen points to this learner's
`pending_points`. Set this to False if you do not
want to modify the state of the learner.
"""
pass
@abc.abstractmethod
def _get_data(self):
pass
@abc.abstractmethod
def _set_data(self):
pass
def copy_from(self, other):
"""Copy over the data from another learner.
Parameters
----------
other : BaseLearner object
The learner from which the data is copied.
"""
self._set_data(other._get_data())
def save(self, fname=None, compress=True):
"""Save the data of the learner into a pickle file.
Parameters
----------
fname : str, optional
The filename of the learner's pickle data file. If None use
the 'fname' attribute, like 'learner.fname = "example.p".
compress : bool, default True
Compress the data upon saving using 'gzip'. When saving
using compression, one must load it with compression too.
Notes
-----
There are **two ways** of naming the files:
1. Using the ``fname`` argument in ``learner.save(fname='example.p')``
2. Setting the ``fname`` attribute, like
``learner.fname = "data/example.p"`` and then ``learner.save()``.
"""
fname = fname or self.fname
data = self._get_data()
save(fname, data, compress)
def load(self, fname=None, compress=True):
"""Load the data of a learner from a pickle file.
Parameters
----------
fname : str, optional
The filename of the saved learner's pickled data file.
If None use the 'fname' attribute, like
'learner.fname = "example.p".
compress : bool, default True
If the data is compressed when saved, one must load it
with compression too.
Notes
-----
See the notes in the `save` doc-string.
"""
fname = fname or self.fname
with suppress(FileNotFoundError, EOFError):
data = load(fname, compress)
self._set_data(data)
def __getstate__(self):
return deepcopy(self.__dict__)
def __setstate__(self, state):
self.__dict__ = state
@property
def fname(self):
"""Filename for the learner when it is saved (or loaded) using
`~adaptive.BaseLearner.save` (or `~adaptive.BaseLearner.load` ).
"""
# This is a property because then it will be availible in the DataSaver
try:
return self._fname
except AttributeError:
raise AttributeError("Set 'learner.fname' or use the 'fname'"
" argument when using 'learner.save' or 'learner.load'.")
@fname.setter
def fname(self, fname):
self._fname = fname
......@@ -2,13 +2,16 @@
from collections import OrderedDict
import functools
from .base_learner import BaseLearner
from ..utils import copy_docstring_from
class DataSaver:
"""Save extra data associated with the values that need to be learned.
Parameters
----------
learner : Learner object
learner : `~adaptive.BaseLearner` instance
The learner that needs to be wrapped.
arg_picker : function
Function that returns the argument that needs to be learned.
......@@ -16,10 +19,11 @@ class DataSaver:
Example
-------
Imagine we have a function that returns a dictionary
of the form: `{'y': y, 'err_est': err_est}`.
of the form: ``{'y': y, 'err_est': err_est}``.
>>> from operator import itemgetter
>>> _learner = Learner1D(f, bounds=(-1.0, 1.0))
>>> learner = DataSaver(_learner, arg_picker=operator.itemgetter('y'))
>>> learner = DataSaver(_learner, arg_picker=itemgetter('y'))
"""
def __init__(self, learner, arg_picker):
......@@ -31,14 +35,35 @@ class DataSaver:
def __getattr__(self, attr):
return getattr(self.learner, attr)
@copy_docstring_from(BaseLearner.tell)
def tell(self, x, result):
y = self.arg_picker(result)
self.extra_data[x] = result
self.learner.tell(x, y)
@copy_docstring_from(BaseLearner.tell_pending)
def tell_pending(self, x):
self.learner.tell_pending(x)
def _get_data(self):
return self.learner._get_data(), self.extra_data
def _set_data(self, data):
learner_data, self.extra_data = data
self.learner._set_data(learner_data)
@copy_docstring_from(BaseLearner.save)
def save(self, fname=None, compress=True):
# We copy this method because the 'DataSaver' is not a
# subclass of the 'BaseLearner'.
BaseLearner.save(self, fname, compress)
@copy_docstring_from(BaseLearner.load)
def load(self, fname=None, compress=True):
# We copy this method because the 'DataSaver' is not a
# subclass of the 'BaseLearner'.
BaseLearner.load(self, fname, compress)
def _ds(learner_type, arg_picker, *args, **kwargs):
args = args[2:] # functools.partial passes the first 2 arguments in 'args'!
......@@ -46,12 +71,12 @@ def _ds(learner_type, arg_picker, *args, **kwargs):
def make_datasaver(learner_type, arg_picker):
"""Create a DataSaver of a `learner_type` that can be instantiated
"""Create a `DataSaver` of a `learner_type` that can be instantiated
with the `learner_type`'s key-word arguments.
Parameters
----------
learner_type : BaseLearner
learner_type : `~adaptive.BaseLearner` type
The learner type that needs to be wrapped.
arg_picker : function
Function that returns the argument that needs to be learned.
......@@ -59,15 +84,16 @@ def make_datasaver(learner_type, arg_picker):
Example
-------
Imagine we have a function that returns a dictionary
of the form: `{'y': y, 'err_est': err_est}`.
of the form: ``{'y': y, 'err_est': err_est}``.
>>> DataSaver = make_datasaver(Learner1D,
... arg_picker=operator.itemgetter('y'))
>>> from operator import itemgetter
>>> DataSaver = make_datasaver(Learner1D, arg_picker=itemgetter('y'))
>>> learner = DataSaver(function=f, bounds=(-1.0, 1.0))
Or when using `BalacingLearner.from_product`:
Or when using `adaptive.BalancingLearner.from_product`:
>>> learner_type = make_datasaver(adaptive.Learner1D,
... arg_picker=operator.itemgetter('y'))
... arg_picker=itemgetter('y'))
>>> learner = adaptive.BalancingLearner.from_product(
... jacobi, learner_type, dict(bounds=(0, 1)), combos)
"""
......
......@@ -330,7 +330,7 @@ class IntegratorLearner(BaseLearner):
The integral value in `self.bounds`.
err : float
The absolute error associated with `self.igral`.
max_ivals : int, default 1000
max_ivals : int, default: 1000
Maximum number of intervals that can be present in the calculation
of the integral. If this amount exceeds max_ivals, the interval
with the smallest error will be discarded.
......@@ -487,6 +487,7 @@ class IntegratorLearner(BaseLearner):
@property
def npoints(self):
"""Number of evaluated points."""
return len(self.done_points)
@property
......@@ -525,3 +526,30 @@ class IntegratorLearner(BaseLearner):
xs, ys = zip(*[(x, y) for ival in ivals
for x, y in sorted(ival.done_points.items())])
return hv.Path((xs, ys))
def _get_data(self):
# Change the defaultdict of SortedSets to a normal dict of sets.
x_mapping = {k: set(v) for k, v in self.x_mapping.items()}
return (self.priority_split,
self.done_points,
self.pending_points,
self._stack,
x_mapping,
self.ivals,
self.first_ival)
def _set_data(self, data):
self.priority_split, self.done_points, self.pending_points, \
self._stack, x_mapping, self.ivals, self.first_ival = data
# Add the pending_points to the _stack such that they are evaluated again
for x in self.pending_points:
if x not in self._stack:
self._stack.append(x)
# x_mapping is a data structure that can't easily be saved
# so we recreate it here
self.x_mapping = defaultdict(lambda: SortedSet([], key=attrgetter('rdepth')))
for k, _set in x_mapping.items():
self.x_mapping[k].update(_set)
......@@ -3,16 +3,72 @@ from copy import deepcopy
import heapq
import itertools
import math
from collections import Iterable
import numpy as np
import sortedcontainers
from .base_learner import BaseLearner
from .learnerND import volume
from .triangulation import simplex_volume_in_embedding
from ..notebook_integration import ensure_holoviews
from ..utils import cache_latest
def uniform_loss(interval, scale, function_values):
def uses_nth_neighbors(n):
"""Decorator to specify how many neighboring intervals the loss function uses.
Wraps loss functions to indicate that they expect intervals together
with ``n`` nearest neighbors
The loss function will then receive the data of the N nearest neighbors
(``nth_neighbors``) aling with the data of the interval itself in a dict.
The `~adaptive.Learner1D` will also make sure that the loss is updated
whenever one of the ``nth_neighbors`` changes.
Examples
--------
The next function is a part of the `curvature_loss_function` function.
>>> @uses_nth_neighbors(1)
... def triangle_loss(xs, ys):
... xs = [x for x in xs if x is not None]
... ys = [y for y in ys if y is not None]
...
... if len(xs) == 2: # we do not have enough points for a triangle
... return xs[1] - xs[0]
...
... N = len(xs) - 2 # number of constructed triangles
... if isinstance(ys[0], Iterable):
... pts = [(x, *y) for x, y in zip(xs, ys)]
... vol = simplex_volume_in_embedding
... else:
... pts = [(x, y) for x, y in zip(xs, ys)]
... vol = volume
... return sum(vol(pts[i:i+3]) for i in range(N)) / N
Or you may define a loss that favours the (local) minima of a function,
assuming that you know your function will have a single float as output.
>>> @uses_nth_neighbors(1)
... def local_minima_resolving_loss(xs, ys):
... dx = xs[2] - xs[1] # the width of the interval of interest
...
... if not ((ys[0] is not None and ys[0] > ys[1])
... or (ys[3] is not None and ys[3] > ys[2])):
... return loss * 100
...
... return loss
"""
def _wrapped(loss_per_interval):
loss_per_interval.nth_neighbors = n
return loss_per_interval
return _wrapped
@uses_nth_neighbors(0)
def uniform_loss(xs, ys):
"""Loss function that samples the domain uniformly.
Works with `~adaptive.Learner1D` only.
......@@ -27,33 +83,58 @@ def uniform_loss(interval, scale, function_values):
... loss_per_interval=uniform_sampling_1d)
>>>
"""
x_left, x_right = interval
x_scale, _ = scale
dx = (x_right - x_left) / x_scale
dx = xs[1] - xs[0]
return dx
def default_loss(interval, scale, function_values):
"""Calculate loss on a single interval
@uses_nth_neighbors(0)
def default_loss(xs, ys):
"""Calculate loss on a single interval.
Currently returns the rescaled length of the interval. If one of the
y-values is missing, returns 0 (so the intervals with missing data are
never touched. This behavior should be improved later.
"""
x_left, x_right = interval
y_right, y_left = function_values[x_right], function_values[x_left]
x_scale, y_scale = scale
dx = (x_right - x_left) / x_scale
if y_scale == 0:
loss = dx
dx = xs[1] - xs[0]
if isinstance(ys[0], Iterable):
dy = [abs(a-b) for a, b in zip(*ys)]
return np.hypot(dx, dy).max()
else:
dy = (y_right - y_left) / y_scale
try:
len(dy)
loss = np.hypot(dx, dy).max()
except TypeError:
loss = math.hypot(dx, dy)
return loss
dy = ys[1] - ys[0]
return np.hypot(dx, dy)
@uses_nth_neighbors(1)
def triangle_loss(xs, ys):
xs = [x for x in xs if x is not None]
ys = [y for y in ys if y is not None]
if len(xs) == 2: # we do not have enough points for a triangle
return xs[1] - xs[0]
N = len(xs) - 2 # number of constructed triangles
if isinstance(ys[0], Iterable):
pts = [(x, *y) for x, y in zip(xs, ys)]
vol = simplex_volume_in_embedding
else:
pts = [(x, y) for x, y in zip(xs, ys)]
vol = volume
return sum(vol(pts[i:i+3]) for i in range(N)) / N
def curvature_loss_function(area_factor=1, euclid_factor=0.02, horizontal_factor=0.02):
@uses_nth_neighbors(1)
def curvature_loss(xs, ys):
xs_middle = xs[1:3]
ys_middle = xs[1:3]
triangle_loss_ = triangle_loss(xs, ys)
default_loss_ = default_loss(xs_middle, ys_middle)
dx = xs_middle[0] - xs_middle[0]
return (area_factor * (triangle_loss_**0.5)
+ euclid_factor * default_loss_
+ horizontal_factor * dx)
return curvature_loss
def linspace(x_left, x_right, n):
......@@ -79,6 +160,15 @@ def _get_neighbors_from_list(xs):
return sortedcontainers.SortedDict(neighbors)
def _get_intervals(x, neighbors, nth_neighbors):
nn = nth_neighbors
i = neighbors.index(x)
start = max(0, i - nn - 1)
end = min(len(neighbors), i + nn + 2)
points = neighbors.keys()[start:end]
return list(zip(points, points[1:]))
class Learner1D(BaseLearner):
"""Learns and predicts a function 'f:ℝ → ℝ^N'.
......@@ -94,23 +184,43 @@ class Learner1D(BaseLearner):
If not provided, then a default is used, which uses the scaled distance
in the x-y plane as the loss. See the notes for more details.
Attributes
----------
data : dict
Sampled points and values.
pending_points : set
Points that still have to be evaluated.
Notes
-----
'loss_per_interval' takes 3 parameters: interval, scale, and function_values,
and returns a scalar; the loss over the interval.
interval : (float, float)
The bounds of the interval.
scale : (float, float)
The x and y scale over all the intervals, useful for rescaling the
interval loss.
function_values : dict(float -> float)
A map containing evaluated function values. It is guaranteed
to have values for both of the points in 'interval'.
`loss_per_interval` takes 2 parameters: ``xs`` and ``ys``, and returns a
scalar; the loss over the interval.
xs : tuple of floats
The x values of the interval, if `nth_neighbors` is greater than zero it
also contains the x-values of the neighbors of the interval, in ascending
order. The interval we want to know the loss of is then the middle
interval. If no neighbor is available (at the edges of the domain) then
`None` will take the place of the x-value of the neighbor.
ys : tuple of function values
The output values of the function when evaluated at the `xs`. This is
either a float or a tuple of floats in the case of vector output.
The `loss_per_interval` function may also have an attribute `nth_neighbors`
that indicates how many of the neighboring intervals to `interval` are used.
If `loss_per_interval` doesn't have such an attribute, it's assumed that is
uses **no** neighboring intervals. Also see the `uses_nth_neighbors`
decorator for more information.
"""
def __init__(self, function, bounds, loss_per_interval=None):
self.function = function
if hasattr(loss_per_interval, 'nth_neighbors'):
self.nth_neighbors = loss_per_interval.nth_neighbors
else:
self.nth_neighbors = 0
self.loss_per_interval = loss_per_interval or default_loss
# A dict storing the loss function for each interval x_n.
......@@ -141,6 +251,12 @@ class Learner1D(BaseLearner):
@property
def vdim(self):
"""Length of the output of ``learner.function``.
If the output is unsized (when it's a scalar)
then `vdim = 1`.
As long as no data is known `vdim = 1`.
"""
if self._vdim is None:
if self.data:
y = next(iter(self.data.values()))
......@@ -155,6 +271,7 @@ class Learner1D(BaseLearner):
@property
def npoints(self):
"""Number of evaluated points."""
return len(self.data)
@cache_latest
......@@ -162,41 +279,77 @@ class Learner1D(BaseLearner):
losses = self.losses if real else self.losses_combined
return max(losses.values()) if len(losses) > 0 else float('inf')
def update_interpolated_loss_in_interval(self, x_left, x_right):
if x_left is not None and x_right is not None:
dx = x_right - x_left
if dx < self._dx_eps:
loss = 0
else:
loss = self.loss_per_interval((x_left, x_right),
self._scale, self.data)
self.losses[x_left, x_right] = loss
# Iterate over all interpolated intervals in between
# x_left and x_right and set the newly interpolated loss.
a, b = x_left, None
while b != x_right:
b = self.neighbors_combined[a][1]
self.losses_combined[a, b] = (b - a) * loss / dx
a = b
def update_losses(self, x, real=True):
def _scale_x(self, x):
if x is None:
return None
return x / self._scale[0]
def _scale_y(self, y):
if y is None:
return None
y_scale = self._scale[1] or 1
return y / y_scale
def _get_point_by_index(self, ind):
if ind < 0 or ind >= len(self.neighbors):
return None
return self.neighbors.keys()[ind]
def _get_loss_in_interval(self, x_left, x_right):
assert x_left is not None and x_right is not None
if x_right - x_left < self._dx_eps:
return 0
nn = self.nth_neighbors
i = self.neighbors.index(x_left)
start = i - nn
end = i + nn + 2
xs = [self._get_point_by_index(i) for i in range(start, end)]
ys = [self.data.get(x, None) for x in xs]
xs_scaled = tuple(self._scale_x(x) for x in xs)
ys_scaled = tuple(self._scale_y(y) for y in ys)
# we need to compute the loss for this interval
return self.loss_per_interval(xs_scaled, ys_scaled)
def _update_interpolated_loss_in_interval(self, x_left, x_right):
if x_left is None or x_right is None:
return
loss = self._get_loss_in_interval(x_left, x_right)
self.losses[x_left, x_right] = loss
# Iterate over all interpolated intervals in between
# x_left and x_right and set the newly interpolated loss.
a, b = x_left, None
dx = x_right - x_left
while b != x_right:
b = self.neighbors_combined[a][1]
self.losses_combined[a, b] = (b - a) * loss / dx
a = b
def _update_losses(self, x, real=True):
"""Update all losses that depend on x"""
# When we add a new point x, we should update the losses
# (x_left, x_right) are the "real" neighbors of 'x'.
x_left, x_right = self.find_neighbors(x, self.neighbors)
x_left, x_right = self._find_neighbors(x, self.neighbors)
# (a, b) are the neighbors of the combined interpolated
# and "real" intervals.
a, b = self.find_neighbors(x, self.neighbors_combined)
a, b = self._find_neighbors(x, self.neighbors_combined)
# (a, b) is splitted into (a, x) and (x, b) so if (a, b) exists
self.losses_combined.pop((a, b), None) # we get rid of (a, b).
if real:
# We need to update all interpolated losses in the interval
# (x_left, x) and (x, x_right). Since the addition of the point
# 'x' could change their loss.
self.update_interpolated_loss_in_interval(x_left, x)
self.update_interpolated_loss_in_interval(x, x_right)
# (x_left, x), (x, x_right) and the nth_neighbors nearest
# neighboring intervals. Since the addition of the
# point 'x' could change their loss.
for ival in _get_intervals(x, self.neighbors, self.nth_neighbors):
self._update_interpolated_loss_in_interval(*ival)
# Since 'x' is in between (x_left, x_right),
# we get rid of the interval.
......@@ -223,7 +376,7 @@ class Learner1D(BaseLearner):
self.losses_combined[x, b] = float('inf')
@staticmethod
def find_neighbors(x, neighbors):
def _find_neighbors(x, neighbors):
if x in neighbors:
return neighbors[x]
pos = neighbors.bisect_left(x)
......@@ -232,14 +385,14 @@ class Learner1D(BaseLearner):
x_right = keys[pos] if pos != len(neighbors) else None
return x_left, x_right
def update_neighbors(self, x, neighbors):
def _update_neighbors(self, x, neighbors):
if x not in neighbors: # The point is new
x_left, x_right = self.find_neighbors(x, neighbors)
x_left, x_right = self._find_neighbors(x, neighbors)
neighbors[x] = [x_left, x_right]
neighbors.get(x_left, [None, None])[1] = x
neighbors.get(x_right, [None, None])[0] = x
def update_scale(self, x, y):
def _update_scale(self, x, y):
"""Update the scale with which the x and y-values are scaled.
For a learner where the function returns a single scalar the scale
......@@ -270,6 +423,9 @@ class Learner1D(BaseLearner):
if x in self.data:
# The point is already evaluated before
return
if y is None:
raise TypeError("Y-value may not be None, use learner.tell_pending(x)"
"to indicate that this value is currently being calculated")
# either it is a float/int, if not, try casting to a np.array
if not isinstance(y, (float, int)):
......@@ -284,16 +440,16 @@ class Learner1D(BaseLearner):
if not self.bounds[0] <= x <= self.bounds[1]:
return
self.update_neighbors(x, self.neighbors_combined)
self.update_neighbors(x, self.neighbors)
self.update_scale(x, y)
self.update_losses(x, real=True)
self._update_neighbors(x, self.neighbors_combined)
self._update_neighbors(x, self.neighbors)
self._update_scale(x, y)
self._update_losses(x, real=True)
# If the scale has increased enough, recompute all losses.
if self._scale[1] > 2 * self._oldscale[1]:
for interval in self.losses:
self.update_interpolated_loss_in_interval(*interval)
self._update_interpolated_loss_in_interval(*interval)
self._oldscale = deepcopy(self._scale)
......@@ -302,8 +458,8 @@ class Learner1D(BaseLearner):
# The point is already evaluated before
return
self.pending_points.add(x)
self.update_neighbors(x, self.neighbors_combined)
self.update_losses(x, real=False)
self._update_neighbors(x, self.neighbors_combined)
self._update_losses(x, real=False)
def tell_many(self, xs, ys, *, force=False):
if not force and not (len(xs) > 0.5 * len(self.data) and len(xs) > 2):
......@@ -342,10 +498,8 @@ class Learner1D(BaseLearner):
# The the losses for the "real" intervals.
self.losses = {}
for x_left, x_right in intervals:
self.losses[x_left, x_right] = (
self.loss_per_interval((x_left, x_right), self._scale, self.data)
if x_right - x_left >= self._dx_eps else 0)
for ival in intervals:
self.losses[ival] = self._get_loss_in_interval(*ival)
# List with "real" intervals that have interpolated intervals inside
to_interpolate = []
......@@ -363,7 +517,7 @@ class Learner1D(BaseLearner):
x_left, x_right = ival
a, b = to_interpolate[-1] if to_interpolate else (None, None)
if b == x_left and (a, b) not in self.losses:
# join (a, b) and (x_left, x_right) --> (a, x_right)
# join (a, b) and (x_left, x_right) (a, x_right)
to_interpolate[-1] = (a, x_right)
else:
to_interpolate.append((x_left, x_right))
......@@ -372,10 +526,10 @@ class Learner1D(BaseLearner):
if ival in self.losses:
# If this interval does not exist it should already
# have an inf loss.
self.update_interpolated_loss_in_interval(*ival)
self._update_interpolated_loss_in_interval(*ival)
def ask(self, n, tell_pending=True):
"""Return n points that are expected to maximally reduce the loss."""
"""Return 'n' points that are expected to maximally reduce the loss."""
points, loss_improvements = self._ask_points_without_adding(n)
if tell_pending:
......@@ -385,7 +539,7 @@ class Learner1D(BaseLearner):
return points, loss_improvements
def _ask_points_without_adding(self, n):
"""Return n points that are expected to maximally reduce the loss.
"""Return 'n' points that are expected to maximally reduce the loss.
Without altering the state of the learner"""
# Find out how to divide the n points over the intervals
# by finding positive integer n_i that minimize max(L_i / n_i) subject
......@@ -459,6 +613,14 @@ class Learner1D(BaseLearner):
return points, loss_improvements
def plot(self):
"""Returns a plot of the evaluated data.
Returns
-------
plot : `holoviews.element.Scatter` (if vdim=1)\
else `holoviews.element.Path`
Plot of the evaluated data.
"""
hv = ensure_holoviews()
if not self.data:
p = hv.Scatter([]) * hv.Path([])
......@@ -478,3 +640,9 @@ class Learner1D(BaseLearner):
self.pending_points = set()
self.losses_combined = deepcopy(self.losses)
self.neighbors_combined = deepcopy(self.neighbors)
def _get_data(self):
return self.data
def _set_data(self, data):
self.tell_many(*zip(*data.items()))
# -*- coding: utf-8 -*-
from collections import OrderedDict
from copy import copy
import itertools
from math import sqrt
......@@ -14,6 +15,19 @@ from ..utils import cache_latest
# Learner2D and helper functions.
def deviations(ip):
"""Returns the deviation of the linear estimate.
Is useful when defining custom loss functions.
Parameters
----------
ip : `scipy.interpolate.LinearNDInterpolator` instance
Returns
-------
numpy array
The deviation per triangle.
"""
values = ip.values / (ip.values.ptp(axis=0).max() or 1)
gradients = interpolate.interpnd.estimate_gradients_2d_global(
ip.tri, values, tol=1e-6)
......@@ -36,6 +50,20 @@ def deviations(ip):
def areas(ip):
"""Returns the area per triangle of the triangulation inside
a `LinearNDInterpolator` instance.
Is useful when defining custom loss functions.
Parameters
----------
ip : `scipy.interpolate.LinearNDInterpolator` instance
Returns
-------
numpy array
The area per triangle in ``ip.tri``.
"""
p = ip.tri.points[ip.tri.vertices]
q = p[:, :-1, :] - p[:, -1, None, :]
areas = abs(q[:, 0, 0] * q[:, 1, 1] - q[:, 0, 1] * q[:, 1, 0]) / 2
......@@ -61,9 +89,9 @@ def uniform_loss(ip):
return np.sqrt(areas(ip))
def resolution_loss(ip, min_distance=0, max_distance=1):
"""Loss function that is similar to the default loss function, but you can
set the maximimum and minimum size of a triangle.
def resolution_loss_function(min_distance=0, max_distance=1):
"""Loss function that is similar to the `default_loss` function, but you
can set the maximimum and minimum size of a triangle.
Works with `~adaptive.Learner2D` only.
......@@ -76,33 +104,31 @@ def resolution_loss(ip, min_distance=0, max_distance=1):
... x, y = xy
... return x**2 + y**2
>>>
>>> from functools import partial
>>> loss = partial(resolution_loss, min_distance=0.01)
>>> loss = resolution_loss_function(min_distance=0.01, max_distance=1)
>>> learner = adaptive.Learner2D(f,
... bounds=[(-1, -1), (1, 1)],
... loss_per_triangle=loss)
>>>
"""
A = areas(ip)
dev = np.sum(deviations(ip), axis=0)
def resolution_loss(ip):
loss = default_loss(ip)
# similar to the default_loss
loss = np.sqrt(A) * dev + A
A = areas(ip)
# Setting areas with a small area to zero such that they won't be chosen again
loss[A < min_distance**2] = 0
# Setting areas with a small area to zero such that they won't be chosen again
loss[A < min_distance**2] = 0
# Setting triangles that have a size larger than max_distance to infinite loss
# such that these triangles will be picked
loss[A > max_distance**2] = np.inf
# Setting triangles that have a size larger than max_distance to infinite loss
# such that these triangles will be picked
loss[A > max_distance**2] = np.inf
return loss
return loss
return resolution_loss
def minimize_triangle_surface_loss(ip):
"""Loss function that is similar to the default loss function in the
`Learner1D`. The loss is the area spanned by the 3D vectors of the
vertices.
`~adaptive.Learner1D`. The loss is the area spanned by the 3D
vectors of the vertices.
Works with `~adaptive.Learner2D` only.
......@@ -206,15 +232,15 @@ class Learner2D(BaseLearner):
pending_points : set
Points that still have to be evaluated and are currently
interpolated, see `data_combined`.
stack_size : int, default 10
stack_size : int, default: 10
The size of the new candidate points stack. Set it to 1
to recalculate the best points at each call to `ask`.
aspect_ratio : float, int, default 1
Average ratio of `x` span over `y` span of a triangle. If
there is more detail in either `x` or `y` the `aspect_ratio`
needs to be adjusted. When `aspect_ratio > 1` the
triangles will be stretched along `x`, otherwise
along `y`.
aspect_ratio : float, int, default: 1
Average ratio of ``x`` span over ``y`` span of a triangle. If
there is more detail in either ``x`` or ``y`` the ``aspect_ratio``
needs to be adjusted. When ``aspect_ratio > 1`` the
triangles will be stretched along ``x``, otherwise
along ``y``.
Methods
-------
......@@ -239,13 +265,13 @@ class Learner2D(BaseLearner):
This sampling procedure is not extremely fast, so to benefit from
it, your function needs to be slow enough to compute.
'loss_per_triangle' takes a single parameter, 'ip', which is a
`loss_per_triangle` takes a single parameter, `ip`, which is a
`scipy.interpolate.LinearNDInterpolator`. You can use the
*undocumented* attributes 'tri' and 'values' of 'ip' to get a
*undocumented* attributes ``tri`` and ``values`` of `ip` to get a
`scipy.spatial.Delaunay` and a vector of function values.
These can be used to compute the loss. The functions
`adaptive.learner.learner2D.areas` and
`adaptive.learner.learner2D.deviations` to calculate the
`~adaptive.learner.learner2D.areas` and
`~adaptive.learner.learner2D.deviations` to calculate the
areas and deviations from a linear interpolation
over each triangle.
"""
......@@ -288,10 +314,17 @@ class Learner2D(BaseLearner):
@property
def npoints(self):
"""Number of evaluated points."""
return len(self.data)
@property
def vdim(self):
"""Length of the output of ``learner.function``.
If the output is unsized (when it's a scalar)
then `vdim = 1`.
As long as no data is known `vdim = 1`.
"""
if self._vdim is None and self.data:
try:
value = next(iter(self.data.values()))
......@@ -336,11 +369,15 @@ class Learner2D(BaseLearner):
return points_combined, values_combined
def data_combined(self):
"""Like `data`, however this includes the points in
`pending_points` for which the values are interpolated."""
# Interpolate the unfinished points
points, values = self._data_combined()
return {tuple(k): v for k, v in zip(points, values)}
def ip(self):
"""A `scipy.interpolate.LinearNDInterpolator` instance
containing the learner's data."""
if self._ip is None:
points, values = self._data_in_bounds()
points = self._scale(points)
......@@ -348,6 +385,9 @@ class Learner2D(BaseLearner):
return self._ip
def ip_combined(self):
"""A `scipy.interpolate.LinearNDInterpolator` instance
containing the learner's data *and* interpolated data of
the `pending_points`."""
if self._ip_combined is None:
points, values = self._data_combined()
points = self._scale(points)
......@@ -464,19 +504,21 @@ class Learner2D(BaseLearner):
Number of points in x and y. If None (default) this number is
evaluated by looking at the size of the smallest triangle.
tri_alpha : float
The opacity (0 <= tri_alpha <= 1) of the triangles overlayed on
top of the image. By default the triangulation is not visible.
The opacity ``(0 <= tri_alpha <= 1)`` of the triangles overlayed
on top of the image. By default the triangulation is not visible.
Returns
-------
plot : holoviews.Overlay or holoviews.HoloMap
A `holoviews.Overlay` of `holoviews.Image * holoviews.EdgePaths`.
If the `learner.function` returns a vector output, a
`holoviews.HoloMap` of the `holoviews.Overlay`s wil be returned.
plot : `holoviews.core.Overlay` or `holoviews.core.HoloMap`
A `holoviews.core.Overlay` of
``holoviews.Image * holoviews.EdgePaths``. If the
`learner.function` returns a vector output, a
`holoviews.core.HoloMap` of the
`holoviews.core.Overlay`\s wil be returned.
Notes
-----
The plot object that is returned if `learner.function` returns a
The plot object that is returned if ``learner.function`` returns a
vector *cannot* be used with the live_plotting functionality.
"""
hv = ensure_holoviews()
......@@ -520,3 +562,13 @@ class Learner2D(BaseLearner):
no_hover = dict(plot=dict(inspection_policy=None, tools=[]))
return im.opts(style=im_opts) * tris.opts(style=tri_opts, **no_hover)
def _get_data(self):
return self.data
def _set_data(self, data):
self.data = data
# Remove points from stack if they already exist
for point in copy(self._stack):
if point in self.data:
self._stack.pop(point)
# -*- coding: utf-8 -*-
from collections import OrderedDict, Iterable
import functools
import heapq
import itertools
import random
......@@ -10,20 +11,21 @@ import scipy.spatial
from .base_learner import BaseLearner
from ..notebook_integration import ensure_holoviews
from ..notebook_integration import ensure_holoviews, ensure_plotly
from .triangulation import (Triangulation, point_in_simplex,
circumsphere, simplex_volume_in_embedding)
circumsphere, simplex_volume_in_embedding,
fast_det)
from ..utils import restore, cache_latest
def volume(simplex, ys=None):
# Notice the parameter ys is there so you can use this volume method as
# as loss function
matrix = np.array(np.subtract(simplex[:-1], simplex[-1]), dtype=float)
dim = len(simplex) - 1
matrix = np.subtract(simplex[:-1], simplex[-1], dtype=float)
# See https://www.jstor.org/stable/2315353
vol = np.abs(np.linalg.det(matrix)) / np.math.factorial(dim)
dim = len(simplex) - 1
vol = np.abs(fast_det(matrix)) / np.math.factorial(dim)
return vol
......@@ -106,9 +108,10 @@ class LearnerND(BaseLearner):
func: callable
The function to learn. Must take a tuple of N real
parameters and return a real number or an arraylike of length M.
bounds : list of 2-tuples
bounds : list of 2-tuples or `scipy.spatial.ConvexHull`
A list ``[(a_1, b_1), (a_2, b_2), ..., (a_n, b_n)]`` containing bounds,
one pair per dimension.
Or a ConvexHull that defines the boundary of the domain.
loss_per_simplex : callable, optional
A function that returns the loss for a simplex.
If not provided, then a default is used, which uses
......@@ -124,15 +127,8 @@ class LearnerND(BaseLearner):
Coordinates of the currently known points
values : numpy array
The values of each of the known points
Methods
-------
plot(n)
If dim == 2, this method will plot the function being learned.
plot_slice(cut_mapping, n)
plot a slice of the function using interpolation of the current data.
the cut_mapping contains the fixed parameters, the other parameters are
used as axes for plotting.
pending_points : set
Points that still have to be evaluated.
Notes
-----
......@@ -156,26 +152,33 @@ class LearnerND(BaseLearner):
"""
def __init__(self, func, bounds, loss_per_simplex=None):
self.ndim = len(bounds)
self._vdim = None
self.loss_per_simplex = loss_per_simplex or default_loss
self.bounds = tuple(tuple(map(float, b)) for b in bounds)
self.data = OrderedDict()
self.pending_points = set()
self._bounds_points = list(map(tuple, itertools.product(*bounds)))
if isinstance(bounds, scipy.spatial.ConvexHull):
hull_points = bounds.points[bounds.vertices]
self._bounds_points = sorted(list(map(tuple, hull_points)))
self._bbox = tuple(zip(hull_points.min(axis=0), hull_points.max(axis=0)))
self._interior = scipy.spatial.Delaunay(self._bounds_points)
else:
self._bounds_points = sorted(list(map(tuple, itertools.product(*bounds))))
self._bbox = tuple(tuple(map(float, b)) for b in bounds)
self.ndim = len(self._bbox)
self.function = func
self._tri = None
self._losses = dict()
self._pending_to_simplex = dict() # vertex -> simplex
self._pending_to_simplex = dict() # vertex simplex
# triangulation of the pending points inside a specific simplex
self._subtriangulations = dict() # simplex -> triangulation
self._subtriangulations = dict() # simplex triangulation
# scale to unit
self._transform = np.linalg.inv(np.diag(np.diff(bounds).flat))
self._transform = np.linalg.inv(np.diag(np.diff(self._bbox).flat))
# create a private random number generator with fixed seed
self._random = random.Random(1)
......@@ -195,10 +198,17 @@ class LearnerND(BaseLearner):
@property
def npoints(self):
"""Number of evaluated points."""
return len(self.data)
@property
def vdim(self):
"""Length of the output of ``learner.function``.
If the output is unsized (when it's a scalar)
then `vdim = 1`.
As long as no data is known `vdim = 1`.
"""
if self._vdim is None and self.data:
try:
value = next(iter(self.data.values()))
......@@ -212,11 +222,15 @@ class LearnerND(BaseLearner):
return all(p in self.data for p in self._bounds_points)
def ip(self):
"""A `scipy.interpolate.LinearNDInterpolator` instance
containing the learner's data."""
# XXX: take our own triangulation into account when generating the ip
return interpolate.LinearNDInterpolator(self.points, self.values)
@property
def tri(self):
"""An `adaptive.learner.triangulation.Triangulation` instance
with all the points of the learner."""
if self._tri is not None:
return self._tri
......@@ -226,16 +240,18 @@ class LearnerND(BaseLearner):
return self._tri
except ValueError:
# A ValueError is raised if we do not have enough points or
# the provided points are coplanar, so we need more points to create
# a valid triangulation
# the provided points are coplanar, so we need more points to
# create a valid triangulation
return None
@property
def values(self):
"""Get the values from `data` as a numpy array."""
return np.array(list(self.data.values()), dtype=float)
@property
def points(self):
"""Get the points from `data` as a numpy array."""
return np.array(list(self.data.keys()), dtype=float)
def tell(self, point, value):
......@@ -267,7 +283,13 @@ class LearnerND(BaseLearner):
return simplex in self.tri.simplices
def inside_bounds(self, point):
return all(mn <= p <= mx for p, (mn, mx) in zip(point, self.bounds))
"""Check whether a point is inside the bounds."""
if hasattr(self, '_interior'):
return self._interior.find_simplex(point, tol=1e-8) >= 0
else:
eps = 1e-8
return all((mn - eps) <= p <= (mx + eps) for p, (mn, mx)
in zip(point, self._bbox))
def tell_pending(self, point, *, simplex=None):
point = tuple(point)
......@@ -341,11 +363,13 @@ class LearnerND(BaseLearner):
assert not self._bounds_available
# pick a random point inside the bounds
# XXX: change this into picking a point based on volume loss
a = np.diff(self.bounds).flat
b = np.array(self.bounds)[:, 0]
r = np.array([self._random.random() for _ in range(self.ndim)])
p = r * a + b
p = tuple(p)
a = np.diff(self._bbox).flat
b = np.array(self._bbox)[:, 0]
p = None
while p is None or not self.inside_bounds(p):
r = np.array([self._random.random() for _ in range(self.ndim)])
p = r * a + b
p = tuple(p)
self.tell_pending(p)
return p, np.inf
......@@ -463,6 +487,10 @@ class LearnerND(BaseLearner):
self._subtriangulations = dict()
self._pending_to_simplex = dict()
##########################
# Plotting related stuff #
##########################
def plot(self, n=None, tri_alpha=0):
"""Plot the function we want to learn, only works in 2D.
......@@ -477,10 +505,10 @@ class LearnerND(BaseLearner):
if self.vdim > 1:
raise NotImplementedError('holoviews currently does not support',
'3D surface plots in bokeh.')
if len(self.bounds) != 2:
if len(self.ndim) != 2:
raise NotImplementedError("Only 2D plots are implemented: You can "
"plot a 2D slice with 'plot_slice'.")
x, y = self.bounds
x, y = self._bbox
lbrt = x[0], y[0], x[1], y[1]
if len(self.data) >= 4:
......@@ -517,13 +545,14 @@ class LearnerND(BaseLearner):
return im.opts(style=im_opts) * tris.opts(style=tri_opts, **no_hover)
def plot_slice(self, cut_mapping, n=None):
"""Plot a 1d or 2d interpolated slice of a N-dimensional function.
"""Plot a 1D or 2D interpolated slice of a N-dimensional function.
Parameters
----------
cut_mapping : dict (int -> float)
cut_mapping : dict (int float)
for each fixed dimension the value, the other dimensions
are interpolated
are interpolated. e.g. ``cut_mapping = {0: 1}``, so from
dimension 0 ('x') to value 1.
n : int
the number of boxes in the interpolation grid along each axis
"""
......@@ -536,7 +565,7 @@ class LearnerND(BaseLearner):
raise NotImplementedError('multidimensional output not yet'
' supported by `plot_slice`')
n = n or 201
values = [cut_mapping.get(i, np.linspace(*self.bounds[i], n))
values = [cut_mapping.get(i, np.linspace(*self._bbox[i], n))
for i in range(self.ndim)]
ind = next(i for i in range(self.ndim) if i not in cut_mapping)
x = values[ind]
......@@ -561,9 +590,9 @@ class LearnerND(BaseLearner):
xys = [xs[:, None], ys[None, :]]
values = [cut_mapping[i] if i in cut_mapping
else xys.pop(0) * (b[1] - b[0]) + b[0]
for i, b in enumerate(self.bounds)]
for i, b in enumerate(self._bbox)]
lbrt = [b for i, b in enumerate(self.bounds)
lbrt = [b for i, b in enumerate(self._bbox)
if i not in cut_mapping]
lbrt = np.reshape(lbrt, (2, 2)).T.flatten().tolist()
......@@ -576,3 +605,254 @@ class LearnerND(BaseLearner):
return im.opts(style=dict(cmap='viridis'))
else:
raise ValueError("Only 1 or 2-dimensional plots can be generated.")
def plot_3D(self, with_triangulation=False):
"""Plot the learner's data in 3D using plotly.
Does *not* work with the
`adaptive.notebook_integration.live_plot` functionality.
Parameters
----------
with_triangulation : bool, default: False
Add the verticices to the plot.
Returns
-------
plot : `plotly.offline.iplot` object
The 3D plot of ``learner.data``.
"""
plotly = ensure_plotly()
plots = []
vertices = self.tri.vertices
if with_triangulation:
Xe, Ye, Ze = [], [], []
for simplex in self.tri.simplices:
for s in itertools.combinations(simplex, 2):
Xe += [vertices[i][0] for i in s] + [None]
Ye += [vertices[i][1] for i in s] + [None]
Ze += [vertices[i][2] for i in s] + [None]
plots.append(plotly.graph_objs.Scatter3d(
x=Xe, y=Ye, z=Ze, mode='lines',
line=dict(color='rgb(125,125,125)', width=1),
hoverinfo='none'
))
Xn, Yn, Zn = zip(*vertices)
colors = [self.data[p] for p in self.tri.vertices]
marker = dict(symbol='circle', size=3, color=colors,
colorscale='Viridis',
line=dict(color='rgb(50,50,50)', width=0.5))
plots.append(plotly.graph_objs.Scatter3d(
x=Xn, y=Yn, z=Zn, mode='markers',
name='actors', marker=marker,
hoverinfo='text'
))
axis = dict(
showbackground=False,
showline=False,
zeroline=False,
showgrid=False,
showticklabels=False,
title='',
)
layout = plotly.graph_objs.Layout(
showlegend=False,
scene=dict(xaxis=axis, yaxis=axis, zaxis=axis),
margin=dict(t=100),
hovermode='closest')
fig = plotly.graph_objs.Figure(data=plots, layout=layout)
return plotly.offline.iplot(fig)
def _get_data(self):
return self.data
def _set_data(self, data):
self.tell_many(*zip(*data.items()))
def _get_iso(self, level=0.0, which='surface'):
if which == 'surface':
if self.ndim != 3 or self.vdim != 1:
raise Exception('Isosurface plotting is only supported'
' for a 3D input and 1D output')
get_surface = True
get_line = False
elif which == 'line':
if self.ndim != 2 or self.vdim != 1:
raise Exception('Isoline plotting is only supported'
' for a 2D input and 1D output')
get_surface = False
get_line = True
vertices = [] # index -> (x,y,z)
faces_or_lines = [] # tuple of indices of the corner points
@functools.lru_cache()
def _get_vertex_index(a, b):
vertex_a = self.tri.vertices[a]
vertex_b = self.tri.vertices[b]
value_a = self.data[vertex_a]
value_b = self.data[vertex_b]
da = abs(value_a - level)
db = abs(value_b - level)
dab = da + db
new_pt = (db / dab * np.array(vertex_a)
+ da / dab * np.array(vertex_b))
new_index = len(vertices)
vertices.append(new_pt)
return new_index
for simplex in self.tri.simplices:
plane_or_line = []
for a, b in itertools.combinations(simplex, 2):
va = self.data[self.tri.vertices[a]]
vb = self.data[self.tri.vertices[b]]
if min(va, vb) < level <= max(va, vb):
vi = _get_vertex_index(a, b)
should_add = True
for pi in plane_or_line:
if np.allclose(vertices[vi], vertices[pi]):
should_add = False
if should_add:
plane_or_line.append(vi)
if get_surface and len(plane_or_line) == 3:
faces_or_lines.append(plane_or_line)
elif get_surface and len(plane_or_line) == 4:
faces_or_lines.append(plane_or_line[:3])
faces_or_lines.append(plane_or_line[1:])
elif get_line and len(plane_or_line) == 2:
faces_or_lines.append(plane_or_line)
if len(faces_or_lines) == 0:
r_min = min(self.data[v] for v in self.tri.vertices)
r_max = max(self.data[v] for v in self.tri.vertices)
raise ValueError(
f"Could not draw isosurface for level={level}, as"
" this value is not inside the function range. Please choose"
f" a level strictly inside interval ({r_min}, {r_max})"
)
return vertices, faces_or_lines
def plot_isoline(self, level=0.0, n=None, tri_alpha=0):
"""Plot the isoline at a specific level, only works in 2D.
Parameters
----------
level : float, default: 0
The value of the function at which you would like to see
the isoline.
n : int
The number of boxes in the interpolation grid along each axis.
This is passed to `plot`.
tri_alpha : float
The opacity of the overlaying triangulation. This is passed
to `plot`.
Returns
-------
`holoviews.core.Overlay`
The plot of the isoline(s). This overlays a `plot` with a
`holoviews.element.Path`.
"""
hv = ensure_holoviews()
if n == -1:
plot = hv.Path([])
else:
plot = self.plot(n=n, tri_alpha=tri_alpha)
if isinstance(level, Iterable):
for l in level:
plot = plot * self.plot_isoline(level=l, n=-1)
return plot
vertices, lines = self.self._get_iso(level, which='line')
paths = [[vertices[i], vertices[j]] for i, j in lines]
contour = hv.Path(paths)
contour_opts = dict(color='black')
contour = contour.opts(style=contour_opts)
return plot * contour
def plot_isosurface(self, level=0.0, hull_opacity=0.2):
"""Plots a linearly interpolated isosurface.
This is the 3D analog of an isoline. Does *not* work with the
`adaptive.notebook_integration.live_plot` functionality.
Parameters
----------
level : float, default: 0.0
the function value which you are interested in.
hull_opacity : float, default: 0.0
the opacity of the hull of the domain.
Returns
-------
plot : `plotly.offline.iplot` object
The plot object of the isosurface.
"""
plotly = ensure_plotly()
vertices, faces = self._get_iso(level, which='surface')
x, y, z = zip(*vertices)
fig = plotly.figure_factory.create_trisurf(
x=x, y=y, z=z, plot_edges=False,
simplices=faces, title="Isosurface")
isosurface = fig.data[0]
isosurface.update(lighting=dict(ambient=1, diffuse=1,
roughness=1, specular=0, fresnel=0))
if hull_opacity < 1e-3:
# Do not compute the hull_mesh.
return plotly.offline.iplot(fig)
hull_mesh = self._get_hull_mesh(opacity=hull_opacity)
return plotly.offline.iplot([isosurface, hull_mesh])
def _get_hull_mesh(self, opacity=0.2):
plotly = ensure_plotly()
hull = scipy.spatial.ConvexHull(self._bounds_points)
# Find the colors of each plane, giving triangles which are coplanar
# the same color, such that a square face has the same color.
color_dict = {}
def _get_plane_color(simplex):
simplex = tuple(simplex)
# If the volume of the two triangles combined is zero then they
# belong to the same plane.
for simplex_key, color in color_dict.items():
points = [hull.points[i] for i in set(simplex_key + simplex)]
points = np.array(points)
if np.linalg.matrix_rank(points[1:] - points[0]) < 3:
return color
if scipy.spatial.ConvexHull(points).volume < 1e-5:
return color
color_dict[simplex] = tuple(random.randint(0, 255)
for _ in range(3))
return color_dict[simplex]
colors = [_get_plane_color(simplex) for simplex in hull.simplices]
x, y, z = zip(*self._bounds_points)
i, j, k = hull.simplices.T
lighting = dict(ambient=1, diffuse=1, roughness=1,
specular=0, fresnel=0)
return plotly.graph_objs.Mesh3d(x=x, y=y, z=z, i=i, j=j, k=k,
facecolor=colors, opacity=opacity,
lighting=lighting)
......@@ -8,18 +8,18 @@ from ..utils import cache_latest
class SKOptLearner(Optimizer, BaseLearner):
"""Learn a function minimum using 'skopt.Optimizer'.
"""Learn a function minimum using ``skopt.Optimizer``.
This is an 'Optimizer' from 'scikit-optimize',
This is an ``Optimizer`` from ``scikit-optimize``,
with the necessary methods added to make it conform
to the 'adaptive' learner interface.
to the ``adaptive`` learner interface.
Parameters
----------
function : callable
The function to learn.
**kwargs :
Arguments to pass to 'skopt.Optimizer'.
Arguments to pass to ``skopt.Optimizer``.
"""
def __init__(self, function, **kwargs):
......@@ -63,6 +63,7 @@ class SKOptLearner(Optimizer, BaseLearner):
@property
def npoints(self):
"""Number of evaluated points."""
return len(self.Xi)
def plot(self, nsamples=200):
......
......@@ -126,6 +126,17 @@ def fast_3d_circumcircle(points):
return center, radius
def fast_det(matrix):
matrix = np.asarray(matrix, dtype=float)
if matrix.shape == (2, 2):
return matrix[0][0] * matrix[1][1] - matrix[1][0] * matrix[0][1]
elif matrix.shape == (3, 3):
a, b, c, d, e, f, g, h, i = matrix.ravel()
return a * (e*i - f*h) - b * (d*i - f*g) + c * (d*h - e*g)
else:
return np.linalg.det(matrix)
def circumsphere(pts):
dim = len(pts) - 1
if dim == 2:
......@@ -140,9 +151,9 @@ def circumsphere(pts):
for i in range(1, len(pts)):
r = np.delete(mat, i, 1)
factor = (-1) ** (i + 1)
center.append(factor * np.linalg.det(r))
center.append(factor * fast_det(r))
a = np.linalg.det(np.delete(mat, 0, 1))
a = fast_det(np.delete(mat, 0, 1))
center = [x / (2 * a) for x in center]
x0 = pts[0]
......@@ -227,9 +238,11 @@ def simplex_volume_in_embedding(vertices) -> float:
sq_dists_mat = scipy.spatial.distance.squareform(bordered)
coeff = - (-2) ** (num_verts-1) * factorial(num_verts-1) ** 2
vol_square = np.linalg.det(sq_dists_mat) / coeff
vol_square = fast_det(sq_dists_mat) / coeff
if vol_square <= 0:
if vol_square < 0:
if abs(vol_square) < 1e-15:
return 0
raise ValueError('Provided vertices do not form a simplex')
return np.sqrt(vol_square)
......@@ -585,7 +598,7 @@ class Triangulation:
prefactor = np.math.factorial(self.dim)
vertices = np.array(self.get_vertices(simplex))
vectors = vertices[1:] - vertices[0]
return float(abs(np.linalg.det(vectors)) / prefactor)
return float(abs(fast_det(vectors)) / prefactor)
def volumes(self):
return [self.volume(sim) for sim in self.simplices]
......@@ -612,7 +625,7 @@ class Triangulation:
Parameters
----------
check : bool, default True
check : bool, default: True
Whether to raise an error if the computed hull is different from
stored.
......
......@@ -3,32 +3,48 @@ import importlib
import asyncio
from contextlib import suppress
import datetime
import random
import warnings
_async_enabled = False
_plotting_enabled = False
_holoviews_enabled = False
_ipywidgets_enabled = False
_plotly_enabled = False
def notebook_extension():
"""Enable ipywidgets, holoviews, and asyncio notebook integration."""
if not in_ipynb():
raise RuntimeError('"adaptive.notebook_extension()" may only be run '
'from a Jupyter notebook.')
global _plotting_enabled
_plotting_enabled = False
global _async_enabled, _holoviews_enabled, _ipywidgets_enabled
# Load holoviews
try:
_holoviews_enabled = False # After closing a notebook the js is gone
if not _holoviews_enabled:
import holoviews
holoviews.notebook_extension('bokeh', logo=False)
_holoviews_enabled = True
except ModuleNotFoundError:
warnings.warn("holoviews is not installed; plotting "
"is disabled.", RuntimeWarning)
# Load ipywidgets
try:
import ipywidgets
import holoviews
holoviews.notebook_extension('bokeh')
_plotting_enabled = True
if not _ipywidgets_enabled:
import ipywidgets
_ipywidgets_enabled = True
except ModuleNotFoundError:
warnings.warn("holoviews and (or) ipywidgets are not installed; plotting "
warnings.warn("ipywidgets is not installed; live_info "
"is disabled.", RuntimeWarning)
global _async_enabled
get_ipython().magic('gui asyncio')
_async_enabled = True
# Enable asyncio integration
if not _async_enabled:
get_ipython().magic('gui asyncio')
_async_enabled = True
def ensure_holoviews():
......@@ -38,6 +54,22 @@ def ensure_holoviews():
raise RuntimeError('holoviews is not installed; plotting is disabled.')
def ensure_plotly():
global _plotly_enabled
try:
import plotly
if not _plotly_enabled:
import plotly.graph_objs
import plotly.figure_factory
import plotly.offline
# This injects javascript and should happen only once
plotly.offline.init_notebook_mode()
_plotly_enabled = True
return plotly
except ModuleNotFoundError:
raise RuntimeError('plotly is not installed; plotting is disabled.')
def in_ipynb():
try:
# If we are running in IPython, then `get_ipython()` is always a global
......@@ -56,23 +88,23 @@ def live_plot(runner, *, plotter=None, update_interval=2, name=None):
Parameters
----------
runner : Runner
runner : `Runner`
plotter : function
A function that takes the learner as a argument and returns a
holoviews object. By default learner.plot() will be called.
holoviews object. By default ``learner.plot()`` will be called.
update_interval : int
Number of second between the updates of the plot.
name : hasable
Name for the `live_plot` task in `adaptive.active_plotting_tasks`.
By default the name is `None` and if another task with the same name
already exists that other live_plot is canceled.
By default the name is None and if another task with the same name
already exists that other `live_plot` is canceled.
Returns
-------
dm : holoviews.DynamicMap
The plot that automatically updates every update_interval.
dm : `holoviews.core.DynamicMap`
The plot that automatically updates every `update_interval`.
"""
if not _plotting_enabled:
if not _holoviews_enabled:
raise RuntimeError("Live plotting is not enabled; did you run "
"'adaptive.notebook_extension()'?")
......@@ -120,13 +152,31 @@ def live_plot(runner, *, plotter=None, update_interval=2, name=None):
return dm
def should_update(status):
try:
# Get the length of the write buffer size
buffer_size = len(status.comm.kernel.iopub_thread._events)
# Make sure to only keep all the messages when the notebook
# is viewed, this means 'buffer_size == 1'. However, when not
# viewing the notebook the buffer fills up. When this happens
# we decide to only add messages to it when a certain probability.
# i.e. we're offline for 12h, with an update_interval of 0.5s,
# and without the reduced probability, we have buffer_size=86400.
# With the correction this is np.log(86400) / np.log(1.1) = 119.2
return 1.1**buffer_size * random.random() < 1
except Exception:
# We catch any Exception because we are using a private API.
return True
def live_info(runner, *, update_interval=0.5):
"""Display live information about the runner.
Returns an interactive ipywidget that can be
visualized in a Jupyter notebook.
"""
if not _plotting_enabled:
if not _holoviews_enabled:
raise RuntimeError("Live plotting is not enabled; did you run "
"'adaptive.notebook_extension()'?")
......@@ -142,7 +192,12 @@ def live_info(runner, *, update_interval=0.5):
async def update():
while not runner.task.done():
await asyncio.sleep(update_interval)
status.value = _info_html(runner)
if should_update(status):
status.value = _info_html(runner)
else:
await asyncio.sleep(0.05)
status.value = _info_html(runner)
cancel.layout.display = 'none'
......@@ -176,7 +231,7 @@ def _info_html(runner):
with suppress(Exception):
info.append(('latest loss', f'{runner.learner._cache["loss"]:.3f}'))
template = '<dt>{}</dt><dd>{}</dd>'
template = '<dt class="ignore-css">{}</dt><dd>{}</dd>'
table = '\n'.join(template.format(k, v) for k, v in info)
return f'''
......
This diff is collapsed.