Basic operators from functions to classes (#496)

* Changed Operator functions to classes

* Updated directionalderivative.py

* Updated basicoperators

* Added name

* Updated docstrings

* Updated docstrings

* Made methods protected

* Made the calc static

* Made the changes

* Updated Gradient

* Updated Gradient

* minor: quick-fix

* minor: quick-fix

* minor: updated second_ddop

Author: rohanbabbar04

pylops/basicoperators/block.py

@@ -1,38 +1,45 @@
 __all__ = ["Block"]
 
+import multiprocessing as mp
+
 from typing import Iterable, Optional
 
 from pylops import LinearOperator
 from pylops.basicoperators import HStack, VStack
-from pylops.utils.typing import DTypeLike
+from pylops.utils.typing import DTypeLike, NDArray
 
 
-def _Block(
-    ops: Iterable[Iterable[LinearOperator]],
-    dtype: Optional[DTypeLike] = None,
-    _HStack=HStack,
-    _VStack=VStack,
-    args_HStack: Optional[dict] = None,
-    args_VStack: Optional[dict] = None,
-) -> LinearOperator:
+class _Block(LinearOperator):
     """Block operator.
 
     Used to be able to provide operators from different libraries to
     Block.
     """
-    if args_HStack is None:
-        args_HStack = {}
-    if args_VStack is None:
-        args_VStack = {}
-    hblocks = [_HStack(hblock, dtype=dtype, **args_HStack) for hblock in ops]
-    return _VStack(hblocks, dtype=dtype, **args_VStack)
-
-
-def Block(
-    ops: Iterable[Iterable[LinearOperator]],
-    nproc: int = 1,
-    dtype: Optional[DTypeLike] = None,
-) -> LinearOperator:
+    def __init__(self, ops: Iterable[Iterable[LinearOperator]],
+                 dtype: Optional[DTypeLike] = None,
+                 _HStack=HStack,
+                 _VStack=VStack,
+                 args_HStack: Optional[dict] = None,
+                 args_VStack: Optional[dict] = None, name: str = 'B'):
+        if args_HStack is None:
+            self.args_HStack = {}
+        else:
+            self.args_HStack = args_HStack
+        if args_VStack is None:
+            self.args_VStack = {}
+        else:
+            self.args_VStack = args_VStack
+        hblocks = [_HStack(hblock, dtype=dtype, **self.args_HStack) for hblock in ops]
+        super().__init__(Op=_VStack(ops=hblocks, dtype=dtype, **self.args_VStack), name=name)
+
+    def _matvec(self, x: NDArray) -> NDArray:
+        return super()._matvec(x)
+
+    def _rmatvec(self, x: NDArray) -> NDArray:
+        return super()._rmatvec(x)
+
+
+class Block(_Block):
     r"""Block operator.
 
     Create a block operator from N lists of M linear operators each.
@@ -116,4 +123,9 @@ def Block(
         \end{bmatrix}
 
     """
-    return _Block(ops, dtype=dtype, args_VStack={"nproc": nproc})
+    def __init__(self, ops: Iterable[Iterable[LinearOperator]],
+                 nproc: int = 1,
+                 dtype: Optional[DTypeLike] = None):
+        if nproc > 1:
+            self.pool = mp.Pool(processes=nproc)
+        super().__init__(ops=ops, dtype=dtype, args_VStack={"nproc": nproc})

pylops/basicoperators/directionalderivative.py

@@ -8,14 +8,7 @@
 from pylops.utils.typing import DTypeLike, InputDimsLike, NDArray
 
 
-def FirstDirectionalDerivative(
-    dims: InputDimsLike,
-    v: NDArray,
-    sampling: int = 1,
-    edge: bool = False,
-    kind: str = "centered",
-    dtype: DTypeLike = "float64",
-) -> LinearOperator:
+class FirstDirectionalDerivative(LinearOperator):
     r"""First Directional derivative.
 
     Apply a directional derivative operator to a multi-dimensional array
@@ -42,11 +35,6 @@ def FirstDirectionalDerivative(
     dtype : :obj:`str`, optional
         Type of elements in input array.
 
-    Returns
-    -------
-    ddop : :obj:`pylops.LinearOperator`
-        First directional derivative linear operator
-
     Notes
     -----
     The FirstDirectionalDerivative applies a first-order derivative
@@ -71,25 +59,39 @@ def FirstDirectionalDerivative(
     operator with :math:`\mathbf{v}` along the main diagonal.
 
     """
-    Gop = Gradient(dims, sampling=sampling, edge=edge, kind=kind, dtype=dtype)
-    if v.ndim == 1:
-        Dop = Diagonal(v, dims=[len(dims)] + list(dims), axis=0, dtype=dtype)
-    else:
-        Dop = Diagonal(v.ravel(), dtype=dtype)
-    Sop = Sum(dims=[len(dims)] + list(dims), axis=0, dtype=dtype)
-    ddop = Sop * Dop * Gop
-    ddop.dims = ddop.dimsd = dims
-    ddop.sampling = sampling
-    return ddop
-
-
-def SecondDirectionalDerivative(
-    dims: InputDimsLike,
-    v: NDArray,
-    sampling: int = 1,
-    edge: bool = False,
-    dtype: DTypeLike = "float64",
-) -> LinearOperator:
+
+    def __init__(self, dims: InputDimsLike,
+                 v: NDArray,
+                 sampling: int = 1,
+                 edge: bool = False,
+                 kind: str = "centered",
+                 dtype: DTypeLike = "float64",
+                 name: str = 'F'):
+        self.sampling = sampling
+        self.edge = edge
+        self.kind = kind
+        self.v = v
+        Op = self._calc_first_ddop(dims=dims, sampling=sampling, edge=edge, kind=kind, dtype=dtype, v=v)
+        super().__init__(Op=Op, name=name)
+
+    def _matvec(self, x: NDArray) -> NDArray:
+        return super()._matvec(x)
+
+    def _rmatvec(self, x: NDArray) -> NDArray:
+        return super()._rmatvec(x)
+
+    @staticmethod
+    def _calc_first_ddop(dims: InputDimsLike, v: NDArray, sampling: int, edge: bool, kind: str, dtype: DTypeLike):
+        Gop = Gradient(dims, sampling=sampling, edge=edge, kind=kind, dtype=dtype)
+        if v.ndim == 1:
+            Dop = Diagonal(v, dims=[len(dims)] + list(dims), axis=0, dtype=dtype)
+        else:
+            Dop = Diagonal(v.ravel(), dtype=dtype)
+        Sop = Sum(dims=[len(dims)] + list(dims), axis=0, dtype=dtype)
+        return Sop * Dop * Gop
+
+
+class SecondDirectionalDerivative(LinearOperator):
     r"""Second Directional derivative.
 
     Apply a second directional derivative operator to a multi-dimensional array
@@ -114,11 +116,6 @@ def SecondDirectionalDerivative(
     dtype : :obj:`str`, optional
         Type of elements in input array.
 
-    Returns
-    -------
-    ddop : :obj:`pylops.LinearOperator`
-        Second directional derivative linear operator
-
     Notes
     -----
     The SecondDirectionalDerivative applies a second-order derivative
@@ -135,8 +132,24 @@ def SecondDirectionalDerivative(
     This operator is sometimes also referred to as directional Laplacian
     in the literature.
     """
-    Dop = FirstDirectionalDerivative(dims, v, sampling=sampling, edge=edge, dtype=dtype)
-    ddop = -Dop.H * Dop
-    ddop.dims = ddop.dimsd = dims
-    ddop.sampling = sampling
-    return ddop
+
+    def __init__(self, dims: InputDimsLike, v: NDArray, sampling: int = 1, edge: bool = False,
+                 dtype: DTypeLike = "float64", name: str = 'S'):
+        self.dims = dims
+        self.v = v
+        self.sampling = sampling
+        self.edge = edge
+        Op = self._calc_second_ddop(dims=dims, v=v, sampling=sampling, edge=edge, dtype=dtype)
+        super().__init__(Op=Op, name=name)
+
+    def _matvec(self, x: NDArray) -> NDArray:
+        return super()._matvec(x)
+
+    def _rmatvec(self, x: NDArray) -> NDArray:
+        return super()._rmatvec(x)
+
+    @staticmethod
+    def _calc_second_ddop(dims: InputDimsLike, v: NDArray, sampling: int, edge: bool, dtype: DTypeLike):
+        Dop = FirstDirectionalDerivative(dims=dims, v=v, sampling=sampling, edge=edge, dtype=dtype)
+        ddop = -Dop.H * Dop
+        return ddop

pylops/basicoperators/gradient.py

@@ -5,16 +5,10 @@
 from pylops import LinearOperator
 from pylops.basicoperators import FirstDerivative, VStack
 from pylops.utils._internal import _value_or_sized_to_tuple
-from pylops.utils.typing import DTypeLike, InputDimsLike
+from pylops.utils.typing import DTypeLike, InputDimsLike, NDArray
 
 
-def Gradient(
-    dims: Union[int, InputDimsLike],
-    sampling: int = 1,
-    edge: bool = False,
-    kind: str = "centered",
-    dtype: DTypeLike = "float64",
-) -> LinearOperator:
+class Gradient(LinearOperator):
     r"""Gradient.
 
     Apply gradient operator to a multi-dimensional array.
@@ -36,11 +30,6 @@ def Gradient(
     dtype : :obj:`str`, optional
         Type of elements in input array.
 
-    Returns
-    -------
-    l2op : :obj:`pylops.LinearOperator`
-        Gradient linear operator
-
     Notes
     -----
     The Gradient operator applies a first-order derivative to each dimension of
@@ -69,26 +58,33 @@ def Gradient(
     axes are instead summed together.
 
     """
-    dims = _value_or_sized_to_tuple(dims)
-    ndims = len(dims)
-    sampling = _value_or_sized_to_tuple(sampling, repeat=ndims)
-
-    gop = VStack(
-        [
-            FirstDerivative(
-                dims,
-                axis=iax,
-                sampling=sampling[iax],
-                edge=edge,
-                kind=kind,
-                dtype=dtype,
-            )
+
+    def __init__(self,
+                 dims: Union[int, InputDimsLike],
+                 sampling: int = 1,
+                 edge: bool = False,
+                 kind: str = "centered",
+                 dtype: DTypeLike = "float64", name: str = 'G'):
+        dims = _value_or_sized_to_tuple(dims)
+        ndims = len(dims)
+        sampling = _value_or_sized_to_tuple(sampling, repeat=ndims)
+        self.sampling = sampling
+        self.edge = edge
+        self.kind = kind
+        Op = VStack([FirstDerivative(
+            dims=dims,
+            axis=iax,
+            sampling=sampling[iax],
+            edge=edge,
+            kind=kind,
+            dtype=dtype,
+        )
             for iax in range(ndims)
-        ]
-    )
-    gop.dims = dims
-    gop.dimsd = (ndims, *gop.dims)
-    gop.sampling = sampling
-    gop.edge = edge
-    gop.kind = kind
-    return gop
+        ])
+        super().__init__(Op=Op, dims=dims, dimsd=(ndims, *dims), name=name)
+
+    def _matvec(self, x: NDArray) -> NDArray:
+        return super()._matvec(x)
+
+    def _rmatvec(self, x: NDArray) -> NDArray:
+        return super()._rmatvec(x)

pylops/basicoperators/laplacian.py

@@ -2,23 +2,16 @@
 
 
 from typing import Tuple
+from pylops.utils.typing import NDArray
 
 from numpy.core.multiarray import normalize_axis_index
 
-from pylops import LinearOperator, aslinearoperator
+from pylops import LinearOperator
 from pylops.basicoperators import SecondDerivative
 from pylops.utils.typing import DTypeLike, InputDimsLike
 
 
-def Laplacian(
-    dims: InputDimsLike,
-    axes: InputDimsLike = (-2, -1),
-    weights: Tuple[float, ...] = (1, 1),
-    sampling: Tuple[float, ...] = (1, 1),
-    edge: bool = False,
-    kind: str = "centered",
-    dtype: DTypeLike = "float64",
-) -> LinearOperator:
+class Laplacian(LinearOperator):
     r"""Laplacian.
 
     Apply second-order centered Laplacian operator to a multi-dimensional array.
@@ -47,11 +40,6 @@ def Laplacian(
     dtype : :obj:`str`, optional
         Type of elements in input array.
 
-    Returns
-    -------
-    l2op : :obj:`pylops.LinearOperator`
-        Laplacian linear operator
-
     Raises
     ------
     ValueError
@@ -69,27 +57,42 @@ def Laplacian(
                   / (\Delta x \Delta y)
 
     """
-    axes = tuple(normalize_axis_index(ax, len(dims)) for ax in axes)
-    if not (len(axes) == len(weights) == len(sampling)):
-        raise ValueError("axes, weights, and sampling have different size")
-
-    l2op = SecondDerivative(
-        dims, axis=axes[0], sampling=sampling[0], edge=edge, kind=kind, dtype=dtype
-    )
-    dims, dimsd = l2op.dims, l2op.dimsd
-
-    l2op *= weights[0]
-    for ax, samp, weight in zip(axes[1:], sampling[1:], weights[1:]):
-        l2op += weight * SecondDerivative(
-            dims, axis=ax, sampling=samp, edge=edge, dtype=dtype
-        )
 
-    l2op = aslinearoperator(l2op)
-    l2op.dims = dims
-    l2op.dimsd = dimsd
-    l2op.axes = axes
-    l2op.weights = weights
-    l2op.sampling = sampling
-    l2op.edge = edge
-    l2op.kind = kind
-    return l2op
+    def __init__(self, dims: InputDimsLike,
+                 axes: InputDimsLike = (-2, -1),
+                 weights: Tuple[float, ...] = (1, 1),
+                 sampling: Tuple[float, ...] = (1, 1),
+                 edge: bool = False,
+                 kind: str = "centered",
+                 dtype: DTypeLike = "float64", name: str = "L"):
+        axes = tuple(normalize_axis_index(ax, len(dims)) for ax in axes)
+        if not (len(axes) == len(weights) == len(sampling)):
+            raise ValueError("axes, weights, and sampling have different size")
+        self.axes = axes
+        self.weights = weights
+        self.sampling = sampling
+        self.edge = edge
+        self.kind = kind
+        Op = self._calc_l2op(dims=dims, axes=axes, sampling=sampling, edge=edge, kind=kind, dtype=dtype,
+                             weights=weights)
+        super().__init__(Op=Op, name=name)
+
+    def _matvec(self, x: NDArray) -> NDArray:
+        return super()._matvec(x)
+
+    def _rmatvec(self, x: NDArray) -> NDArray:
+        return super()._rmatvec(x)
+
+    @staticmethod
+    def _calc_l2op(dims: InputDimsLike, axes: InputDimsLike, weights: Tuple[float, ...], sampling: Tuple[float, ...],
+                   edge: bool, kind: str, dtype: DTypeLike):
+        l2op = SecondDerivative(
+            dims, axis=axes[0], sampling=sampling[0], edge=edge, kind=kind, dtype=dtype
+        )
+        dims = l2op.dims
+        l2op *= weights[0]
+        for ax, samp, weight in zip(axes[1:], sampling[1:], weights[1:]):
+            l2op += weight * SecondDerivative(
+                dims, axis=ax, sampling=samp, edge=edge, dtype=dtype
+            )
+        return l2op