diff --git a/nifty5/domain_tuple.py b/nifty5/domain_tuple.py
index bbe7cff14fe0ccfe5f2072dd4fff8216984cb458..ef6bdd996b3b94b9be8580c40b7a5e5619104a5b 100644
--- a/nifty5/domain_tuple.py
+++ b/nifty5/domain_tuple.py
@@ -104,6 +104,16 @@ class DomainTuple(object):
         """
         return self._shape
 
+    @property
+    def local_shape(self):
+        """tuple of int: number of pixels along each axis on the local task
+
+        The shape of the array-like object required to store information
+        living on part of the domain which is stored on the local MPI task.
+        """
+        from .dobj import local_shape
+        return local_shape(self._shape)
+
     @property
     def size(self):
         """int : total number of pixels.
diff --git a/nifty5/domains/domain.py b/nifty5/domains/domain.py
index 38db023bd752b25501e8b2029bfdc85e77e971c3..75002a8dae88204446258c9fb77ffcc103b782f7 100644
--- a/nifty5/domains/domain.py
+++ b/nifty5/domains/domain.py
@@ -88,6 +88,16 @@ class Domain(NiftyMetaBase()):
         """
         raise NotImplementedError
 
+    @property
+    def local_shape(self):
+        """tuple of int: number of pixels along each axis on the local task
+
+        The shape of the array-like object required to store information
+        living on part of the domain which is stored on the local MPI task.
+        """
+        from ..dobj import local_shape
+        return local_shape(self.shape)
+
     @abc.abstractproperty
     def size(self):
         """int: total number of pixels.
diff --git a/nifty5/domains/log_rg_space.py b/nifty5/domains/log_rg_space.py
index 16c644cdb922d8f416b2a29dd5827c126a6779fc..511d990ef4d983c9658a09c66bcd05afacb4fd59 100644
--- a/nifty5/domains/log_rg_space.py
+++ b/nifty5/domains/log_rg_space.py
@@ -3,7 +3,7 @@ from ..sugar import exp
 
 import numpy as np
 
-from ..dobj import ibegin
+from .. import dobj
 from ..field import Field
 from .structured_domain import StructuredDomain
 
@@ -62,26 +62,22 @@ class LogRGSpace(StructuredDomain):
                           np.zeros(len(self.shape)), True)
 
     def get_k_length_array(self):
-        out = Field(self, dtype=np.float64)
-        oloc = out.local_data
-        ib = ibegin(out.val)
-        res = np.arange(oloc.shape[0], dtype=np.float64) + ib[0]
+        ib = dobj.ibegin_from_shape(self._shape)
+        res = np.arange(self.local_shape[0], dtype=np.float64) + ib[0]
         res = np.minimum(res, self.shape[0]-res)*self.bindistances[0]
         if len(self.shape) == 1:
-            oloc[()] = res
-            return out
+            return Field.from_local_data(self, res)
         res *= res
         for i in range(1, len(self.shape)):
-            tmp = np.arange(oloc.shape[i], dtype=np.float64) + ib[i]
+            tmp = np.arange(self.local_shape[i], dtype=np.float64) + ib[i]
             tmp = np.minimum(tmp, self.shape[i]-tmp)*self.bindistances[i]
             tmp *= tmp
             res = np.add.outer(res, tmp)
-        oloc[()] = np.sqrt(res)
-        return out
+        return Field.from_local_data(self, np.sqrt(res))
 
     def get_expk_length_array(self):
         # FIXME This is a hack! Only for plotting. Seems not to be the final version.
-        out = exp(self.get_k_length_array())
-        out.val[1:] = out.val[:-1]
-        out.val[0] = 0
-        return out
+        out = exp(self.get_k_length_array()).to_global_data().copy()
+        out[1:] = out[:-1]
+        out[0] = 0
+        return Field.from_global_data(self, out)
diff --git a/nifty5/domains/rg_space.py b/nifty5/domains/rg_space.py
index 118ee47dbeb95ffd2aefa5811ac3e8ea34220496..d3db6a98f08910c18ef5351f2a0637805887cb01 100644
--- a/nifty5/domains/rg_space.py
+++ b/nifty5/domains/rg_space.py
@@ -95,22 +95,18 @@ class RGSpace(StructuredDomain):
     def get_k_length_array(self):
         if (not self.harmonic):
             raise NotImplementedError
-        out = Field(self, dtype=np.float64)
-        oloc = out.local_data
-        ibegin = dobj.ibegin(out.val)
-        res = np.arange(oloc.shape[0], dtype=np.float64) + ibegin[0]
+        ibegin = dobj.ibegin_from_shape(self._shape)
+        res = np.arange(self.local_shape[0], dtype=np.float64) + ibegin[0]
         res = np.minimum(res, self.shape[0]-res)*self.distances[0]
         if len(self.shape) == 1:
-            oloc[()] = res
-            return out
+            return Field.from_local_data(self, res)
         res *= res
         for i in range(1, len(self.shape)):
-            tmp = np.arange(oloc.shape[i], dtype=np.float64) + ibegin[i]
+            tmp = np.arange(self.local_shape[i], dtype=np.float64) + ibegin[i]
             tmp = np.minimum(tmp, self.shape[i]-tmp)*self.distances[i]
             tmp *= tmp
             res = np.add.outer(res, tmp)
-        oloc[()] = np.sqrt(res)
-        return out
+        return Field.from_local_data(self, np.sqrt(res))
 
     def get_unique_k_lengths(self):
         if (not self.harmonic):
@@ -147,8 +143,7 @@ class RGSpace(StructuredDomain):
         from ..sugar import exp
         tmp = x*x
         tmp *= -2.*np.pi*np.pi*sigma*sigma
-        exp(tmp, out=tmp)
-        return tmp
+        return exp(tmp)
 
     def get_fft_smoothing_kernel_function(self, sigma):
         if (not self.harmonic):
diff --git a/nifty5/extra/operator_tests.py b/nifty5/extra/operator_tests.py
index 45b8437e5d81c69f527c4eefdbbde9c46ab2611c..39db35a4489ce5e8171682d30eddd4e561e32a8a 100644
--- a/nifty5/extra/operator_tests.py
+++ b/nifty5/extra/operator_tests.py
@@ -35,9 +35,9 @@ def adjoint_implementation(op, domain_dtype, target_dtype, atol, rtol):
     needed_cap = op.TIMES | op.ADJOINT_TIMES
     if (op.capability & needed_cap) != needed_cap:
         return
-    f1 = from_random("normal", op.domain, dtype=domain_dtype).lock()
-    f2 = from_random("normal", op.target, dtype=target_dtype).lock()
-    res1 = f1.vdot(op.adjoint_times(f2).lock())
+    f1 = from_random("normal", op.domain, dtype=domain_dtype)
+    f2 = from_random("normal", op.target, dtype=target_dtype)
+    res1 = f1.vdot(op.adjoint_times(f2))
     res2 = op.times(f1).vdot(f2)
     np.testing.assert_allclose(res1, res2, atol=atol, rtol=rtol)
 
@@ -46,12 +46,12 @@ def inverse_implementation(op, domain_dtype, target_dtype, atol, rtol):
     needed_cap = op.TIMES | op.INVERSE_TIMES
     if (op.capability & needed_cap) != needed_cap:
         return
-    foo = from_random("normal", op.target, dtype=target_dtype).lock()
-    res = op(op.inverse_times(foo).lock())
+    foo = from_random("normal", op.target, dtype=target_dtype)
+    res = op(op.inverse_times(foo))
     _assert_allclose(res, foo, atol=atol, rtol=rtol)
 
-    foo = from_random("normal", op.domain, dtype=domain_dtype).lock()
-    res = op.inverse_times(op(foo).lock())
+    foo = from_random("normal", op.domain, dtype=domain_dtype)
+    res = op.inverse_times(op(foo))
     _assert_allclose(res, foo, atol=atol, rtol=rtol)
 
 
diff --git a/nifty5/field.py b/nifty5/field.py
index 980019f857d5650098969c4d648dbf52fae25193..2c8aa84d5cee3a1b088af2fa9c95e37cb3e5b7ea 100644
--- a/nifty5/field.py
+++ b/nifty5/field.py
@@ -35,7 +35,7 @@ class Field(object):
     ----------
     domain : None, DomainTuple, tuple of Domain, or Domain
 
-    val : None, Field, data_object, or scalar
+    val : Field, data_object or scalar
         The values the array should contain after init. A scalar input will
         fill the whole array with this scalar. If a data_object is provided,
         its dimensions must match the domain's.
@@ -49,32 +49,30 @@ class Field(object):
     many convenience functions for Field conatruction!
     """
 
-    def __init__(self, domain=None, val=None, dtype=None, copy=False,
-                 locked=False):
+    def __init__(self, domain=None, val=None, dtype=None):
         self._domain = self._infer_domain(domain=domain, val=val)
 
         dtype = self._infer_dtype(dtype=dtype, val=val)
         if isinstance(val, Field):
             if self._domain != val._domain:
                 raise ValueError("Domain mismatch")
-            self._val = dobj.from_object(val.val, dtype=dtype, copy=copy,
-                                         set_locked=locked)
+            self._val = val._val
+
         elif (np.isscalar(val)):
             self._val = dobj.full(self._domain.shape, dtype=dtype,
                                   fill_value=val)
         elif isinstance(val, dobj.data_object):
             if self._domain.shape == val.shape:
-                self._val = dobj.from_object(val, dtype=dtype, copy=copy,
-                                             set_locked=locked)
+                if dtype == val.dtype:
+                    self._val = val
+                else:
+                    self._val = dobj.from_object(val, dtype, True, True)
             else:
                 raise ValueError("Shape mismatch")
-        elif val is None:
-            self._val = dobj.empty(self._domain.shape, dtype=dtype)
         else:
             raise TypeError("unknown source type")
 
-        if locked:
-            dobj.lock(self._val)
+        dobj.lock(self._val)
 
     # prevent implicit conversion to bool
     def __nonzero__(self):
@@ -84,7 +82,7 @@ class Field(object):
         raise TypeError("Field does not support implicit conversion to bool")
 
     @staticmethod
-    def full(domain, val, dtype=None):
+    def full(domain, val):
         """Creates a Field with a given domain, filled with a constant value.
 
         Parameters
@@ -101,11 +99,7 @@ class Field(object):
         """
         if not np.isscalar(val):
             raise TypeError("val must be a scalar")
-        return Field(DomainTuple.make(domain), val, dtype)
-
-    @staticmethod
-    def empty(domain, dtype=None):
-        return Field(DomainTuple.make(domain), None, dtype)
+        return Field(DomainTuple.make(domain), val)
 
     @staticmethod
     def from_global_data(domain, arr, sum_up=False):
@@ -152,11 +146,6 @@ class Field(object):
 
         Returns a handle to the part of the array data residing on the local
         task (or to the entore array if MPI is not active).
-
-        Notes
-        -----
-        If the field is not locked, the array data can be modified.
-        Use with care!
         """
         return dobj.local_data(self._val)
 
@@ -196,8 +185,6 @@ class Field(object):
             return dtype
         if val is None:
             raise ValueError("could not infer dtype")
-        if isinstance(val, Field):
-            return val.dtype
         return np.result_type(val)
 
     @staticmethod
@@ -223,41 +210,6 @@ class Field(object):
                      val=dobj.from_random(random_type, dtype=dtype,
                                           shape=domain.shape, **kwargs))
 
-    def fill(self, fill_value):
-        """Fill `self` uniformly with `fill_value`
-
-        Parameters
-        ----------
-        fill_value: float or complex or int
-            The value to fill the field with.
-        """
-        self._val.fill(fill_value)
-        return self
-
-    def lock(self):
-        """Write-protect the data content of `self`.
-
-        After this call, it will no longer be possible to change the data
-        entries of `self`. This is convenient if, for example, a
-        DiagonalOperator wants to ensure that its diagonal cannot be modified
-        inadvertently, without making copies.
-
-        Notes
-        -----
-        This will not only prohibit modifications to the entries of `self`, but
-        also to the entries of any other Field or numpy array pointing to the
-        same data. If an unlocked instance is needed, use copy().
-
-        The fact that there is no `unlock()` method is deliberate.
-        """
-        dobj.lock(self._val)
-        return self
-
-    @property
-    def locked(self):
-        """bool : True iff the field's data content has been locked"""
-        return dobj.locked(self._val)
-
     @property
     def val(self):
         """dobj.data_object : the data object storing the field's entries
@@ -303,43 +255,6 @@ class Field(object):
             raise ValueError(".imag called on a non-complex Field")
         return Field(self._domain, self.val.imag)
 
-    def copy(self):
-        """ Returns a full copy of the Field.
-
-        The returned object will be an identical copy of the original Field.
-        The copy will be writeable, even if `self` was locked.
-
-        Returns
-        -------
-        Field
-            An identical, but unlocked copy of 'self'.
-        """
-        return Field(val=self, copy=True)
-
-    def empty_copy(self):
-        """ Returns a Field with identical domain and data type, but
-        uninitialized data.
-
-        Returns
-        -------
-        Field
-            A copy of 'self', with uninitialized data.
-        """
-        return Field(self._domain, dtype=self.dtype)
-
-    def locked_copy(self):
-        """ Returns a read-only version of the Field.
-
-        If `self` is locked, returns `self`. Otherwise returns a locked copy
-        of `self`.
-
-        Returns
-        -------
-        Field
-            A read-only version of `self`.
-        """
-        return self if self.locked else Field(val=self, copy=True, locked=True)
-
     def scalar_weight(self, spaces=None):
         """Returns the uniform volume element for a sub-domain of `self`.
 
@@ -392,7 +307,7 @@ class Field(object):
             res *= self._domain[i].total_volume
         return res
 
-    def weight(self, power=1, spaces=None, out=None):
+    def weight(self, power=1, spaces=None):
         """ Weights the pixels of `self` with their invidual pixel-volume.
 
         Parameters
@@ -404,21 +319,12 @@ class Field(object):
             Determines on which sub-domain the operation takes place.
             If None, the entire domain is used.
 
-        out : Field or None
-            if not None, the result is returned in a new Field
-            otherwise the contents of "out" are overwritten with the result.
-            "out" may be identical to "self"!
-
         Returns
         -------
         Field
             The weighted field.
         """
-        if out is None:
-            out = self.copy()
-        else:
-            if out is not self:
-                out.copy_content_from(self)
+        aout = self.local_data.copy()
 
         spaces = utilities.parse_spaces(spaces, len(self._domain))
 
@@ -435,12 +341,12 @@ class Field(object):
                 if dobj.distaxis(self._val) >= 0 and ind == 0:
                     # we need to distribute the weights along axis 0
                     wgt = dobj.local_data(dobj.from_global_data(wgt))
-                out.local_data[()] *= wgt**power
+                aout *= wgt**power
         fct = fct**power
         if fct != 1.:
-            out *= fct
+            aout *= fct
 
-        return out
+        return Field.from_local_data(self._domain, aout)
 
     def vdot(self, x=None, spaces=None):
         """ Computes the dot product of 'self' with x.
@@ -508,7 +414,7 @@ class Field(object):
     # ---General unary/contraction methods---
 
     def __pos__(self):
-        return self.copy()
+        return self
 
     def __neg__(self):
         return Field(self._domain, -self.val)
@@ -538,7 +444,7 @@ class Field(object):
                                   for i, dom in enumerate(self._domain)
                                   if i not in spaces)
 
-            return Field(domain=return_domain, val=data, copy=False)
+            return Field(domain=return_domain, val=data)
 
     def sum(self, spaces=None):
         """Sums up over the sub-domains given by `spaces`.
@@ -713,13 +619,6 @@ class Field(object):
             return self._contraction_helper('std', spaces)
         return sqrt(self.var(spaces))
 
-    def copy_content_from(self, other):
-        if not isinstance(other, Field):
-            raise TypeError("argument must be a Field")
-        if other._domain != self._domain:
-            raise ValueError("domains are incompatible.")
-        self.local_data[()] = other.local_data[()]
-
     def __repr__(self):
         return "<nifty5.Field>"
 
@@ -745,13 +644,13 @@ class Field(object):
         return self.isEquivalentTo(other)
 
 
-for op in ["__add__", "__radd__", "__iadd__",
-           "__sub__", "__rsub__", "__isub__",
-           "__mul__", "__rmul__", "__imul__",
-           "__div__", "__rdiv__", "__idiv__",
-           "__truediv__", "__rtruediv__", "__itruediv__",
-           "__floordiv__", "__rfloordiv__", "__ifloordiv__",
-           "__pow__", "__rpow__", "__ipow__",
+for op in ["__add__", "__radd__",
+           "__sub__", "__rsub__",
+           "__mul__", "__rmul__",
+           "__div__", "__rdiv__",
+           "__truediv__", "__rtruediv__",
+           "__floordiv__", "__rfloordiv__",
+           "__pow__", "__rpow__",
            "__lt__", "__le__", "__gt__", "__ge__", "__eq__", "__ne__"]:
     def func(op):
         def func2(self, other):
@@ -761,11 +660,11 @@ for op in ["__add__", "__radd__", "__iadd__",
                 if other._domain != self._domain:
                     raise ValueError("domains are incompatible.")
                 tval = getattr(self.val, op)(other.val)
-                return self if tval is self.val else Field(self._domain, tval)
+                return Field(self._domain, tval)
 
             if np.isscalar(other) or isinstance(other, dobj.data_object):
                 tval = getattr(self.val, op)(other)
-                return self if tval is self.val else Field(self._domain, tval)
+                return Field(self._domain, tval)
 
             return NotImplemented
         return func2
diff --git a/nifty5/minimization/conjugate_gradient.py b/nifty5/minimization/conjugate_gradient.py
index 0c3d828b6469b10a0636f6acd6ef89804b2b311d..94cd8991e9f74bf908df26f923f79d6c2ed03251 100644
--- a/nifty5/minimization/conjugate_gradient.py
+++ b/nifty5/minimization/conjugate_gradient.py
@@ -66,7 +66,7 @@ class ConjugateGradient(Minimizer):
             return energy, status
 
         r = energy.gradient
-        d = r.copy() if preconditioner is None else preconditioner(r)
+        d = r if preconditioner is None else preconditioner(r)
 
         previous_gamma = r.vdot(d).real
         if previous_gamma == 0:
diff --git a/nifty5/minimization/energy.py b/nifty5/minimization/energy.py
index 2edf844e0209ee62d7f2c133e401f9ff6496b65e..07df8fb6c5a876a2f09c1e680ca0b9be0f725b5d 100644
--- a/nifty5/minimization/energy.py
+++ b/nifty5/minimization/energy.py
@@ -52,7 +52,7 @@ class Energy(NiftyMetaBase()):
 
     def __init__(self, position):
         super(Energy, self).__init__()
-        self._position = position.lock()
+        self._position = position
 
     def at(self, position):
         """ Returns a new Energy object, initialized at `position`.
diff --git a/nifty5/minimization/energy_sum.py b/nifty5/minimization/energy_sum.py
index d90e9d134eca40cc9d3cd85a889d39b6f7e192e1..6a5f14af8408b7c4f4da363d6bfe56db38d46074 100644
--- a/nifty5/minimization/energy_sum.py
+++ b/nifty5/minimization/energy_sum.py
@@ -63,7 +63,7 @@ class EnergySum(Energy):
     @memo
     def gradient(self):
         return my_lincomb(map(lambda v: v.gradient, self._energies),
-                          self._factors).lock()
+                          self._factors)
 
     @property
     @memo
diff --git a/nifty5/minimization/line_energy.py b/nifty5/minimization/line_energy.py
index 2ef3518fe66012954f2b8163ea5a05873f9bb8bd..477ad01df6bbc5bbddfc841569420929774a6717 100644
--- a/nifty5/minimization/line_energy.py
+++ b/nifty5/minimization/line_energy.py
@@ -48,7 +48,7 @@ class LineEnergy(object):
     def __init__(self, line_position, energy, line_direction, offset=0.):
         super(LineEnergy, self).__init__()
         self._line_position = float(line_position)
-        self._line_direction = line_direction.lock()
+        self._line_direction = line_direction
 
         if self._line_position == float(offset):
             self._energy = energy
diff --git a/nifty5/minimization/quadratic_energy.py b/nifty5/minimization/quadratic_energy.py
index 87f7e63fb14c690cfaaba28ad00fee57052934cd..2467d6c5d554ddf1cf027c2d19c6cb1ec31160f1 100644
--- a/nifty5/minimization/quadratic_energy.py
+++ b/nifty5/minimization/quadratic_energy.py
@@ -35,7 +35,6 @@ class QuadraticEnergy(Energy):
         else:
             Ax = self._A(self.position)
             self._grad = Ax if b is None else Ax - b
-        self._grad.lock()
         self._value = 0.5*self.position.vdot(Ax)
         if b is not None:
             self._value -= b.vdot(self.position)
diff --git a/nifty5/minimization/scipy_minimizer.py b/nifty5/minimization/scipy_minimizer.py
index f3b712f916e9a8d13e3b5c0f5cd8049c479c845f..3cf3c1d6e74d3de888bc146f8d0729d705e7d3bb 100644
--- a/nifty5/minimization/scipy_minimizer.py
+++ b/nifty5/minimization/scipy_minimizer.py
@@ -33,7 +33,7 @@ def _toFlatNdarray(fld):
 
 
 def _toField(arr, dom):
-    return Field.from_global_data(dom, arr.reshape(dom.shape))
+    return Field.from_global_data(dom, arr.reshape(dom.shape).copy())
 
 
 class _MinHelper(object):
@@ -44,7 +44,7 @@ class _MinHelper(object):
     def _update(self, x):
         pos = _toField(x, self._domain)
         if (pos != self._energy.position).any():
-            self._energy = self._energy.at(pos.locked_copy())
+            self._energy = self._energy.at(pos)
 
     def fun(self, x):
         self._update(x)
diff --git a/nifty5/minimization/vl_bfgs.py b/nifty5/minimization/vl_bfgs.py
index 501fba2bc5867d7f9186e8df8f9548e81e8648cd..67e67ca925dba209b8124805a40b52cd606f07b7 100644
--- a/nifty5/minimization/vl_bfgs.py
+++ b/nifty5/minimization/vl_bfgs.py
@@ -109,8 +109,8 @@ class _InformationStore(object):
         self.max_history_length = max_history_length
         self.s = [None]*max_history_length
         self.y = [None]*max_history_length
-        self.last_x = x0.copy()
-        self.last_gradient = gradient.copy()
+        self.last_x = x0
+        self.last_gradient = gradient
         self.k = 0
 
         mmax = max_history_length
@@ -233,7 +233,7 @@ class _InformationStore(object):
         self.s[self.k % mmax] = x - self.last_x
         self.y[self.k % mmax] = gradient - self.last_gradient
 
-        self.last_x = x.copy()
-        self.last_gradient = gradient.copy()
+        self.last_x = x
+        self.last_gradient = gradient
 
         self.k += 1
diff --git a/nifty5/multi/multi_field.py b/nifty5/multi/multi_field.py
index f8074b03c15eaa75b4ae75b8ae0d015bce4a4cce..18011add3443c9d73c16ebc94b4dadc5bcb95c7e 100644
--- a/nifty5/multi/multi_field.py
+++ b/nifty5/multi/multi_field.py
@@ -69,18 +69,6 @@ class MultiField(object):
                                                   dtype[key], **kwargs)
                            for key in sorted(domain.keys())})
 
-    def fill(self, fill_value):
-        """Fill `self` uniformly with `fill_value`
-
-        Parameters
-        ----------
-        fill_value: float or complex or int
-            The value to fill the field with.
-        """
-        for val in self._val.values():
-            val.fill(fill_value)
-        return self
-
     def _check_domain(self, other):
         if other._domain != self._domain:
             raise ValueError("domains are incompatible.")
@@ -92,27 +80,6 @@ class MultiField(object):
             result += sub_field.vdot(x[key])
         return result
 
-    def lock(self):
-        for v in self.values():
-            v.lock()
-        return self
-
-    @property
-    def locked(self):
-        return all(v.locked for v in self.values())
-
-    def copy(self):
-        return MultiField({key: val.copy() for key, val in self.items()})
-
-    def locked_copy(self):
-        if self.locked:
-            return self
-        return MultiField({key: val.locked_copy()
-                          for key, val in self.items()})
-
-    def empty_copy(self):
-        return MultiField({key: val.empty_copy() for key, val in self.items()})
-
     @staticmethod
     def build_dtype(dtype, domain):
         if isinstance(dtype, dict):
@@ -121,12 +88,6 @@ class MultiField(object):
             dtype = np.float64
         return {key: dtype for key in domain.keys()}
 
-    @staticmethod
-    def empty(domain, dtype=None):
-        dtype = MultiField.build_dtype(dtype, domain)
-        return MultiField({key: Field.empty(dom, dtype=dtype[key])
-                           for key, dom in domain.items()})
-
     @staticmethod
     def full(domain, val):
         return MultiField({key: Field.full(dom, val)
@@ -241,9 +202,9 @@ for op in ["__add__", "__radd__",
                         result_val[key] = getattr(self[key], op)(other[key])
                     if op in ("__add__", "__radd__"):
                         for key in only_self_keys:
-                            result_val[key] = self[key].copy()
+                            result_val[key] = self[key]
                         for key in only_other_keys:
-                            result_val[key] = other[key].copy()
+                            result_val[key] = other[key]
                     elif op in ("__mul__", "__rmul__"):
                         pass
                     else:
diff --git a/nifty5/operators/diagonal_operator.py b/nifty5/operators/diagonal_operator.py
index ec6df7ce2f2a85b15fd57c2c9d21188ed1e5940c..b619083e9bec89fd91609c9996590bc1732c42bd 100644
--- a/nifty5/operators/diagonal_operator.py
+++ b/nifty5/operators/diagonal_operator.py
@@ -185,7 +185,7 @@ class DiagonalOperator(EndomorphicOperator):
         res = Field.from_random(random_type="normal", domain=self._domain,
                                 dtype=dtype)
         if from_inverse:
-            res.local_data[()] /= np.sqrt(self._ldiag)
+            res /= np.sqrt(self._ldiag)
         else:
-            res.local_data[()] *= np.sqrt(self._ldiag)
+            res *= np.sqrt(self._ldiag)
         return res
diff --git a/nifty5/operators/dof_distributor.py b/nifty5/operators/dof_distributor.py
index 10418a1b521fd38717282b4276e45782028f47f5..c53ac3362ab41c15b66cba1fa140070117a4a512 100644
--- a/nifty5/operators/dof_distributor.py
+++ b/nifty5/operators/dof_distributor.py
@@ -120,15 +120,14 @@ class DOFDistributor(LinearOperator):
         return res
 
     def _times(self, x):
-        res = Field.empty(self._target, dtype=x.dtype)
         if dobj.distaxis(x.val) in x.domain.axes[self._space]:
             arr = x.to_global_data()
         else:
             arr = x.local_data
         arr = arr.reshape(self._hshape)
+        oarr = arr[(slice(None), self._dofdex, slice(None))]
+        return Field.from_local_data(self._target, oarr.reshape(self._target.local_shape))
         oarr = res.local_data.reshape(self._pshape)
-        oarr[()] = arr[(slice(None), self._dofdex, slice(None))]
-        return res
 
     def apply(self, x, mode):
         self._check_input(x, mode)
diff --git a/nifty5/operators/inversion_enabler.py b/nifty5/operators/inversion_enabler.py
index 9980d8492f0781d436232712f63d5f38ae8eeb66..5efd5848b57248acaaa9ac705db856713c885d19 100644
--- a/nifty5/operators/inversion_enabler.py
+++ b/nifty5/operators/inversion_enabler.py
@@ -23,6 +23,7 @@ from ..minimization.conjugate_gradient import ConjugateGradient
 from ..minimization.iteration_controller import IterationController
 from ..minimization.quadratic_energy import QuadraticEnergy
 from .endomorphic_operator import EndomorphicOperator
+from ..sugar import full
 
 
 class InversionEnabler(EndomorphicOperator):
@@ -65,7 +66,7 @@ class InversionEnabler(EndomorphicOperator):
         if self._op.capability & mode:
             return self._op.apply(x, mode)
 
-        x0 = x.empty_copy().fill(0.)
+        x0 = full(x.domain, 0.)
         invmode = self._modeTable[self.INVERSE_BIT][self._ilog[mode]]
         invop = self._op._flip_modes(self._ilog[invmode])
         prec = self._approximation
diff --git a/nifty5/operators/scaling_operator.py b/nifty5/operators/scaling_operator.py
index b0aacb7d75c5c872cd3c40ec1ed692479ea47fbf..7f74dab68ab8eab5aa32b4ff80602083dddcc363 100644
--- a/nifty5/operators/scaling_operator.py
+++ b/nifty5/operators/scaling_operator.py
@@ -22,6 +22,7 @@ from ..field import Field
 from ..multi.multi_field import MultiField
 from .endomorphic_operator import EndomorphicOperator
 from ..domain_tuple import DomainTuple
+from ..sugar import full
 
 
 class ScalingOperator(EndomorphicOperator):
@@ -61,9 +62,9 @@ class ScalingOperator(EndomorphicOperator):
         self._check_input(x, mode)
 
         if self._factor == 1.:
-            return x.copy()
+            return x
         if self._factor == 0.:
-            return x.empty_copy().fill(0.)
+            return full(self.domain, 0.)
 
         if mode == self.TIMES:
             return x*self._factor
diff --git a/nifty5/operators/selection_operator.py b/nifty5/operators/selection_operator.py
index 8723d027bd14c7b17b72ca53790ac8b2a5002413..953dfae5dd269282821fe138c5ae421997fa6e85 100644
--- a/nifty5/operators/selection_operator.py
+++ b/nifty5/operators/selection_operator.py
@@ -50,10 +50,9 @@ class SelectionOperator(LinearOperator):
         return self.TIMES | self.ADJOINT_TIMES
 
     def apply(self, x, mode):
-        # FIXME Is the copying necessary?
         self._check_input(x, mode)
         if mode == self.TIMES:
-            return x[self._key].copy()
+            return x[self._key]
         else:
             from ..multi.multi_field import MultiField
-            return MultiField({self._key: x.copy()})
+            return MultiField({self._key: x})
diff --git a/nifty5/operators/symmetrizing_operator.py b/nifty5/operators/symmetrizing_operator.py
index a0e429f802671a2f7eeae89901f4eb7f3e7d9184..98d28d34e3223bf09b8c8f765d7c44af5e48a580 100644
--- a/nifty5/operators/symmetrizing_operator.py
+++ b/nifty5/operators/symmetrizing_operator.py
@@ -15,7 +15,7 @@ class SymmetrizingOperator(EndomorphicOperator):
 
     def apply(self, x, mode):
         self._check_input(x, mode)
-        tmp = x.copy().val
+        tmp = x.val.copy()
         ax = dobj.distaxis(tmp)
         globshape = tmp.shape
         for i in range(self._ndim):
diff --git a/nifty5/sugar.py b/nifty5/sugar.py
index 4afa3a98dcc0c0b21530dd2efd3b4105a16cd856..b5069e2bac77a17713a19bcaa87a8848fcdf6102 100644
--- a/nifty5/sugar.py
+++ b/nifty5/sugar.py
@@ -31,7 +31,7 @@ from .logger import logger
 
 __all__ = ['PS_field', 'power_analyze', 'create_power_operator',
            'create_harmonic_smoothing_operator', 'from_random',
-           'full', 'empty', 'from_global_data', 'from_local_data',
+           'full', 'from_global_data', 'from_local_data',
            'makeDomain', 'sqrt', 'exp', 'log', 'tanh', 'conjugate',
            'get_signal_variance', 'makeOp']
 
@@ -203,12 +203,6 @@ def full(domain, val):
     return Field.full(domain, val)
 
 
-def empty(domain, dtype):
-    if isinstance(domain, (dict, MultiDomain)):
-        return MultiField.empty(domain, dtype)
-    return Field.empty(domain, dtype)
-
-
 def from_random(random_type, domain, dtype=np.float64, **kwargs):
     if isinstance(domain, (dict, MultiDomain)):
         return MultiField.from_random(random_type, domain, dtype, **kwargs)
@@ -248,26 +242,13 @@ _current_module = sys.modules[__name__]
 
 for f in ["sqrt", "exp", "log", "tanh", "conjugate"]:
     def func(f):
-        def func2(x, out=None):
+        def func2(x):
             if isinstance(x, MultiField):
-                if out is not None:
-                    if (not isinstance(out, MultiField) or
-                            x._domain != out._domain):
-                        raise ValueError("Bad 'out' argument")
-                    for key, value in x.items():
-                        func2(value, out=out[key])
-                    return out
                 return MultiField({key: func2(val) for key, val in x.items()})
             elif isinstance(x, Field):
                 fu = getattr(dobj, f)
-                if out is not None:
-                    if not isinstance(out, Field) or x._domain != out._domain:
-                        raise ValueError("Bad 'out' argument")
-                    fu(x.val, out=out.val)
-                    return out
-                else:
-                    return Field(domain=x._domain, val=fu(x.val))
+                return Field(domain=x._domain, val=fu(x.val))
             else:
-                return getattr(np, f)(x, out)
+                return getattr(np, f)(x)
         return func2
     setattr(_current_module, f, func(f))
diff --git a/test/test_field.py b/test/test_field.py
index 1bfad9e90de2ad89484ba85e3a735d6f9c6d9ba2..c5687e2babe6a45ca807a09dfa258333d4bdd5f7 100644
--- a/test/test_field.py
+++ b/test/test_field.py
@@ -124,21 +124,6 @@ class Test_Functionality(unittest.TestCase):
         res = m.vdot(m, spaces=1)
         assert_allclose(res.local_data, 37.5)
 
-    def test_lock(self):
-        s1 = ift.RGSpace((10,))
-        f1 = ift.Field(s1, 27)
-        assert_equal(f1.locked, False)
-        f1.lock()
-        assert_equal(f1.locked, True)
-        with assert_raises(ValueError):
-            f1 += f1
-        assert_equal(f1.locked_copy() is f1, True)
-
-    def test_fill(self):
-        s1 = ift.RGSpace((10,))
-        f1 = ift.Field(s1, 27)
-        assert_equal(f1.fill(10).local_data, 10)
-
     def test_dataconv(self):
         s1 = ift.RGSpace((10,))
         ld = np.arange(ift.dobj.local_shape(s1.shape)[0])
@@ -158,12 +143,6 @@ class Test_Functionality(unittest.TestCase):
         assert_equal(f.local_data, 5)
         f = ift.Field(None, 5)
         assert_equal(f.local_data, 5)
-        assert_equal(f.empty_copy().domain, f.domain)
-        assert_equal(f.empty_copy().dtype, f.dtype)
-        assert_equal(f.copy().domain, f.domain)
-        assert_equal(f.copy().dtype, f.dtype)
-        assert_equal(f.copy().local_data, f.local_data)
-        assert_equal(f.copy() is f, False)
 
     def test_trivialities(self):
         s1 = ift.RGSpace((10,))
@@ -182,8 +161,7 @@ class Test_Functionality(unittest.TestCase):
     def test_weight(self):
         s1 = ift.RGSpace((10,))
         f = ift.Field(s1, 10.)
-        f2 = f.copy()
-        f.weight(1, out=f2)
+        f2 = f.weight(1)
         assert_equal(f.weight(1).local_data, f2.local_data)
         assert_equal(f.total_volume(), 1)
         assert_equal(f.total_volume(0), 1)
@@ -233,10 +211,6 @@ class Test_Functionality(unittest.TestCase):
             f1.vdot(42)
         with assert_raises(ValueError):
             f1.vdot(ift.Field(s2, 1.))
-        with assert_raises(TypeError):
-            f1.copy_content_from(1)
-        with assert_raises(ValueError):
-            f1.copy_content_from(ift.Field(s2, 1.))
         with assert_raises(TypeError):
             ift.full(s1, [2, 3])
 
@@ -246,9 +220,6 @@ class Test_Functionality(unittest.TestCase):
         assert_equal(f.local_data, 27)
         assert_equal(f.shape, (200,))
         assert_equal(f.dtype, np.int)
-        fx = ift.empty(f.domain, f.dtype)
-        assert_equal(f.dtype, fx.dtype)
-        assert_equal(f.shape, fx.shape)
         fx = ift.full(f.domain, 0)
         assert_equal(f.dtype, fx.dtype)
         assert_equal(f.shape, fx.shape)
diff --git a/test/test_multi_field.py b/test/test_multi_field.py
index 79fff128bb0df7a1d8fe7def0ed2f28e49608742..9409b3f6a4b94950531f4fa296baeb77d13ec0bd 100644
--- a/test/test_multi_field.py
+++ b/test/test_multi_field.py
@@ -32,19 +32,6 @@ class Test_Functionality(unittest.TestCase):
         f2 = ift.from_random("normal", domain=dom, dtype=np.complex128)
         assert_allclose(f1.vdot(f2), np.conj(f2.vdot(f1)))
 
-    def test_lock(self):
-        f1 = ift.full(dom, 27)
-        assert_equal(f1.locked, False)
-        f1.lock()
-        assert_equal(f1.locked, True)
-        assert_equal(f1.locked_copy() is f1, True)
-
-    def test_fill(self):
-        f1 = ift.full(dom, 27)
-        f1.fill(10)
-        for val in f1.values():
-            assert_equal((val == 10).all(), True)
-
     def test_dataconv(self):
         f1 = ift.full(dom, 27)
         f2 = ift.from_global_data(dom, f1.to_global_data())