diff --git a/doc/source/pre/whatsnew/1.3.rst b/doc/source/pre/whatsnew/1.3.rst
index 707608cb65e2d54c7854f01b0ded189c2ce2d815..f7dbd8fc0ef4373499d976fc3b908f04548e629e 100644
--- a/doc/source/pre/whatsnew/1.3.rst
+++ b/doc/source/pre/whatsnew/1.3.rst
@@ -219,3 +219,10 @@ image in the fundamental domain.
 
 This change is documented here for completeness.  We expect that the vast
 majority of users of Kwant will not be affected by it.
+
+ .. _whatsnew13-params-api-change:
+
+API change that affects low-level systems
+-----------------------------------------
+The `~kwant.system.System.hamiltonian` method of low-level systems must now accept a
+`params` keyword parameter.
diff --git a/doc/source/pre/whatsnew/1.4.rst b/doc/source/pre/whatsnew/1.4.rst
index 0efda7e2d05d814838a4a26f2caf44a3246e4373..c962fb9741f8f6ac8709aa8f8bd62b8ced65f377 100644
--- a/doc/source/pre/whatsnew/1.4.rst
+++ b/doc/source/pre/whatsnew/1.4.rst
@@ -5,78 +5,63 @@ This article explains the user-visible changes in Kwant 1.4.0.  Subsequently,
 the user-visible changes for each maintenance release of the 1.4.x series are
 listed (if there were any).
 
-Value functions may no longer take unnamed arguments
-----------------------------------------------------
-The system parameter (``args`` and ``params``) handling of Kwant has been
-rewritten following the principle that each value function must take a clearly
-specified set of named parameters.  This allows to make the parameter handling
-less error-prone (as explained in the following section) and faster.  For this
-reason, value functions may no longer accept ``*args`` or ``**kwargs``.
-
-If you are using such functions, we recommend that you replace them by named
-parameters.
-
-Value functions may no longer have default values for parameters
-----------------------------------------------------------------
-Using value functions with default values for parameters can be
-problematic, especially when re-using value functions between simulations.
-When parameters have default values it is easy to forget that such a
-parameter exists at all, because it is not necessary to provide them explicitly
-to functions that use the Kwant system. This means that other value functions
-might be introduced that also depend on the same parameter,
-but in an inconsistent way (e.g. a parameter 'phi' that is a superconducting
-phase in one value function, but a peierls phase in another). This leads
-to bugs that are confusing and hard to track down.
-
-For this reason value functions may no longer have default values for paramters.
-Concretely this means that the following no longer works::
-
-  syst = kwant.Builder()
+Summary: release highlights
+---------------------------
+* Adding magnetic field to systems, even in complicated cases, is now specially
+  :ref:`supported <whatsnew14-magnetic>`.
+* The :ref:`KPM module can now calculate conductivities
+  <whatsnew14-kpm-conductivity>`.
+* The `Qsymm library <https://gitlab.kwant-project.org/qt/qsymm>`_ for
+  Hamiltonian symmetry analysis has been :ref:`integrated <whatsnew14-qsymm>`.
+* The handling of system parameters has been :ref:`improved
+  <whatsnew14-parameters>` and optimized.
+* Plotting has been improved, most notably through the addition of a :ref:`routine
+  that plots densities with interpolation <whatsnew14-density-plots>`.
+* :ref:`Installing Kwant on Windows <whatsnew14-windows>` is now much easier
+  thanks to Conda packages.
+
+Backwards-incompatible changes:
+
+* `Restrictions on value functions when named parameters are given`_
+
+.. _whatsnew14-magnetic:
 
-  # Parameter 't' has a default value of 1
-  def onsite(site, V, t=1):
-    return V = 2 * t
+Automatic Peierls phase calculation
+-----------------------------------
+When defining systems with orbital magnetic fields it is often cumbersome to
+manually calculate the phases required by the Peierls substitution, and to
+ensure that the chosen gauge is consistent across the whole system
+(this is especially true for systems with leads that point in different
+directions). This release introduces `kwant.physics.magnetic_gauge`,
+which calculates the Peierls phases for you::
 
-  def hopping(site_a, site_b, t=1):
-    return -t
+  import numpy as np
+  import kwant
 
-  syst[...] = onsite
-  syst[...] = hopping
+  def hopping(a, b, t, peierls):
+      return -t * peierls(a, b)
 
-  # Raises ValueError
+  syst = make_system(hopping)
+  lead = make_lead(hopping).substituted(peierls='peierls_lead')
+  syst.attach_lead(lead)
   syst = syst.finalized()
 
-As a solution, simply remove the default values and always provide ``t``.
-To deal with many parameters, the following idiom may be useful::
-
-  defaults = dict(a=0, b=1, c=2, d=3)
-  ...
-  smatrix = kwant.smatrix(syst, E, params=dict(defaults, d=4, e=5))
+  gauge = kwant.physics.magnetic_gauge(syst)
 
-Note that this allows to override defaults as well as to add additional
-parameters.
+  def B_syst(pos):
+     return np.exp(-np.sum(pos * pos))
 
-Installation on Microsoft Windows is available via Conda
---------------------------------------------------------
-Kwant is now packaged for the Conda package manager on Windows, and using
-Conda is the preferred method for installing Kwant on that platform.
-Please refer to the
-`installation section <https://kwant-project.org/install#microsoft-windows>`_
-of the Kwant website for details.
-Currently the MUMPS solver is not available for the Windows version of the
-Conda package; we hope to include MUMPS support in a later patch release.
+  # B_syst in scattering region, 0 in lead.
+  # Ensure that the fields match at the system/lead interface!
+  peierls_syst, peierls_lead = gauge(B_syst, 0)
 
-Minimum required versions for some dependencies have increased
---------------------------------------------------------------
-Kwant now requires at least the following versions:
+  params = dict(t=1, peierls=peierls_syst, peierls_lead=peierls_lead)
+  kwant.hamiltonian_submatrix(syst, params=params)
 
-+ Python 3.5
-+ numpy 0.11.0
-+ scipy 0.17.0
-+ matplotlib 1.5.1
+Note that the API for this functionality is provisional, and may be
+revised in a future version of Kwant.
 
-These versions (or newer) are available in the latest stable releases
-of Ubuntu and Debian GNU/Linux.
+.. _whatsnew14-kpm-conductivity:
 
 Conductivity calculations using `kwant.kpm.conductivity`
 --------------------------------------------------------
@@ -90,17 +75,7 @@ potentials at finite temperature::
   conductivities = [sigma_xy(mu=mu, temperature=0.1)
                     for mu in np.linspace(0, 4)]
 
-`kwant.physics.Bands` can optionally return eigenvectors and velocities
------------------------------------------------------------------------
-`kwant.physics.Bands` now takes extra parameters that allow it to
-return the mode eigenvectors, and also the derivatives of the dispersion
-relation (up to second order) using the Hellman-Feynman relation::
-
-  syst = make_system().finalized()
-
-  bands = kwant.physics.Bands(syst)
-  (energies, velocities, vectors) = bands(k=0, derivative_order=1,
-                                          return_eigenvectors=True)
+.. _whatsnew14-qsymm:
 
 Integration with Qsymm library
 ------------------------------
@@ -136,40 +111,7 @@ and vice versa, and for working with continuum Hamiltonians (such as would be us
 This integration requires separately installing Qsymm, which is available on the
 `Python Package Index <https://pypi.org/project/qsymm/>`_.
 
-Automatic Peierls phase calculation
------------------------------------
-When defining systems with orbital magnetic fields it is often cumbersome to
-manually calculate the phases required by the Peierls substitution, and to
-ensure that the chosen gauge is consistent across the whole system
-(this is especially true for systems with leads that point in different
-directions). This release introduces `kwant.physics.magnetic_gauge`,
-which calculates the Peierls phases for you::
-
-  import numpy as np
-  import kwant
-
-  def hopping(a, b, t, peierls):
-      return -t * peierls(a, b)
-
-  syst = make_system(hopping)
-  lead = make_lead(hopping).substituted(peierls='peierls_lead')
-  syst.attach_lead(lead)
-  syst = syst.finalized()
-
-  gauge = kwant.physics.magnetic_gauge(syst)
-
-  def B_syst(pos):
-     return np.exp(-np.sum(pos * pos))
-
-  # B_syst in scattering region, 0 in lead.
-  # Ensure that the fields match at the system/lead interface!
-  peierls_syst, peierls_lead = gauge(B_syst, 0)
-
-  params = dict(t=1, peierls=peierls_syst, peierls_lead=peierls_lead)
-  kwant.hamiltonian_submatrix(syst, params=params)
-
-Note that the API for this functionality is provisional, and may be
-revised in a future version of Kwant.
+.. _whatsnew14-parameters:
 
 System parameter substitution
 -----------------------------
@@ -248,27 +190,7 @@ a dictionary using ``params``::
 
 Providing ``args`` will be removed in a future Kwant version.
 
-Finalized Builders keep track of which sites were added when attaching leads
-----------------------------------------------------------------------------
-When attaching leads to an irregularly shaped scattering region, Kwant adds
-sites in order to make the interface with the leads "smooth". Previously,
-the information of which sites were added was not inspectable after finalization.
-Now the sites that were added from each lead are available in the ``lead_paddings``
-attribute. See the documentation for `~kwant.builder.FiniteSystem` for details.
-
-`kwant.continuum.discretize` can be used with rectangular lattices
-------------------------------------------------------------------
-Previously the discretizer could only be used with lattices with the same
-lattice constant in all directions. Now it is possible to pass rectangular
-lattices to the discretizer::
-
-  kwant.continuum.discretize(
-    'k_x**2 + k_y**2',
-    grid=kwant.lattice.general([(1, 0), (0, 2]),
-  )
-
-This is useful when you need a finer discretization step in some spatial
-directions, and a coarser one in others.
+.. _whatsnew14-density-plots:
 
 Interpolated density plots
 --------------------------
@@ -304,8 +226,119 @@ different families are plotted in different colors, which improves the
 default plotting style. You can still customize the site coloring using
 the ``site_color`` parameter, as before.
 
+`kwant.physics.Bands` can optionally return eigenvectors and velocities
+-----------------------------------------------------------------------
+`kwant.physics.Bands` now takes extra parameters that allow it to
+return the mode eigenvectors, and also the derivatives of the dispersion
+relation (up to second order) using the Hellman-Feynman relation::
+
+  syst = make_system().finalized()
+
+  bands = kwant.physics.Bands(syst)
+  (energies, velocities, vectors) = bands(k=0, derivative_order=1,
+                                          return_eigenvectors=True)
+
+Finalized Builders keep track of which sites were added when attaching leads
+----------------------------------------------------------------------------
+When attaching leads to an irregularly shaped scattering region, Kwant adds
+sites in order to make the interface with the leads "smooth". Previously,
+the information of which sites were added was not inspectable after finalization.
+Now the sites that were added from each lead are available in the ``lead_paddings``
+attribute. See the documentation for `~kwant.builder.FiniteSystem` for details.
+
+`kwant.continuum.discretize` can be used with rectangular lattices
+------------------------------------------------------------------
+Previously the discretizer could only be used with lattices with the same
+lattice constant in all directions. Now it is possible to pass rectangular
+lattices to the discretizer::
+
+  kwant.continuum.discretize(
+    'k_x**2 + k_y**2',
+    grid=kwant.lattice.general([(1, 0), (0, 2]),
+  )
+
+This is useful when you need a finer discretization step in some spatial
+directions, and a coarser one in others.
+
+Restrictions on value functions when named parameters are given
+---------------------------------------------------------------
+New restrictions apply to how value functions may accept arguments, when named
+parameters are given through ``params``.  (Nothing changes when the now
+deprcated ``args`` mechanism is used).  The restrictions follow the principle
+that each value function must take a clearly specified set of named parameters.
+This allows to make the parameter handling less error-prone and faster.
+
+In particular, when ``params`` is used, it is no longer possible for value
+functions to
+- take ``*args`` or ``**kwargs``,
+- take keyword-only parameters,
+- have default parameters for arguments.
+
+As an example, the following snippet no longer works because it uses default
+values::
+
+  syst = kwant.Builder()
+
+  # Parameter 't' has a default value of 1
+  def onsite(site, V, t=1):
+    return V = 2 * t
+
+  def hopping(site_a, site_b, t=1):
+    return -t
+
+  syst[...] = onsite
+  syst[...] = hopping
+
+  # Raises ValueError
+  syst = syst.finalized()
+
+As a solution, simply remove the default values and always provide ``t``.
+To deal with many parameters, the following idiom may be useful::
+
+  defaults = dict(a=0, b=1, c=2, d=3)
+  ...
+  smatrix = kwant.smatrix(syst, E, params=dict(defaults, d=4, e=5))
+
+Note that this allows to override defaults as well as to add additional
+parameters.
+
+.. _whatsnew14-windows:
+
+Installation on Microsoft Windows is available via Conda
+--------------------------------------------------------
+Kwant is now packaged for the Conda package manager on Windows, and using
+Conda is the preferred method for installing Kwant on that platform.
+Please refer to the
+`installation section <https://kwant-project.org/install#microsoft-windows>`_
+of the Kwant website for details.
+Currently the MUMPS solver is not available for the Windows version of the
+Conda package; we hope to include MUMPS support in a later patch release.
+
+Minimum required versions for some dependencies have increased
+--------------------------------------------------------------
+Kwant now requires at least the following versions:
+
++ Python 3.5
++ numpy 0.11.0
++ scipy 0.17.0
++ matplotlib 1.5.1
+
+These versions (or newer) are available in the latest stable releases
+of Ubuntu and Debian GNU/Linux.
+
 Changes in Kwant 1.4.1
 ----------------------
 - The list of user-visible changes was rearranged to emphasize
-  backwards-incompatible changes by moving them to the top of the list and
-  adding the entry `Value functions may no longer take unnamed arguments`_.
+  backwards-incompatible changes by moving them to the top of the list.
+- Restrictions on value functions no longer apply when the old ``args``
+  mechanism is used, this restores most of the backwards compatibility with
+  previous Kwant versions: `Restrictions on value functions when named
+  parameters are given`_.
+- The ``args`` parameter passing mechanism works again with
+  `~kwant.wraparound`-treated systems.  Some restriction continue to appply,
+  notably it is not possible to use ``wraparound`` with value functions that
+  take ``*args`` or ``*kwargs``.
+- Kwant no longer requires the existence of a `parameters` attribute for
+  low-level systems.
+- A note about an :ref:`whatsnew13-params-api-change` that ocurred in Kwant
+  1.3 was added.
diff --git a/kwant/_system.pyx b/kwant/_system.pyx
index 45e0da438e6b9a7322307331982ea76d8cf272c7..71860ed674dbabfade12a69c9eb2df9f31b6f65e 100644
--- a/kwant/_system.pyx
+++ b/kwant/_system.pyx
@@ -243,17 +243,6 @@ def make_dense_full(ham, args, params, CGraph gr, diag,
     return h_sub
 
 
-def _check_parameters_match(expected_parameters, params):
-    if params is None:
-        params = {}
-    missing = set(expected_parameters) - set(params)
-
-    if missing:
-        msg = ('System is missing required parameters: ',
-               ', '.join(map('"{}"'.format, missing)))
-        raise TypeError(''.join(msg))
-
-
 @deprecate_args
 @cython.binding(True)
 @cython.embedsignature(True)
@@ -301,9 +290,6 @@ def hamiltonian_submatrix(self, args=(), to_sites=None, from_sites=None,
     n = self.graph.num_nodes
     matrix = ta.matrix
 
-    if not args:  # Then perhaps parameters
-        _check_parameters_match(self.parameters, params)
-
     if from_sites is None:
         diag = n * [None]
         from_norb = np.empty(n, gint_dtype)
diff --git a/kwant/builder.py b/kwant/builder.py
index 8eb11c43652f60d6b9ae54358476cddccff9ab3d..ad482a75c838438b709296704ef3ff8ee2419022 100644
--- a/kwant/builder.py
+++ b/kwant/builder.py
@@ -1868,6 +1868,9 @@ class _FinalizedBuilderMixin:
             if param_names is not None:  # 'value' is callable
                 site = self.symmetry.to_fd(self.sites[i])
                 if params:
+                    # See body of _value_params_pair_cache().
+                    if isinstance(param_names, Exception):
+                        raise param_names
                     args = map(params.__getitem__, param_names)
                 try:
                     value = value(site, *args)
@@ -1891,6 +1894,9 @@ class _FinalizedBuilderMixin:
                 sites = self.sites
                 site_i, site_j = self.symmetry.to_fd(sites[i], sites[j])
                 if params:
+                    # See body of _value_params_pair_cache().
+                    if isinstance(param_names, Exception):
+                        raise param_names
                     args = map(params.__getitem__, param_names)
                 try:
                     value = value(site_i, site_j, *args)
@@ -1940,7 +1946,18 @@ def _value_params_pair_cache(nstrip):
             if isinstance(value, _Substituted):
                 entry = value.func, value.params[nstrip:]
             elif callable(value):
-                entry = value, get_parameters(value)[nstrip:]
+                try:
+                    param_names = get_parameters(value)
+                except ValueError as ex:
+                    # The parameter names are determined and stored in advance
+                    # for future use.  This has failed, but it will only turn
+                    # into a problem if user code ever uses the 'params'
+                    # mechanism.  To maintain backwards compatibility, we catch
+                    # and store the exception so that it can be raised whenever
+                    # appropriate.
+                    entry = value, ex
+                else:
+                    entry = value, param_names[nstrip:]
             else:
                 # None means: value is not callable. (That's faster to check.)
                 entry = value, None
@@ -2036,13 +2053,23 @@ class FiniteSystem(_FinalizedBuilderMixin, system.FiniteSystem):
         hoppings = [cache(builder._get_edge(sites[tail], sites[head]))
                     for tail, head in g]
 
-        # System parameters are the union of the parameters
-        # of onsites and hoppings.
-        # Here 'onsites' and 'hoppings' are pairs whos second element
-        # is a tuple of parameter names when matrix element is a function,
-        # and None otherwise.
-        parameters = frozenset(chain.from_iterable(
-            p for _, p in chain(onsites, hoppings) if p))
+        # Compute the union of the parameters of onsites and hoppings.  Here,
+        # 'onsites' and 'hoppings' are pairs whose second element is one of
+        # three things:
+        #
+        # * a tuple of parameter names when the matrix element is a function,
+        # * 'None' when it is a constant,
+        # * an exception when the parameter names could not have been
+        #   determined (See body of _value_params_pair_cache()).
+        parameters = []
+        for _, names in chain(onsites, hoppings):
+            if isinstance(names, Exception):
+                parameters = None
+                break
+            if names:
+                parameters.extend(names)
+        else:
+            parameters = frozenset(parameters)
 
         self.graph = g
         self.sites = sites
@@ -2205,13 +2232,23 @@ class InfiniteSystem(_FinalizedBuilderMixin, system.InfiniteSystem):
                 tail, head = sym.to_fd(tail, head)
             hoppings.append(cache(builder._get_edge(tail, head)))
 
-        # System parameters are the union of the parameters
-        # of onsites and hoppings.
-        # Here 'onsites' and 'hoppings' are pairs whos second element
-        # is a tuple of parameter names when matrix element is a function,
-        # and None otherwise.
-        parameters = frozenset(chain.from_iterable(
-            p for _, p in chain(onsites, hoppings) if p))
+        # Compute the union of the parameters of onsites and hoppings.  Here,
+        # 'onsites' and 'hoppings' are pairs whose second element is one of
+        # three things:
+        #
+        # * a tuple of parameter names when the matrix element is a function,
+        # * 'None' when it is a constant,
+        # * an exception when the parameter names could not have been
+        #   determined (See body of _value_params_pair_cache()).
+        parameters = []
+        for _, names in chain(onsites, hoppings):
+            if isinstance(names, Exception):
+                parameters = None
+                break
+            if names:
+                parameters.extend(names)
+        else:
+            parameters = frozenset(parameters)
 
         self.graph = g
         self.sites = sites
diff --git a/kwant/tests/test_wraparound.py b/kwant/tests/test_wraparound.py
index f9b536ec3ac696e955c179d0a236f5d4d08c68c2..aa53e1c3268043e8334bf8081754b5d24ef3f469 100644
--- a/kwant/tests/test_wraparound.py
+++ b/kwant/tests/test_wraparound.py
@@ -9,6 +9,7 @@
 import tempfile
 import itertools
 import numpy as np
+from numpy.testing import assert_equal
 import tinyarray as ta
 import pytest
 
@@ -380,3 +381,18 @@ def test_fd_mismatch():
                 for k in np.linspace(-np.pi, np.pi, 5)]
 
     assert np.allclose(spectrum1, spectrum2)
+
+
+# There seems no more specific way to only filter KwantDeprecationWarning.
+@pytest.mark.filterwarnings('ignore')
+def test_args_params_equivalence():
+    for lat in [kwant.lattice.square(), kwant.lattice.honeycomb(),
+                kwant.lattice.kagome()]:
+        syst = kwant.Builder(kwant.TranslationalSymmetry(*lat.prim_vecs))
+        syst[lat.shape((lambda pos: True), (0, 0))] = 1
+        syst[lat.neighbors(1)] = 0.1
+        syst[lat.neighbors(2)] = lambda a, b, param: 0.01
+        syst = wraparound(syst).finalized()
+        shs = syst.hamiltonian_submatrix
+        assert_equal(shs(args=["bla", 1, 2]),
+                     shs(params=dict(param="bla", k_x=1, k_y=2)))
diff --git a/kwant/wraparound.py b/kwant/wraparound.py
index fdfcc9f2598714ad5189ceff5e50d6fc5440ecd7..c3c5a38cccb9f3a8562bb894df52d7c791d3a955 100644
--- a/kwant/wraparound.py
+++ b/kwant/wraparound.py
@@ -165,42 +165,45 @@ def wraparound(builder, keep=None, *, coordinate_names='xyz'):
             params[name] = inspect.Parameter(
                 name, inspect.Parameter.POSITIONAL_ONLY)
 
-        # Add all the other parameters, except for the momenta.
+        # Add all the other parameters (except for the momenta).  Setup the
+        # 'selections'.
         selections = []
         for val in vals:
             if not callable(val):
                 selections.append(())
                 continue
             val_params = get_parameters(val)[num_sites:]
+            assert val_params[mnp:] == momenta
+            val_params = val_params[:mnp]
             selections.append((*site_params, *val_params))
             for p in val_params:
-                # Skip parameters that exist in previously added functions,
-                # and the momenta, which will be placed at the end.
-                if p in params or p in momenta:
+                # Skip parameters that exist in previously added functions.
+                if p in params:
                     continue
                 params[p] = inspect.Parameter(
                     p, inspect.Parameter.POSITIONAL_ONLY)
 
-        # Finally, add the momenta.
-        for k in momenta:
-            params[k] = inspect.Parameter(
-                k, inspect.Parameter.POSITIONAL_ONLY)
-
         # Sort values such that ones with the same arguments are bunched.
         # Prepare 'val_selection_pairs' that is used in the function 'f' above.
         params_keys = list(params.keys())
         val_selection_pairs = []
         prev_selection = None
         argsort = sorted(range(len(selections)), key=selections.__getitem__)
+        momenta_sel = tuple(range(mnp, 0, 1))
         for i in argsort:
             selection = selections[i]
             if selection and selection != prev_selection:
                 prev_selection = selection = tuple(
-                    params_keys.index(s) for s in selection)
+                    params_keys.index(s) for s in selection) + momenta_sel
             else:
                 selection = ()
             val_selection_pairs.append((vals[i], selection))
 
+        # Finally, add the momenta.
+        for k in momenta:
+            params[k] = inspect.Parameter(
+                k, inspect.Parameter.POSITIONAL_ONLY)
+
         f.__signature__ = inspect.Signature(params.values())
         return f
 
@@ -220,7 +223,7 @@ def wraparound(builder, keep=None, *, coordinate_names='xyz'):
         sym = TranslationalSymmetry(*periods)
         momenta.pop(keep)
     momenta = tuple(momenta)
-    mnp = -len(sym.periods)      # Used by the bound functions above.
+    mnp = -len(momenta)         # Used by the bound functions above.
 
     # Store the names of the momentum parameters and the symmetry of the
     # old Builder (this will be needed for band structure plotting)