Repository: apytypes/apytypes (language: python)

Issue:
Array creation methods
### From shape or value
- ~`empty`~
- ~`empty_like`~
- [x] `eye` (#491)
- [x] `identity` (#491)
- [x] `ones` (#491)
- [x] `ones_like` (#491)
- [x] `zeros` (#491)
- [x] `zeros_like` (#491)
- [x] `full` (#491)
- [x] `full_like` (#491)

`ones` and `zeros` should be static methods on each array class. `ones_like` and `zeros_like` can be methods on the array class and have Python methods. `full` and `full_like` should take both floats/ints and the corresponding APy-scalar and should probably warn if the APy-scalar is quantized (has a longer word length than the resulting array).

### From existing data
- ~`array`~
- ~`asarray`~
- ~`asanyarray`~
- ~`ascontiguousarray`~ ?
- ~`asmatrix`~
- ~`astype`~
- [ ] `copy` can have a Python function.
- ~`frombuffer`~
- ~`from_dlpack`~ ?
- [ ] `fromfile` 
- [ ] `fromfunction`
- [ ] `fromiter`
- [ ] `fromstring`
- [ ] `loadtxt`

### ~Creating record arrays~

### ~Creating character arrays~

### Numerical ranges
- [x] `arange` (#491,  #513)
- [ ] `linspace`
- [ ] `logspace`
- [ ] `geomspace`
- [ ] `meshgrid`
- [ ] `mgrid`
- [ ] `ogrid`

### Building matrices
- [ ] `diag`
- [ ] `diagflat`
- [ ] `tri`
- [ ] `tril`
- [ ] `triu`
- [ ] `vander`

### ~The matrix class~


Discussion / hints:
How should we handle bit sizes when receiving a Python integer or float? Should the user always specify the number of integer, fractional, exponent, and mantissa bits, or should we try to determine appropriate values dynamically? If the input value is an APy-Scalar, we could use the same bit sizes as the scalar for the array.

Fix patch:
diff --git a/CHANGELOG.md b/CHANGELOG.md
index 6c38e0f3b..e6a2d008d 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -42,6 +42,15 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
   - `min()`
   - `nanmax()`
   - `nanmin()`
+- Array creation functions added:
+  - `eye()`
+  - `identity()`
+  - `ones()`
+  - `ones_like()`
+  - `zeros()`
+  - `zeros_like()`
+  - `full()`
+  - `full_like()`
 - Added `array_utils.h` to contain shared functions between float and fixed arrays.
 - Added support for higher dimensions (ndim > 2) in `transpose()`.
 
diff --git a/docs/api/apyfixedarray.rst b/docs/api/apyfixedarray.rst
index ad7351a9e..c122b20ce 100644
--- a/docs/api/apyfixedarray.rst
+++ b/docs/api/apyfixedarray.rst
@@ -15,6 +15,19 @@
 
    .. automethod:: from_float
 
+   Other creation functions
+   ------------------------
+
+   .. automethod:: zeros
+
+   .. automethod:: ones
+
+   .. automethod:: eye
+
+   .. automethod:: identity
+
+   .. automethod:: full
+
    Change word length
    ------------------
 
diff --git a/docs/api/apyfloatarray.rst b/docs/api/apyfloatarray.rst
index 5028f64e5..26388c021 100644
--- a/docs/api/apyfloatarray.rst
+++ b/docs/api/apyfloatarray.rst
@@ -15,6 +15,19 @@
 
    .. automethod:: from_float
 
+   Other creation functions
+   ------------------------
+
+   .. automethod:: zeros
+
+   .. automethod:: ones
+
+   .. automethod:: eye
+
+   .. automethod:: identity
+
+   .. automethod:: full
+
    Change word length
    ------------------
 
diff --git a/docs/api/arrayfunctions.rst b/docs/api/arrayfunctions.rst
index b0aa0b651..0619e5618 100644
--- a/docs/api/arrayfunctions.rst
+++ b/docs/api/arrayfunctions.rst
@@ -14,3 +14,19 @@ Array functions
 .. autofunction:: apytypes.swapaxes
 
 .. autofunction:: apytypes.expand_dims
+
+.. autofunction:: apytypes.zeros
+
+.. autofunction:: apytypes.ones
+
+.. autofunction:: apytypes.eye
+
+.. autofunction:: apytypes.identity
+
+.. autofunction:: apytypes.full
+
+.. autofunction:: apytypes.zeros_like
+
+.. autofunction:: apytypes.ones_like
+
+.. autofunction:: apytypes.full_like
diff --git a/lib/apytypes/__init__.py b/lib/apytypes/__init__.py
index 601bbae4a..689061dfe 100644
--- a/lib/apytypes/__init__.py
+++ b/lib/apytypes/__init__.py
@@ -28,6 +28,15 @@
     moveaxis,
     swapaxes,
     expand_dims,
+    zeros,
+    zeros_like,
+    ones,
+    ones_like,
+    eye,
+    identity,
+    full,
+    full_like,
+    arange,
 )
 
 from apytypes._version import version as __version__
@@ -58,6 +67,15 @@
     "moveaxis",
     "swapaxes",
     "expand_dims",
+    "zeros",
+    "zeros_like",
+    "ones",
+    "ones_like",
+    "eye",
+    "identity",
+    "full",
+    "full_like",
+    "arange",
 ]
 
 APyFloat.__doc__ = r"""
diff --git a/lib/apytypes/__init__.pyi b/lib/apytypes/__init__.pyi
index e2497c146..dbdee15c0 100644
--- a/lib/apytypes/__init__.pyi
+++ b/lib/apytypes/__init__.pyi
@@ -23,15 +23,24 @@ from ._apytypes import (
     set_float_quantization_seed as set_float_quantization_seed,
 )
 from ._array_functions import (
+    arange as arange,
     convolve as convolve,
     expand_dims as expand_dims,
+    eye as eye,
+    full as full,
+    full_like as full_like,
+    identity as identity,
     moveaxis as moveaxis,
+    ones as ones,
+    ones_like as ones_like,
     ravel as ravel,
     reshape as reshape,
     shape as shape,
     squeeze as squeeze,
     swapaxes as swapaxes,
     transpose as transpose,
+    zeros as zeros,
+    zeros_like as zeros_like,
 )
 
 __all__: list = [
@@ -60,6 +69,15 @@ __all__: list = [
     "moveaxis",
     "swapaxes",
     "expand_dims",
+    "zeros",
+    "zeros_like",
+    "ones",
+    "ones_like",
+    "eye",
+    "identity",
+    "full",
+    "full_like",
+    "arange",
 ]
 
 annotations: __future__._Feature = ...
diff --git a/lib/apytypes/_apytypes.pyi b/lib/apytypes/_apytypes.pyi
index 65171b11a..346aa3343 100644
--- a/lib/apytypes/_apytypes.pyi
+++ b/lib/apytypes/_apytypes.pyi
@@ -1330,6 +1330,135 @@ class APyFixedArray:
         :class:`APyFixedArray`
         """
 
+    @staticmethod
+    def zeros(
+        shape: tuple,
+        int_bits: int | None = None,
+        frac_bits: int | None = None,
+        bits: int | None = None,
+    ) -> APyFixedArray:
+        """
+        Initializes an array with zeros.
+
+        Parameters
+        ----------
+        shape : tuple
+            Shape of the array.
+        int_bits : int, optional
+            Number of integer bits. Default is None.
+        frac_bits : int, optional
+            Number of fractional bits. Default is None.
+        bits : int, optional
+            Total number of bits. Default is None.
+
+        Returns
+        -------
+        APyFixedArray
+            An array initialized with zeros.
+        """
+
+    @staticmethod
+    def ones(
+        shape: tuple,
+        int_bits: int | None = None,
+        frac_bits: int | None = None,
+        bits: int | None = None,
+    ) -> APyFixedArray:
+        """
+        Initializes an array with ones.
+
+        Parameters
+        ----------
+        shape : tuple
+            Shape of the array.
+        int_bits : int, optional
+            Number of integer bits. Default is None.
+        frac_bits : int, optional
+            Number of fractional bits. Default is None.
+        bits : int, optional
+            Total number of bits. Default is None.
+
+        Returns
+        -------
+        APyFixedArray
+            An array initialized with ones.
+        """
+
+    @staticmethod
+    def eye(
+        n: int,
+        m: int | None = None,
+        int_bits: int | None = None,
+        frac_bits: int | None = None,
+        bits: int | None = None,
+    ) -> APyFixedArray:
+        """
+        Initializes an array with ones on the diagonal.
+
+        Parameters
+        ----------
+        n : int
+            Number of rows.
+        m : int, optional
+            Number of columns. Default is None.
+        int_bits : int, optional
+            Number of integer bits. Default is None.
+        frac_bits : int, optional
+            Number of fractional bits. Default is None.
+        bits : int, optional
+            Total number of bits. Default is None.
+
+        Returns
+        -------
+        APyFixedArray
+            An array with the specified value on the diagonal.
+        """
+
+    @staticmethod
+    def identity(
+        n: int,
+        int_bits: int | None = None,
+        frac_bits: int | None = None,
+        bits: int | None = None,
+    ) -> APyFixedArray:
+        """
+        Initializes an identity matrix with ones on the diagonal.
+
+        Parameters
+        ----------
+        n : int
+            Number of rows (and columns) in n x n output.
+        int_bits : int, optional
+            Number of integer bits. Default is None.
+        frac_bits : int, optional
+            Number of fractional bits. Default is None.
+        bits : int, optional
+            Total number of bits. Default is None.
+
+        Returns
+        -------
+        APyFixedArray
+            An identity matrix with ones on the diagonal.
+        """
+
+    @staticmethod
+    def full(shape: tuple, fill_value: APyFixed) -> APyFixedArray:
+        """
+        Initializes an array with the specified value.
+
+        Parameters
+        ----------
+        shape : tuple
+            Shape of the array.
+        fill_value : APyFixed
+            Value to fill the array.
+
+        Returns
+        -------
+        APyFixedArray
+            An array filled with the specified value.
+        """
+
     def __lshift__(self, shift_amnt: int) -> APyFixedArray: ...
     def __matmul__(self, rhs: APyFixedArray) -> APyFixedArray: ...
     def __repr__(self) -> str: ...
@@ -2390,6 +2519,126 @@ class APyFloatArray:
         :class:`APyFloatArray`
         """
 
+    @staticmethod
+    def zeros(
+        shape: tuple, exp_bits: int, man_bits: int, bias: int | None = None
+    ) -> APyFloatArray:
+        """
+        Initializes an array with zeros.
+
+        Parameters
+        ----------
+        shape : tuple
+            Shape of the array.
+        exp_bits : int
+            Number of exponent bits.
+        man_bits : int
+            Number of mantissa bits.
+        bias : optional
+            Set bias. Default is None.
+
+        Returns
+        -------
+        APyFloatArray
+            An array filled with zeros.
+        """
+
+    @staticmethod
+    def ones(
+        shape: tuple, exp_bits: int, man_bits: int, bias: int | None = None
+    ) -> APyFloatArray:
+        """
+        Initializes an array with ones.
+
+        Parameters
+        ----------
+        shape : tuple
+            Shape of the array.
+        exp_bits : int
+            Number of exponent bits.
+        man_bits : int
+            Number of mantissa bits.
+        bias : optional
+            Set bias. Default is None.
+
+        Returns
+        -------
+        APyFloatArray
+            An array filled with ones.
+        """
+
+    @staticmethod
+    def eye(
+        n: int,
+        exp_bits: int,
+        man_bits: int,
+        m: int | None = None,
+        bias: int | None = None,
+    ) -> APyFloatArray:
+        """
+        Initializes an array with ones on the diagonal.
+
+        Parameters
+        ----------
+        n : int
+            Number of rows (and columns) in the n x n output.
+        exp_bits : int
+            Number of exponent bits.
+        man_bits : int
+            Number of mantissa bits.
+        m : int, optional
+            Number of columns. Default is None, which results in an n x n output.
+        bias : optional
+            Set bias. Default is None.
+
+        Returns
+        -------
+        APyFloatArray
+            An array with the specified value on the diagonal.
+        """
+
+    @staticmethod
+    def identity(
+        n: int, exp_bits: int, man_bits: int, bias: int | None = None
+    ) -> APyFloatArray:
+        """
+        Initializes an identity matrix with ones on the diagonal.
+
+        Parameters
+        ----------
+        n : int
+            Number of rows (and columns) in the n x n output.
+        exp_bits : int
+            Number of exponent bits.
+        man_bits : int
+            Number of mantissa bits.
+        bias : optional
+            Set bias. Default is None.
+
+        Returns
+        -------
+        APyFloatArray
+            An identity matrix with ones on the diagonal.
+        """
+
+    @staticmethod
+    def full(shape: tuple, fill_value: APyFloat) -> APyFloatArray:
+        """
+        Initializes an array filled with the specified value.
+
+        Parameters
+        ----------
+        shape : tuple
+            Shape of the array.
+        fill_value : APyFloat
+            Value to fill the array.
+
+        Returns
+        -------
+        APyFloatArray
+            An array filled with the specified value.
+        """
+
     def __matmul__(self, rhs: APyFloatArray) -> APyFloatArray | APyFloat: ...
     def __repr__(self) -> str: ...
     def __len__(self) -> int: ...
diff --git a/lib/apytypes/_array_functions.py b/lib/apytypes/_array_functions.py
index 00b9ded37..f9e6bf467 100644
--- a/lib/apytypes/_array_functions.py
+++ b/lib/apytypes/_array_functions.py
@@ -411,9 +411,478 @@ def expand_dims(a, axis):
     return a.reshape(shape)
 
 
+# =============================================================================
+# array creation
+# =============================================================================
+
+
+def zeros(shape, int_bits=None, frac_bits=None, exp_bits=None, mantissa_bits=None):
+    """
+    Initializes an array with zeros.
+
+    Word lengths need to be specified and decide the return type of the array.
+
+    Parameters
+    ----------
+    shape : tuple
+        Shape of the array.
+    int_bits : int, optional
+        Number of integer bits.
+    frac_bits : int, optional
+        Number of fractional bits.
+    exp_bits : int, optional
+        Number of exponential bits.
+    mantissa_bits : int, optional
+        Number of mantissa bits.
+
+    Returns
+    -------
+    result : :class:`APyFloatArray` or :class:`APyFixedArray`
+        The initialized array with zeros.
+    """
+
+    from apytypes import APyFixedArray, APyFloatArray
+
+    a_type = _determine_array_type(int_bits, frac_bits, exp_bits, mantissa_bits)
+    try:
+        zeros = a_type.zeros
+    except AttributeError:
+        raise TypeError(f"Cannot make zeros array of type {type(a_type)}")
+    if a_type is APyFixedArray:
+        return zeros(shape=shape, int_bits=int_bits, frac_bits=frac_bits)
+
+    if a_type is APyFloatArray:
+        return zeros(shape=shape, exp_bits=exp_bits, man_bits=mantissa_bits)
+
+    raise ValueError("Only 'fixed' and 'float' array_types are defined")
+
+
+def ones(shape, int_bits=None, frac_bits=None, exp_bits=None, mantissa_bits=None):
+    """
+    Initializes an array with ones.
+
+    Word lengths need to be specified and decide the return type of the array.
+
+    Parameters
+    ----------
+    shape : tuple
+        Shape of the array.
+    int_bits : int, optional
+        Number of integer bits.
+    frac_bits : int, optional
+        Number of fractional bits.
+    exp_bits : int, optional
+        Number of exponential bits.
+    mantissa_bits : int, optional
+        Number of mantissa bits.
+
+    Returns
+    -------
+    result : :class:`APyFloatArray` or :class:`APyFixedArray`
+        The array initialized filled with ones.
+    """
+
+    from apytypes import APyFixedArray, APyFloatArray
+
+    a_type = _determine_array_type(int_bits, frac_bits, exp_bits, mantissa_bits)
+    try:
+        ones = a_type.ones
+    except AttributeError:
+        raise TypeError(f"Cannot make ones array of type {type(a_type)}")
+    if a_type is APyFixedArray or isinstance(a_type, APyFixedArray):
+        return ones(shape=shape, int_bits=int_bits, frac_bits=frac_bits)
+
+    if a_type is APyFloatArray or isinstance(a_type, APyFloatArray):
+        return ones(shape=shape, exp_bits=exp_bits, man_bits=mantissa_bits)
+
+    raise ValueError("Only 'fixed' and 'float' array_types are defined")
+
+
+def eye(
+    n: int, m=None, int_bits=None, frac_bits=None, exp_bits=None, mantissa_bits=None
+):
+    """
+    Return a 2-D array with ones on the diagonal and zeros elsewhere.
+
+    Word lengths need to be specified and decide the return type of the array.
+
+    Parameters
+    ----------
+    a_type : :class:`APyFloatArray` or :class:`APyFixedArray`, optional
+        The type of array to initialize.
+    n : int
+        Number of rows in the output.
+    m : int, optional
+        Number of columns in the output. If None, defaults to N.
+    int_bits : int, optional
+        Number of integer bits.
+    frac_bits : int, optional
+        Number of fractional bits.
+    exp_bits : int, optional
+        Number of exponential bits.
+    mantissa_bits : int, optional
+        Number of mantissa bits.
+
+    Returns
+    -------
+    result : :class:`APyFloatArray` or :class:`APyFixedArray`
+        An array where all elements are equal to zero, except for the k-th diagonal, whose values are equal to one.
+    """
+
+    from apytypes import APyFixedArray, APyFloatArray
+
+    a_type = _determine_array_type(int_bits, frac_bits, exp_bits, mantissa_bits)
+    try:
+        eye = a_type.eye
+    except AttributeError:
+        raise TypeError(f"Cannot make eye array of type {type(a_type)}")
+    if a_type is APyFixedArray or isinstance(a_type, APyFixedArray):
+        return eye(n=n, m=m, int_bits=int_bits, frac_bits=frac_bits)
+
+    if a_type is APyFloatArray or isinstance(a_type, APyFloatArray):
+        return eye(n=n, m=m, exp_bits=exp_bits, man_bits=mantissa_bits)
+
+    raise ValueError("Only 'fixed' and 'float' array_types are defined")
+
+
+def identity(n, int_bits=None, frac_bits=None, exp_bits=None, mantissa_bits=None):
+    """
+    Return the identity array.
+
+    Word lengths need to be specified and decides return type of Array.
+
+    Parameters
+    ----------
+    a_type : :class:`APyFloatArray` or :class:`APyFixedArray`, optional
+        The type of array to initialize.
+    n : int
+        Number of rows (and columns) in n x n output.
+    int_bits : int, optional
+        Number of integer bits.
+    frac_bits : int, optional
+        Number of fractional bits.
+    exp_bits : int, optional
+        Number of exponential bits.
+    mantissa_bits : int, optional
+        Number of mantissa bits.
+
+    Returns
+    -------
+    result : :class:`APyFloatArray` or :class:`APyFixedArray`
+        n x n array with its main diagonal set to one, and all other elements 0.
+    """
+
+    from apytypes import APyFixedArray, APyFloatArray
+
+    a_type = _determine_array_type(int_bits, frac_bits, exp_bits, mantissa_bits)
+    try:
+        identity = a_type.identity
+    except AttributeError:
+        raise TypeError(f"Cannot make identity array of type {type(a_type)}")
+
+    if a_type is APyFixedArray or isinstance(a_type, APyFixedArray):
+        return identity(n=n, int_bits=int_bits, frac_bits=frac_bits)
+
+    if a_type is APyFloatArray or isinstance(a_type, APyFloatArray):
+        return identity(n=n, exp_bits=exp_bits, man_bits=mantissa_bits)
+
+    raise ValueError("Only 'fixed' and 'float' array_types are defined")
+
+
+def full(
+    shape,
+    fill_value,
+    int_bits=None,
+    frac_bits=None,
+    exp_bits=None,
+    mantissa_bits=None,
+):
+    """
+    Return a new array of given shape and type, filled with fill_value.
+
+    If fill_value is an int or float, you must specify the word lengths (int_bits, frac_bits or exp_bits, mantissa_bits).
+    If fill_value is an APyFloat or APyFixed, the array will use the provided word lengths if specified.
+    If no word lengths are specified, the resulting array will inherit the word lengths from fill_value.
+
+    Parameters
+    ----------
+    shape : tuple
+        Shape of the array.
+    fill_value : :class:`APyFloat` or :class:`APyFixedArray` or int or float
+        Fill value.
+    int_bits : int, optional
+        Number of integer bits.
+    frac_bits : int, optional
+        Number of fractional bits.
+    exp_bits : int, optional
+        Number of exponential bits.
+    mantissa_bits : int, optional
+        Number of mantissa bits.
+
+    Returns
+    -------
+    result : :class:`APyFloatArray` or :class:`APyFixedArray`
+        Array of fill_value with the given shape.
+    """
+
+    from apytypes import APyFixedArray, APyFixed, APyFloat, APyFloatArray
+
+    fill_value = _normalize_fill_value(
+        fill_value, int_bits, frac_bits, exp_bits, mantissa_bits
+    )
+    try:
+        if isinstance(fill_value, APyFixed):
+            full = APyFixedArray.full
+        elif isinstance(fill_value, APyFloat):
+            full = APyFloatArray.full
+        else:
+            raise (AttributeError)
+    except AttributeError:
+        raise TypeError(f"Cannot make full array of type {type(fill_value)}")
+    return full(shape, fill_value)
+
+
+def zeros_like(a, int_bits=None, frac_bits=None, exp_bits=None, mantissa_bits=None):
+    """
+    Return an array of zeros with the same shape and type as a given array.
+    Defaults to `a` wordlength
+
+    Parameters
+    ----------
+    a : :class:`APyFloatArray` or :class:`APyFixedArray`
+        The shape and data-type define these same attributes of the returned array.
+    int_bits : int, optional
+        Number of integer bits.
+    frac_bits : int, optional
+        Number of fractional bits.
+    exp_bits : int, optional
+        Number of exponential bits.
+    mantissa_bits : int, optional
+        Number of mantissa bits.
+
+    Returns
+    -------
+    result : :class:`APyFloatArray` or :class:`APyFixedArray`
+        The initialized array with zeros.
+    """
+    try:
+        zeros = a.zeros
+    except AttributeError:
+        raise TypeError(f"Cannot make zeros array of type {type(a)}")
+
+    from apytypes import APyFixedArray, APyFloatArray
+
+    if isinstance(a, APyFixedArray):
+        int_bits = int_bits if int_bits is not None else a.int_bits
+        frac_bits = frac_bits if int_bits is not None else a.frac_bits
+        return zeros(shape=a.shape, int_bits=int_bits, frac_bits=frac_bits)
+
+    if isinstance(a, APyFloatArray):
+        exp_bits = exp_bits if exp_bits is not None else a.exp_bits
+        mantissa_bits = mantissa_bits if mantissa_bits is not None else a.man_bits
+        return zeros(shape=a.shape, exp_bits=exp_bits, man_bits=mantissa_bits)
+
+    raise ValueError("Only 'fixed' and 'float' array_types are defined")
+
+
+def ones_like(a, int_bits=None, frac_bits=None, exp_bits=None, mantissa_bits=None):
+    """
+    Return an array of ones with the same shape and type as a given array.
+    Defaults to `a` wordlength
+
+    Parameters
+    ----------
+    a : :class:`APyFloatArray` or :class:`APyFixedArray`
+        The shape and data-type of a define these same attributes of the returned array.
+    int_bits : int, optional
+        Number of integer bits.
+    frac_bits : int, optional
+        Number of fractional bits.
+    exp_bits : int, optional
+        Number of exponential bits.
+    mantissa_bits : int, optional
+        Number of mantissa bits.
+
+    Returns
+    -------
+    result : :class:`APyFloatArray` or :class:`APyFixedArray`
+        The array initialized filled with ones.
+    """
+    try:
+        ones = a.ones
+    except AttributeError:
+        raise TypeError(f"Cannot make ones array of type {type(a)}")
+
+    from apytypes import APyFixedArray, APyFloatArray
+
+    if isinstance(a, APyFixedArray):
+        int_bits = int_bits if int_bits is not None else a.int_bits
+        frac_bits = frac_bits if frac_bits is not None else a.frac_bits
+        return ones(shape=a.shape, int_bits=int_bits, frac_bits=frac_bits)
+
+    if isinstance(a, APyFloatArray):
+        exp_bits = exp_bits if exp_bits is not None else a.exp_bits
+        mantissa_bits = mantissa_bits if mantissa_bits is not None else a.man_bits
+        return ones(shape=a.shape, exp_bits=exp_bits, man_bits=mantissa_bits)
+
+    raise ValueError("Only 'fixed' and 'float' array_types are defined")
+
+
+def full_like(
+    a,
+    fill_value,
+    int_bits=None,
+    frac_bits=None,
+    exp_bits=None,
+    mantissa_bits=None,
+):
+    """
+    Return a full array with the same shape and type as a given array.
+
+    If fill_value is an int or float, you must specify the word lengths (int_bits, frac_bits or exp_bits, mantissa_bits).
+    If fill_value is an APyFloat or APyFixed, the array will use the provided word lengths if specified.
+    If no word lengths are specified, the resulting array will inherit the word lengths from fill_value.
+
+    Parameters
+    ----------
+    a : :class:`APyFloatArray` or :class:`APyFixedArray`
+        The shape and array type of the returned array.
+    fill_value : :class:`APyFloat` or :class:`APyFixed` or int or float
+        The value to fill the array with.
+    int_bits : int, optional
+        Number of integer bits for APyFixed.
+    frac_bits : int, optional
+        Number of fractional bits for APyFixed.
+    exp_bits : int, optional
+        Number of exponent bits for APyFloat.
+    mantissa_bits : int, optional
+        Number of mantissa bits for APyFloat.
+
+    Returns
+    -------
+    result : :class:`APyFloatArray` or :class:`APyFixedArray`
+        Array filled with fill_value, having the same shape and type as `a`.
+    """
+    try:
+        full = a.full
+    except AttributeError:
+        raise TypeError(f"Cannot make full array of type {type(a)}")
+
+    fill_value = _normalize_fill_value(
+        fill_value, int_bits, frac_bits, exp_bits, mantissa_bits
+    )
+    return full(a.shape, fill_value)
+
+
+def arange(
+    start,
+    stop=None,
+    step=1,
+    int_bits=None,
+    frac_bits=None,
+    exp_bits=None,
+    mantissa_bits=None,
+):
+    """
+    Simple and naive version of arange.
+
+    Parameters
+    ----------
+    start : int, float
+        start number
+    stop : int, optional
+        stop number
+    step : int, float, optional
+        step size in range
+    int_bits : int, optional
+        Number of integer bits for APyFixed.
+    frac_bits : int, optional
+        Number of fractional bits for APyFixed.
+    exp_bits : int, optional
+        Number of exponent bits for APyFloat.
+    mantissa_bits : int, optional
+        Number of mantissa bits for APyFloat.
+
+    Returns
+    -------
+    result : :class:`APyFloatArray` or :class:`APyFixedArray`
+        Array filled with fill_value, having the same shape and type as `a`.
+    """
+    if stop is None:
+        start, stop = 0, start
+    if step == 0:
+        raise ValueError("Step must not be zero")
+
+    floats = []
+    value = start
+    while (step > 0 and value < stop) or (step < 0 and value > stop):
+        floats.append(value)
+        value += step
+
+    a_type = _determine_array_type(int_bits, frac_bits, exp_bits, mantissa_bits)
+    from apytypes import APyFixedArray, APyFloatArray
+
+    if a_type is APyFixedArray:
+        return APyFixedArray.from_float(floats, int_bits=int_bits, frac_bits=frac_bits)
+    elif a_type is APyFloatArray:
+        return APyFloatArray.from_float(
+            floats, exp_bits=exp_bits, man_bits=mantissa_bits
+        )
+    raise ValueError("Undefined inputs")
+
+
 # =============================================================================
 # Helpers
 # =============================================================================
+def _determine_array_type(int_bits, frac_bits, exp_bits, mantissa_bits):
+    from apytypes import APyFixedArray, APyFloatArray
+
+    if int_bits is not None and frac_bits is not None:
+        return APyFixedArray
+    elif exp_bits is not None and mantissa_bits is not None:
+        return APyFloatArray
+    raise ValueError("You need to specify wordlengths in this function")
+
+
+def _normalize_fill_value(
+    fill_value,
+    int_bits=None,
+    frac_bits=None,
+    exp_bits=None,
+    mantissa_bits=None,
+):
+    from apytypes import APyFixed, APyFloat
+
+    if isinstance(fill_value, (int, float)):
+        if int_bits is not None and frac_bits is not None:
+            return APyFixed.from_float(
+                fill_value, int_bits=int_bits, frac_bits=frac_bits
+            )
+        elif exp_bits is not None and mantissa_bits is not None:
+            return APyFloat.from_float(
+                fill_value, exp_bits=exp_bits, man_bits=mantissa_bits
+            )
+        else:
+            raise ValueError(
+                "You need to specify wordlengths if input is a Python float or integer"
+            )
+
+    elif isinstance(fill_value, APyFixed):
+        int_bits = int_bits if int_bits is not None else fill_value.int_bits
+        frac_bits = frac_bits if frac_bits is not None else fill_value.frac_bits
+        return fill_value.cast(int_bits, frac_bits)
+
+    elif isinstance(fill_value, APyFloat):
+        exp_bits = exp_bits if exp_bits is not None else fill_value.exp_bits
+        mantissa_bits = (
+            mantissa_bits if mantissa_bits is not None else fill_value.man_bits
+        )
+        return fill_value.cast(exp_bits, mantissa_bits, fill_value.bias)
+
+    else:
+        raise ValueError(
+            "Only int, float, APyFloat and APyFixed are supported for fill_value"
+        )
 
 
 def _normalize_axis(axis: int, ndim: int) -> int:
diff --git a/lib/apytypes/_array_functions.pyi b/lib/apytypes/_array_functions.pyi
index 47f5887fe..8b1aa65cb 100644
--- a/lib/apytypes/_array_functions.pyi
+++ b/lib/apytypes/_array_functions.pyi
@@ -1,6 +1,6 @@
-from typing import overload, Sequence
-from apytypes import APyFloatArray, APyFixedArray
-from apytypes.typing import APyArray
+from typing import overload, Sequence, Union
+from apytypes import APyFloatArray, APyFixedArray, APyFixed, APyFloat
+from apytypes.typing import APyArray, APyScalar
 
 @overload
 def squeeze(
@@ -61,3 +61,220 @@ def expand_dims(a: APyFloatArray, axis: int | Sequence[int]) -> APyFloatArray: .
 @overload
 def expand_dims(a: APyFixedArray, axis: int | Sequence[int]) -> APyFixedArray: ...
 def expand_dims(a: APyArray, axis: int | Sequence[int]) -> APyArray: ...
+@overload
+def identity(
+    n: int,
+    int_bits: int = None,
+    frac_bits: int = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyFixedArray: ...
+def identity(
+    n: int,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: int = None,
+    mantissa_bits: int = None,
+) -> APyFloatArray: ...
+@overload
+def eye(
+    n: int,
+    m: int = None,
+    int_bits: int = None,
+    frac_bits: int = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyFixedArray: ...
+def eye(
+    n: int,
+    m: int = None,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: int = None,
+    mantissa_bits: int = None,
+) -> APyFloatArray: ...
+@overload
+def zeros(
+    shape: tuple[int, ...],
+    int_bits: int = None,
+    frac_bits: int = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyFixedArray: ...
+def zeros(
+    shape: tuple[int, ...],
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: int = None,
+    mantissa_bits: int = None,
+) -> APyFloatArray: ...
+@overload
+def ones(
+    shape: tuple[int, ...],
+    int_bits: int = None,
+    frac_bits: int = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyFixedArray: ...
+def ones(
+    shape: tuple[int, ...],
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: int = None,
+    mantissa_bits: int = None,
+) -> APyFloatArray: ...
+@overload
+def full(
+    shape: tuple[int, ...],
+    fill_value: APyScalar = None,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyArray: ...
+@overload
+def full(
+    shape: tuple[int, ...],
+    fill_value: APyFixed = None,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyFixedArray: ...
+@overload
+def full(
+    shape: tuple[int, ...],
+    fill_value: APyFloat = None,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyFloatArray: ...
+@overload
+def full(
+    shape: tuple[int, ...],
+    fill_value: Union[int, float] = None,
+    int_bits: int = None,
+    frac_bits: int = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyFixedArray: ...
+@overload
+def full(
+    shape: tuple[int, ...],
+    fill_value: Union[int, float] = None,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: int = None,
+    mantissa_bits: int = None,
+) -> APyFloatArray: ...
+@overload
+def zeros_like(
+    a: APyArray,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyArray: ...
+@overload
+def zeros_like(
+    a: APyFixedArray,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyFixedArray: ...
+def zeros_like(
+    a: APyFloatArray,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyFloatArray: ...
+@overload
+def ones_like(
+    a: APyArray,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyArray: ...
+@overload
+def ones_like(
+    a: APyFixedArray,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyFixedArray: ...
+def ones_like(
+    a: APyFloatArray,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyFloatArray: ...
+@overload
+def full_like(
+    a: APyArray,
+    fill_value: APyScalar = None,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyArray: ...
+@overload
+def full_like(
+    a: APyArray,
+    fill_value: APyFixed = None,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyFixedArray: ...
+@overload
+def full_like(
+    a: APyArray,
+    fill_value: APyFloat = None,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyFloatArray: ...
+@overload
+def full_like(
+    a: APyArray,
+    fill_value: Union[int, float] = None,
+    int_bits: int = None,
+    frac_bits: int = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyFixedArray: ...
+def full_like(
+    a: APyArray,
+    fill_value: Union[int, float] = None,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: int = None,
+    mantissa_bits: int = None,
+) -> APyFloatArray: ...
+@overload
+def arange(
+    start: Union[int, float],
+    stop: Union[int, float] = None,
+    step: Union[int, float, None] = None,
+    int_bits: int = None,
+    frac_bits: int = None,
+    exp_bits: None = None,
+    mantissa_bits: None = None,
+) -> APyFixedArray: ...
+@overload
+def arange(
+    start: Union[int, float],
+    stop: Union[int, float] = None,
+    step: Union[int, float, None] = None,
+    int_bits: None = None,
+    frac_bits: None = None,
+    exp_bits: int = None,
+    mantissa_bits: int = None,
+) -> APyFloatArray: ...
diff --git a/lib/apytypes/_typing.py b/lib/apytypes/_typing.py
index 139444794..82d42d8cd 100644
--- a/lib/apytypes/_typing.py
+++ b/lib/apytypes/_typing.py
@@ -1,6 +1,7 @@
 from typing import Union
-from apytypes import APyFloatArray, APyFixedArray
+from apytypes import APyFloatArray, APyFixedArray, APyFloat, APyFixed
 # The following are type aliases. Once Python 3.9 is dropped, they should be annotated
 # using ``typing.TypeAlias`` and Unions should be converted to using ``|`` syntax.
 
 APyArray = Union[APyFixedArray, APyFloatArray]
+APyScalar = Union[APyFixed, APyFloat]
diff --git a/src/apyfixed_util.h b/src/apyfixed_util.h
index 0cab319dd..6ca2e9e29 100644
--- a/src/apyfixed_util.h
+++ b/src/apyfixed_util.h
@@ -20,6 +20,22 @@ mp_limb_t get_data_from_double(double value, int bits, int frac_bits, int shift_
 //! Fast integer power by squaring.
 APyFixed ipow(APyFixed base, unsigned int n);
 
+//! Get bit pattern for the value one
+static APY_INLINE APyFixed one(int bits, int int_bits)
+{
+    std::size_t bit_index = bits - int_bits;
+    std::size_t limb_bits = sizeof(mp_limb_t) * 8;
+    std::size_t limb_index = bit_index / limb_bits;
+    std::size_t bit_offset = bit_index % limb_bits;
+
+    std::size_t num_limbs = limb_index + 1;
+    std::vector<mp_limb_t> data(num_limbs, static_cast<mp_limb_t>(0));
+
+    // Set the specified bit to 1
+    data[limb_index] |= static_cast<mp_limb_t>(1) << bit_offset;
+
+    return APyFixed(bits, int_bits, data);
+}
 /* ********************************************************************************** *
  * *    Fixed-point iterator based in-place quantization with multi-limb support    * *
  * ********************************************************************************** */
diff --git a/src/apyfixedarray.cc b/src/apyfixedarray.cc
index 9d2cf8d1e..63435c283 100644
--- a/src/apyfixedarray.cc
+++ b/src/apyfixedarray.cc
@@ -1,5 +1,4 @@
-// Python object access through Nanobind
-#include <nanobind/nanobind.h>
+#include "nanobind/nanobind.h"
 #include <nanobind/ndarray.h>
 #include <nanobind/stl/variant.h> // std::variant (with nanobind support)
 #include <variant>
@@ -1762,6 +1761,82 @@ APyFixedArray APyFixedArray::from_array(
     return result;
 }
 
+/* ****************************************************************************** *
+ *                    Static methods for creating arrays                          *
+ * ****************************************************************************** */
+
+APyFixedArray APyFixedArray::zeros(
+    const nb::tuple& shape,
+    std::optional<int> int_bits,
+    std::optional<int> frac_bits,
+    std::optional<int> bits
+)
+{
+    std::vector<std::size_t> new_shape = ::shape_from_tuple(shape);
+    APyFixedArray result(new_shape, int_bits, frac_bits, bits);
+    return result;
+}
+
+APyFixedArray APyFixedArray::ones(
+    const nb::tuple& shape,
+    std::optional<int> int_bits,
+    std::optional<int> frac_bits,
+    std::optional<int> bits
+)
+{
+    const int final_int_bits
+        = int_bits ? int_bits.value() : bits.value() - frac_bits.value();
+    const int final_bits = bits ? bits.value() : int_bits.value() + frac_bits.value();
+
+    APyFixed fixed_one = one(final_bits, final_int_bits);
+    return full(shape, fixed_one);
+}
+
+APyFixedArray APyFixedArray::eye(
+    const nb::int_& N,
+    std::optional<nb::int_> M,
+    std::optional<int> int_bits,
+    std::optional<int> frac_bits,
+    std::optional<int> bits
+)
+{
+    const int final_int_bits
+        = int_bits ? int_bits.value() : bits.value() - frac_bits.value();
+    const int final_bits = bits ? bits.value() : int_bits.value() + frac_bits.value();
+    APyFixed fixed_one = one(final_bits, final_int_bits);
+
+    // Use N for both dimensions if M is not provided
+    nb::tuple shape = M ? nb::make_tuple(N, M.value()) : nb::make_tuple(N, N);
+
+    return diagonal(shape, fixed_one);
+}
+
+APyFixedArray APyFixedArray::identity(
+    const nb::int_& N,
+    std::optional<int> int_bits,
+    std::optional<int> frac_bits,
+    std::optional<int> bits
+)
+{
+    return eye(N, std::nullopt, int_bits, frac_bits, bits);
+}
+
+APyFixedArray APyFixedArray::full(const nb::tuple& shape, const APyFixed& fill_value)
+{
+    std::vector<std::size_t> new_shape = ::shape_from_tuple(shape);
+    APyFixedArray result(new_shape, fill_value.bits(), fill_value.int_bits());
+
+    std::size_t num_elem = fold_shape(new_shape);
+    std::size_t itemsize = result._itemsize;
+
+    for (std::size_t index = 0; index < num_elem; ++index) {
+        std::copy_n(
+            fill_value._data.begin(), itemsize, result._data.begin() + index * itemsize
+        );
+    }
+    return result;
+}
+
 /* ********************************************************************************** *
  * *                            Private member functions                            * *
  * ********************************************************************************** */
@@ -2326,3 +2401,27 @@ void APyFixedArray::_set_values_from_ndarray(const nb::ndarray<nb::c_contig>& nd
         "unsupported `dtype` expecting integer/float"
     );
 }
+
+APyFixedArray
+APyFixedArray::diagonal(const nb::tuple& shape, const APyFixed& fill_value)
+{
+    std::vector<std::size_t> new_shape = ::shape_from_tuple(shape);
+    if (new_shape.size() > 2) {
+        throw nb::value_error(
+            "Creating higher dimensional diagonal arrays are not yet defined"
+        );
+    }
+    APyFixedArray result(new_shape, fill_value.bits(), fill_value.int_bits());
+    std::size_t itemsize = result._itemsize;
+
+    std::size_t min_dim = *std::min_element(new_shape.begin(), new_shape.end());
+    std::vector<std::size_t> strides = ::strides_from_shape(new_shape);
+
+    for (std::size_t i = 0; i < min_dim; ++i) {
+        std::size_t index = i * std::accumulate(strides.begin(), strides.end(), 0);
+        std::copy_n(
+            fill_value._data.begin(), itemsize, result._data.begin() + index * itemsize
+        );
+    }
+    return result;
+}
diff --git a/src/apyfixedarray.h b/src/apyfixedarray.h
index 36725bec3..7f2b0b1f9 100644
--- a/src/apyfixedarray.h
+++ b/src/apyfixedarray.h
@@ -318,6 +318,47 @@ class APyFixedArray : public APyBuffer<mp_limb_t> {
         std::optional<int> bits = std::nullopt
     );
 
+    /* ****************************************************************************** *
+     *                     Static methods for array initialization                    *
+     * ****************************************************************************** */
+
+    //! Create an `APyFixedArray` initialized with zeros
+    static APyFixedArray zeros(
+        const nb::tuple& shape,
+        std::optional<int> int_bits = std::nullopt,
+        std::optional<int> frac_bits = std::nullopt,
+        std::optional<int> bits = std::nullopt
+    );
+    //! Create an `APyFixedArray` initialized with ones
+    static APyFixedArray ones(
+        const nb::tuple& shape,
+        std::optional<int> int_bits = std::nullopt,
+        std::optional<int> frac_bits = std::nullopt,
+        std::optional<int> bits = std::nullopt
+    );
+
+    //! Create an `APyFixedArray` with ones on the diagonal and zeros elsewhere
+    static APyFixedArray
+    eye(const nb::int_& N,
+        std::optional<nb::int_> M = std::nullopt,
+        std::optional<int> int_bits = std::nullopt,
+        std::optional<int> frac_bits = std::nullopt,
+        std::optional<int> bits = std::nullopt);
+
+    //! Create a square `APyFixedArray` with ones on the diagonal and zeros elsewhere
+    static APyFixedArray identity(
+        const nb::int_& N,
+        std::optional<int> int_bits = std::nullopt,
+        std::optional<int> frac_bits = std::nullopt,
+        std::optional<int> bits = std::nullopt
+    );
+
+    //! Create an `APyFixedArray` initialized with a specified fill value
+    static APyFixedArray full(const nb::tuple& shape, const APyFixed& fill_value);
+
+    //! Create an `APyFixedArray` with a specified diagonal value
+    static APyFixedArray diagonal(const nb::tuple& shape, const APyFixed& fill_value);
+
 private:
     /* ****************************************************************************** *
      * *                          Private member functions                          * *
diff --git a/src/apyfixedarray_wrapper.cc b/src/apyfixedarray_wrapper.cc
index 6cb778c2b..0a3a09c74 100644
--- a/src/apyfixedarray_wrapper.cc
+++ b/src/apyfixedarray_wrapper.cc
@@ -1083,6 +1083,138 @@ void bind_fixed_array(nb::module_& m)
             :class:`APyFixedArray`
             )pbdoc"
         )
+        .def_static(
+            "zeros",
+            &APyFixedArray::zeros,
+            nb::arg("shape"),
+            nb::arg("int_bits") = nb::none(),
+            nb::arg("frac_bits") = nb::none(),
+            nb::arg("bits") = nb::none(),
+            R"pbdoc(
+        Initializes an array with zeros.
+
+        Parameters
+        ----------
+        shape : tuple
+            Shape of the array.
+        int_bits : int, optional
+            Number of integer bits. Default is None.
+        frac_bits : int, optional
+            Number of fractional bits. Default is None.
+        bits : int, optional
+            Total number of bits. Default is None.
+
+        Returns
+        -------
+        APyFixedArray
+            An array initialized with zeros.
+    )pbdoc"
+        )
+        .def_static(
+            "ones",
+            &APyFixedArray::ones,
+            nb::arg("shape"),
+            nb::arg("int_bits") = nb::none(),
+            nb::arg("frac_bits") = nb::none(),
+            nb::arg("bits") = nb::none(),
+            R"pbdoc(
+        Initializes an array with ones.
+
+        Parameters
+        ----------
+        shape : tuple
+            Shape of the array.
+        int_bits : int, optional
+            Number of integer bits. Default is None.
+        frac_bits : int, optional
+            Number of fractional bits. Default is None.
+        bits : int, optional
+            Total number of bits. Default is None.
+
+        Returns
+        -------
+        APyFixedArray
+            An array initialized with ones.
+    )pbdoc"
+        )
+        .def_static(
+            "eye",
+            &APyFixedArray::eye,
+            nb::arg("n"),
+            nb::arg("m") = nb::none(),
+            nb::arg("int_bits") = nb::none(),
+            nb::arg("frac_bits") = nb::none(),
+            nb::arg("bits") = nb::none(),
+            R"pbdoc(
+        Initializes an array with ones on the diagonal.
+
+        Parameters
+        ----------
+        n : int
+            Number of rows.
+        m : int, optional
+            Number of columns. Default is None.
+        int_bits : int, optional
+            Number of integer bits. Default is None.
+        frac_bits : int, optional
+            Number of fractional bits. Default is None.
+        bits : int, optional
+            Total number of bits. Default is None.
+
+        Returns
+        -------
+        APyFixedArray
+            An array with the specified value on the diagonal.
+    )pbdoc"
+        )
+        .def_static(
+            "identity",
+            &APyFixedArray::identity,
+            nb::arg("n"),
+            nb::arg("int_bits") = nb::none(),
+            nb::arg("frac_bits") = nb::none(),
+            nb::arg("bits") = nb::none(),
+            R"pbdoc(
+        Initializes an identity matrix with ones on the diagonal.
+
+        Parameters
+        ----------
+        n : int
+            Number of rows (and columns) in n x n output.
+        int_bits : int, optional
+            Number of integer bits. Default is None.
+        frac_bits : int, optional
+            Number of fractional bits. Default is None.
+        bits : int, optional
+            Total number of bits. Default is None.
+
+        Returns
+        -------
+        APyFixedArray
+            An identity matrix with ones on the diagonal.
+    )pbdoc"
+        )
+        .def_static(
+            "full",
+            &APyFixedArray::full,
+            nb::arg("shape"),
+            nb::arg("fill_value"),
+            R"pbdoc(
+        Initializes an array with the specified value.
+
+        Parameters
+        ----------
+        shape : tuple
+            Shape of the array.
+        fill_value : APyFixed
+            Value to fill the array.
+
+        Returns
+        -------
+        APyFixedArray
+            An array filled with the specified value.
+    )pbdoc"
+        )
 
         /*
          * Dunder methods
diff --git a/src/apyfloat_util.h b/src/apyfloat_util.h
index 549d47ed7..c1f42b074 100644
--- a/src/apyfloat_util.h
+++ b/src/apyfloat_util.h
@@ -5,6 +5,8 @@
 #include "apyfixed_util.h"
 #include "apyfloat.h"
 #include "apytypes_common.h"
+#include <optional>
+#include <vector>
 
 /*!
  * Sizes of APyFloat datatypes
@@ -53,6 +55,18 @@ exp_t calc_bias_general(
     int new_exp_bits, int exp_bits1, exp_t bias1, int exp_bits2, exp_t bias2
 );
 
+//! Get bit pattern for the value one
+static APY_INLINE APyFloat
+one(std::uint8_t exp_bits,
+    std::uint8_t man_bits,
+    std::optional<exp_t> bias = std::nullopt)
+{
+    APyFloat result(exp_bits, man_bits, bias.value_or(APyFloat::ieee_bias(exp_bits)));
+    APyFloatData data = { 0, result.get_bias(), 0 };
+    result.set_data(data);
+    return result;
+}
+
 //! Quantize mantissa
 APY_INLINE void quantize_mantissa(
     man_t& man,
diff --git a/src/apyfloatarray.cc b/src/apyfloatarray.cc
index 3a3b86b4b..e1eafbf54 100644
--- a/src/apyfloatarray.cc
+++ b/src/apyfloatarray.cc
@@ -995,6 +995,102 @@ std::variant<APyFloatArray, APyFloat> APyFloatArray::matmul(const APyFloatArray&
     ));
 }
 
+/* ****************************************************************************** *
+ *                     Static methods for array initialization                    *
+ * ****************************************************************************** */
+
+APyFloatArray APyFloatArray::zeros(
+    const nb::tuple& shape,
+    std::uint8_t exp_bits,
+    std::uint8_t man_bits,
+    std::optional<exp_t> bias
+)
+{
+    std::vector<std::size_t> new_shape = ::shape_from_tuple(shape);
+    APyFloatArray result(new_shape, exp_bits, man_bits, bias);
+    return result;
+}
+
+APyFloatArray APyFloatArray::ones(
+    const nb::tuple& shape,
+    std::uint8_t exp_bits,
+    std::uint8_t man_bits,
+    std::optional<exp_t> bias
+)
+{
+    APyFloat float_one = one(exp_bits, man_bits, bias);
+    return full(shape, float_one);
+}
+
+APyFloatArray APyFloatArray::eye(
+    const nb::int_& N,
+    std::uint8_t exp_bits,
+    std::uint8_t man_bits,
+    std::optional<nb::int_> M,
+    std::optional<exp_t> bias
+)
+{
+    APyFloat float_one = one(exp_bits, man_bits, bias);
+    // Use N for both dimensions if M is not provided
+    nb::tuple shape = M ? nb::make_tuple(N, M.value()) : nb::make_tuple(N, N);
+    return diagonal(shape, float_one);
+}
+
+APyFloatArray APyFloatArray::identity(
+    const nb::int_& N,
+    std::uint8_t exp_bits,
+    std::uint8_t man_bits,
+    std::optional<exp_t> bias
+)
+{
+    return eye(N, exp_bits, man_bits, std::nullopt, bias);
+}
+
+APyFloatArray APyFloatArray::full(const nb::tuple& shape, const APyFloat& fill_value)
+{
+    std::vector<std::size_t> new_shape = ::shape_from_tuple(shape);
+    APyFloatArray result(
+        new_shape,
+        fill_value.get_exp_bits(),
+        fill_value.get_man_bits(),
+        fill_value.get_bias()
+    );
+
+    std::size_t num_elem = ::fold_shape(new_shape);
+    result.data = std::vector<APyFloatData>(num_elem, fill_value.get_data());
+    return result;
+}
+
+APyFloatArray
+APyFloatArray::diagonal(const nb::tuple& shape, const APyFloat& fill_value)
+{
+    std::vector<std::size_t> new_shape = ::shape_from_tuple(shape);
+    if (new_shape.size() > 2) {
+        throw nb::value_error(
+            "Creating higher dimensional diagonal arrays are not yet defined"
+        );
+    }
+    APyFloatArray result(
+        new_shape,
+        fill_value.get_exp_bits(),
+        fill_value.get_man_bits(),
+        fill_value.get_bias()
+    );
+
+    std::size_t min_dim = *std::min_element(new_shape.begin(), new_shape.end());
+    std::vector<std::size_t> strides = ::strides_from_shape(new_shape);
+    std::size_t multiplier = std::accumulate(strides.begin(), strides.end(), 0);
+    for (std::size_t i = 0; i < min_dim; ++i) {
+        std::size_t index = i * multiplier;
+        result.data[index] = fill_value.get_data();
+    }
+    return result;
+}
+
+/* ****************************************************************************** *
+ *                            Methods for array modification                      *
+ * ****************************************************************************** */
+
 APyFloatArray
 APyFloatArray::squeeze(std::optional<std::variant<nb::int_, nb::tuple>> axis) const
 {
diff --git a/src/apyfloatarray.h b/src/apyfloatarray.h
index 3cad0addb..f596243c8 100644
--- a/src/apyfloatarray.h
+++ b/src/apyfloatarray.h
@@ -95,6 +95,43 @@ class APyFloatArray {
     //! Set data fields based on an and-array of doubles
     void _set_values_from_ndarray(const nanobind::ndarray<nanobind::c_contig>& ndarray);
 
+    /* ****************************************************************************** *
+     *                     Static methods for array initialization                    *
+     * ****************************************************************************** */
+
+    static APyFloatArray zeros(
+        const nb::tuple& shape,
+        std::uint8_t exp_bits,
+        std::uint8_t man_bits,
+        std::optional<exp_t> bias = std::nullopt
+    );
+
+    static APyFloatArray ones(
+        const nb::tuple& shape,
+        std::uint8_t exp_bits,
+        std::uint8_t man_bits,
+        std::optional<exp_t> bias = std::nullopt
+    );
+
+    static APyFloatArray
+    eye(const nb::int_& N,
+        std::uint8_t exp_bits,
+        std::uint8_t man_bits,
+        std::optional<nb::int_> M,
+        std::optional<exp_t> bias = std::nullopt);
+
+    static APyFloatArray identity(
+        const nb::int_& N,
+        std::uint8_t exp_bits,
+        std::uint8_t man_bits,
+        std::optional<exp_t> bias = std::nullopt
+    );
+
+    //! Create an `APyFloatArray` initialized with a specified fill value
+    static APyFloatArray full(const nb::tuple& shape, const APyFloat& fill_value);
+
+    static APyFloatArray diagonal(const nb::tuple& shape, const APyFloat& fill_value);
+
     /* ****************************************************************************** *
      * *                          Public member functions                           * *
      * ****************************************************************************** */
diff --git a/src/apyfloatarray_wrapper.cc b/src/apyfloatarray_wrapper.cc
index cd0634d6f..77951bf9d 100644
--- a/src/apyfloatarray_wrapper.cc
+++ b/src/apyfloatarray_wrapper.cc
@@ -505,7 +505,139 @@ void bind_float_array(nb::module_& m)
             :class:`APyFloatArray`
             )pbdoc"
         )
+        .def_static(
+            "zeros",
+            &APyFloatArray::zeros,
+            nb::arg("shape"),
+            nb::arg("exp_bits"),
+            nb::arg("man_bits"),
+            nb::arg("bias") = std::nullopt,
+            R"pbdoc(
+
+            Initializes an array with zeros.
+
+            Parameters
+            ----------
+            shape : tuple
+                Shape of the array.
+            exp_bits : int
+                Number of exponent bits.
+            man_bits : int
+                Number of mantissa bits.
+            bias : optional
+                Set bias. Default is None.
+
+            Returns
+            -------
+            APyFloatArray
+                An array filled with zeros.
+    )pbdoc"
+        )
+        .def_static(
+            "ones",
+            &APyFloatArray::ones,
+            nb::arg("shape"),
+            nb::arg("exp_bits"),
+            nb::arg("man_bits"),
+            nb::arg("bias") = std::nullopt,
+            R"pbdoc(
+        Initializes an array with ones.
+
+        Parameters
+        ----------
+        shape : tuple
+            Shape of the array.
+        exp_bits : int
+            Number of exponent bits.
+        man_bits : int
+            Number of mantissa bits.
+        bias : optional
+            Set bias. Default is None.
+
+        Returns
+        -------
+        APyFloatArray
+            An array filled with ones.
+    )pbdoc"
+        )
+        .def_static(
+            "eye",
+            &APyFloatArray::eye,
+            nb::arg("n"),
+            nb::arg("exp_bits"),
+            nb::arg("man_bits"),
+            nb::arg("m") = nb::none(),
+            nb::arg("bias") = std::nullopt,
+            R"pbdoc(
+        Initializes an array with ones on the diagonal.
 
+        Parameters
+        ----------
+        n : int
+            Number of rows (and columns) in the n x n output.
+        exp_bits : int
+            Number of exponent bits.
+        man_bits : int
+            Number of mantissa bits.
+        m : int, optional
+            Number of columns. Default is None, which results in an n x n output.
+        bias : optional
+            Set bias. Default is None.
+
+        Returns
+        -------
+        APyFloatArray
+            An array with the specified value on the diagonal.
+    )pbdoc"
+        )
+        .def_static(
+            "identity",
+            &APyFloatArray::identity,
+            nb::arg("n"),
+            nb::arg("exp_bits"),
+            nb::arg("man_bits"),
+            nb::arg("bias") = std::nullopt,
+            R"pbdoc(
+        Initializes an identity matrix with ones on the diagonal.
+
+        Parameters
+        ----------
+        n : int
+            Number of rows (and columns) in the n x n output.
+        exp_bits : int
+            Number of exponent bits.
+        man_bits : int
+            Number of mantissa bits.
+        bias : optional
+            Set bias. Default is None.
+
+        Returns
+        -------
+        APyFloatArray
+            An identity matrix with ones on the diagonal.
+    )pbdoc"
+        )
+        .def_static(
+            "full",
+            &APyFloatArray::full,
+            nb::arg("shape"),
+            nb::arg("fill_value"),
+            R"pbdoc(
+        Initializes an array filled with the specified value.
+
+        Parameters
+        ----------
+        shape : tuple
+            Shape of the array.
+        fill_value : APyFloat
+            Value to fill the array.
+
+        Returns
+        -------
+        APyFloatArray
+            An array filled with the specified value.
+    )pbdoc"
+        )
         /*
          * Dunder methods
          */
diff --git a/src/array_utils.h b/src/array_utils.h
index 93f3164a4..ed51c2f39 100644
--- a/src/array_utils.h
+++ b/src/array_utils.h
@@ -63,6 +63,20 @@ shape_from_tuple(nb::tuple new_shape, size_t elem_count)
     return new_shape_vec;
 }
 
+static APY_INLINE std::vector<std::size_t> shape_from_tuple(nb::tuple new_shape)
+{
+    std::vector<std::size_t> new_shape_vec;
+    for (auto it = new_shape.begin(); it != new_shape.end(); ++it) {
+        int current_value = nb::cast<int>(*it);
+        if (current_value < 0) {
+            throw nb::value_error("Negative dimensions or are not allowed.");
+        } else {
+            new_shape_vec.push_back(static_cast<std::size_t>(current_value));
+        }
+    }
+    return new_shape_vec;
+}
+
 /**
  * @brief Converts the provided axes into a vector of positive indices.
  *


Test patch:
diff --git a/lib/test/test_array_functions.py b/lib/test/test_array_functions.py
index 85d305677..c9681a5db 100644
--- a/lib/test/test_array_functions.py
+++ b/lib/test/test_array_functions.py
@@ -7,6 +7,17 @@
     moveaxis,
     expand_dims,
     swapaxes,
+    eye,
+    identity,
+    zeros,
+    zeros_like,
+    ones,
+    ones_like,
+    full,
+    full_like,
+    arange,
+    APyFloat,
+    APyFixed,
     APyFixedArray,
     APyFloatArray,
 )
@@ -132,3 +143,432 @@ def test_expanddims():
     b = _floatArrange(10).reshape((1, 2, 5))
     c = expand_dims(a, 0)
     assert c.is_identical(b), "expanddims failed for float array"
+
+
+@pytest.mark.parametrize(
+    "n, m, nums",
+    [
+        (1, 1, [1]),
+        (1, 2, [1, 0]),
+        (2, None, [1, 0, 0, 1]),
+        (2, 3, [1, 0, 0, 0, 1, 0]),
+        (4, None, [1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1]),
+    ],
+)
+def test_eye(n, m, nums):
+    def check_eye(int_bits=None, frac_bits=None, exp_bits=None, mantissa_bits=None):
+        a = eye(
+            n=n,
+            m=m,
+            int_bits=int_bits,
+            frac_bits=frac_bits,
+            exp_bits=exp_bits,
+            mantissa_bits=mantissa_bits,
+        )
+        if isinstance(a, APyFixedArray):
+            b = APyFixedArray.from_float(
+                nums, int_bits=int_bits, frac_bits=frac_bits
+            ).reshape((n, m if m is not None else n))
+        elif isinstance(a, APyFloatArray):
+            b = APyFloatArray.from_float(
+                nums, exp_bits=exp_bits, man_bits=mantissa_bits
+            ).reshape((n, m if m is not None else n))
+        assert a.is_identical(
+            b
+        ), f"eye on {type(a).__name__} didn't work when n={n}, m={m}. Expected result was {nums} but got \n {a}"
+
+    # Test cases for APyFixedArray
+    check_eye(int_bits=5, frac_bits=5)
+    check_eye(int_bits=12314, frac_bits=1832)
+
+    # Test cases for APyFloatArray
+    check_eye(exp_bits=13, mantissa_bits=28)
+    check_eye(exp_bits=16, mantissa_bits=5)
+
+
+@pytest.mark.parametrize(
+    "n, nums",
+    [
+        (1, [1]),
+        (2, [1, 0, 0, 1]),
+        (3, [1, 0, 0, 0, 1, 0, 0, 0, 1]),
+        (4, [1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1]),
+        (
+            5,
+            [1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1],
+        ),
+    ],
+)
+def test_identity(n, nums):
+    def check_identity(
+        int_bits=None, frac_bits=None, exp_bits=None, mantissa_bits=None
+    ):
+        a = identity(
+            n=n,
+            int_bits=int_bits,
+            frac_bits=frac_bits,
+            exp_bits=exp_bits,
+            mantissa_bits=mantissa_bits,
+        )
+        if isinstance(a, APyFixedArray):
+            b = APyFixedArray.from_float(
+                nums, int_bits=int_bits, frac_bits=frac_bits
+            ).reshape((n, n))
+        elif isinstance(a, APyFloatArray):
+            b = APyFloatArray.from_float(
+                nums, exp_bits=exp_bits, man_bits=mantissa_bits
+            ).reshape((n, n))
+        assert a.is_identical(
+            b
+        ), f"identity on {type(a).__name__} didn't work when n={n}. Expected result was {nums} but got \n {a}"
+
+    # Test cases for APyFixedArray
+    check_identity(int_bits=5, frac_bits=5)
+    check_identity(int_bits=38231, frac_bits=1237)
+
+    # Test cases for APyFloatArray
+    check_identity(exp_bits=13, mantissa_bits=28)
+    check_identity(exp_bits=16, mantissa_bits=5)
+
+
+@pytest.mark.parametrize("shape", [(i, j) for i in range(1, 5) for j in range(1, 5)])
+def test_zeros(shape):
+    def check_zeros(int_bits=None, frac_bits=None, exp_bits=None, mantissa_bits=None):
+        a = zeros(
+            shape,
+            int_bits=int_bits,
+            frac_bits=frac_bits,
+            exp_bits=exp_bits,
+            mantissa_bits=mantissa_bits,
+        )
+        if isinstance(a, APyFixedArray):
+            b = APyFixedArray.from_float(
+                [0] * (shape[0] * shape[1]), int_bits=int_bits, frac_bits=frac_bits
+            ).reshape(shape)
+        if isinstance(a, APyFloatArray):
+            b = APyFloatArray.from_float(
+                [0] * (shape[0] * shape[1]), exp_bits=exp_bits, man_bits=mantissa_bits
+            ).reshape(shape)
+        assert a.is_identical(
+            b
+        ), f"zeros on {a.__name__} didn't work when shape={shape}."
+
+    # Test cases for APyFixedArray
+    check_zeros(int_bits=5, frac_bits=5)
+
+    check_zeros(int_bits=12314, frac_bits=1832)
+
+    check_zeros(int_bits=0, frac_bits=5)
+    # Test cases for APyFloatArray
+    check_zeros(exp_bits=13, mantissa_bits=28)
+
+    check_zeros(exp_bits=13, mantissa_bits=0)
+    check_zeros(exp_bits=16, mantissa_bits=5)
+
+
+@pytest.mark.parametrize("shape", [(i, j) for i in range(1, 5) for j in range(1, 5)])
+def test_zeros_like(shape):
+    def check_zeros_like(
+        ArrayType, int_bits=None, frac_bits=None, exp_bits=None, mantissa_bits=None
+    ):
+        if ArrayType is APyFixedArray:
+            b = ArrayType.from_float(
+                [0] * (shape[0] * shape[1]), int_bits=int_bits, frac_bits=frac_bits
+            ).reshape(shape)
+            a = zeros_like(b, int_bits=int_bits, frac_bits=frac_bits)
+        elif ArrayType is APyFloatArray:
+            b = ArrayType.from_float(
+                [0] * (shape[0] * shape[1]), exp_bits=exp_bits, man_bits=mantissa_bits
+            ).reshape(shape)
+            a = zeros_like(b, exp_bits=exp_bits, mantissa_bits=mantissa_bits)
+        assert a.is_identical(
+            b
+        ), f"zeros_like on {ArrayType.__name__} didn't work when shape={shape}."
+
+    # Test cases for APyFixedArray
+    check_zeros_like(APyFixedArray, int_bits=5, frac_bits=5)
+    check_zeros_like(APyFixedArray, int_bits=12314, frac_bits=1832)
+
+    # Test cases for APyFloatArray
+    check_zeros_like(APyFloatArray, exp_bits=13, mantissa_bits=28)
+    check_zeros_like(APyFloatArray, exp_bits=16, mantissa_bits=5)
+
+
+@pytest.mark.parametrize("shape", [(i, j) for i in range(1, 5) for j in range(1, 5)])
+def test_ones(shape):
+    def check_ones(int_bits=None, frac_bits=None, exp_bits=None, mantissa_bits=None):
+        a = ones(
+            shape,
+            int_bits=int_bits,
+            frac_bits=frac_bits,
+            exp_bits=exp_bits,
+            mantissa_bits=mantissa_bits,
+        )
+        if isinstance(a, APyFixedArray):
+            b = APyFixedArray.from_float(
+                [1] * (shape[0] * shape[1]), int_bits=int_bits, frac_bits=frac_bits
+            ).reshape(shape)
+        elif isinstance(a, APyFloatArray):
+            b = APyFloatArray.from_float(
+                [1] * (shape[0] * shape[1]), exp_bits=exp_bits, man_bits=mantissa_bits
+            ).reshape(shape)
+        assert a.is_identical(
+            b
+        ), f"ones on {type(a).__name__} didn't work when shape={shape}."
+
+    # Test cases for APyFixedArray
+    check_ones(int_bits=5, frac_bits=5)
+    check_ones(int_bits=12314, frac_bits=1832)
+
+    # Test cases for APyFloatArray
+    check_ones(exp_bits=13, mantissa_bits=28)
+    check_ones(exp_bits=16, mantissa_bits=5)
+
+
+@pytest.mark.parametrize("shape", [(i, j) for i in range(1, 5) for j in range(1, 5)])
+def test_ones_like(shape):
+    def check_ones_like(
+        ArrayType, int_bits=None, frac_bits=None, exp_bits=None, mantissa_bits=None
+    ):
+        if ArrayType is APyFixedArray:
+            b = ArrayType.from_float(
+                [1] * (shape[0] * shape[1]), int_bits=int_bits, frac_bits=frac_bits
+            ).reshape(shape)
+            a = ones_like(b, int_bits=int_bits, frac_bits=frac_bits)
+        elif ArrayType is APyFloatArray:
+            b = ArrayType.from_float(
+                [1] * (shape[0] * shape[1]), exp_bits=exp_bits, man_bits=mantissa_bits
+            ).reshape(shape)
+            a = ones_like(b, exp_bits=exp_bits, mantissa_bits=mantissa_bits)
+        assert a.is_identical(
+            b
+        ), f"ones_like on {ArrayType.__name__} didn't work when shape={shape}."
+
+    # Test cases for APyFixedArray
+    check_ones_like(APyFixedArray, int_bits=5, frac_bits=5)
+    check_ones_like(APyFixedArray, int_bits=12314, frac_bits=1832)
+
+    # Test cases for APyFloatArray
+    check_ones_like(APyFloatArray, exp_bits=13, mantissa_bits=28)
+    check_ones_like(APyFloatArray, exp_bits=16, mantissa_bits=5)
+
+
+@pytest.mark.parametrize("shape", [(i, j) for i in range(1, 5) for j in range(1, 5)])
+def test_full(shape):
+    def check_full(
+        ArrayType,
+        fill_value,
+        int_bits=None,
+        frac_bits=None,
+        exp_bits=None,
+        mantissa_bits=None,
+    ):
+        if ArrayType is APyFixedArray:
+            b = ArrayType.from_float(
+                [fill_value] * (shape[0] * shape[1]),
+                int_bits=int_bits,
+                frac_bits=frac_bits,
+            ).reshape(shape)
+            num = APyFixed.from_float(
+                fill_value, int_bits=int_bits, frac_bits=frac_bits
+            )
+            a = full(shape, num, int_bits=int_bits, frac_bits=frac_bits)
+        elif ArrayType is APyFloatArray:
+            b = ArrayType.from_float(
+                [fill_value] * (shape[0] * shape[1]),
+                exp_bits=exp_bits,
+                man_bits=mantissa_bits,
+            ).reshape(shape)
+            num = APyFloat.from_float(
+                fill_value, exp_bits=exp_bits, man_bits=mantissa_bits
+            )
+            a = full(shape, num, exp_bits=exp_bits, mantissa_bits=mantissa_bits)
+        assert a.is_identical(
+            b
+        ), f"full on {ArrayType.__name__} didn't work when shape={shape}."
+
+    # Test cases for APyFixedArray
+    check_full(APyFixedArray, fill_value=shape[0] + shape[1], int_bits=5, frac_bits=5)
+    check_full(
+        APyFixedArray, fill_value=shape[0] + shape[1], int_bits=12314, frac_bits=1832
+    )
+
+    check_full(
+        APyFixedArray, fill_value=shape[0] + shape[1], int_bits=0, frac_bits=1832
+    )
+
+    # Test cases for APyFloatArray
+    check_full(
+        APyFloatArray, fill_value=shape[0] + shape[1], exp_bits=13, mantissa_bits=28
+    )
+    check_full(
+        APyFloatArray, fill_value=shape[0] + shape[1], exp_bits=16, mantissa_bits=5
+    )
+    check_full(
+        APyFloatArray, fill_value=shape[0] + shape[1], exp_bits=0, mantissa_bits=5
+    )
+
+
+@pytest.mark.parametrize("shape", [(i, j) for i in range(1, 5) for j in range(1, 5)])
+def test_full_like(shape):
+    def check_full_like(
+        ArrayType,
+        fill_value,
+        int_bits=None,
+        frac_bits=None,
+        exp_bits=None,
+        mantissa_bits=None,
+    ):
+        if ArrayType is APyFixedArray:
+            b = ArrayType.from_float(
+                [fill_value] * (shape[0] * shape[1]),
+                int_bits=int_bits,
+                frac_bits=frac_bits,
+            ).reshape(shape)
+            num = APyFixed.from_float(
+                fill_value, int_bits=int_bits, frac_bits=frac_bits
+            )
+            a = full_like(b, num, int_bits=int_bits, frac_bits=frac_bits)
+            assert a.is_identical(
+                b
+            ), f"full_like on {ArrayType.__name__} didn't work when shape={shape}."
+            a = full_like(b, fill_value, int_bits=int_bits, frac_bits=frac_bits)
+            assert a.is_identical(
+                b
+            ), f"full_like on {ArrayType.__name__} didn't work when shape={shape}."
+        elif ArrayType is APyFloatArray:
+            b = ArrayType.from_float(
+                [fill_value] * (shape[0] * shape[1]),
+                exp_bits=exp_bits,
+                man_bits=mantissa_bits,
+            ).reshape(shape)
+            num = APyFloat.from_float(
+                fill_value, exp_bits=exp_bits, man_bits=mantissa_bits
+            )
+            a = full_like(b, num, exp_bits=exp_bits, mantissa_bits=mantissa_bits)
+            assert a.is_identical(
+                b
+            ), f"full_like on {ArrayType.__name__} didn't work when shape={shape}."
+            a = full_like(b, fill_value, exp_bits=exp_bits, mantissa_bits=mantissa_bits)
+            assert a.is_identical(
+                b
+            ), f"full_like on {ArrayType.__name__} didn't work when shape={shape}."
+
+    # Test cases for APyFixedArray
+    check_full_like(
+        APyFixedArray, fill_value=shape[0] + shape[1], int_bits=5, frac_bits=5
+    )
+    check_full_like(
+        APyFixedArray, fill_value=shape[0] + shape[1], int_bits=12314, frac_bits=1832
+    )
+    check_full_like(
+        APyFixedArray, fill_value=shape[0] + shape[1], int_bits=0, frac_bits=5
+    )
+
+    check_full_like(
+        APyFixedArray, fill_value=shape[0] + shape[1], int_bits=5, frac_bits=0
+    )
+
+    check_full_like(
+        APyFixedArray, fill_value=shape[0] + shape[1], int_bits=5, frac_bits=0
+    )
+
+    # Test cases for APyFloatArray
+
+    check_full_like(
+        APyFloatArray, fill_value=shape[0] + shape[1], exp_bits=13, mantissa_bits=28
+    )
+    check_full_like(
+        APyFloatArray, fill_value=shape[0] + shape[1], exp_bits=16, mantissa_bits=5
+    )
+    check_full_like(
+        APyFloatArray, fill_value=shape[0] + shape[1], exp_bits=16, mantissa_bits=0
+    )
+
+
+def test_arrange():
+    a = arange(10, int_bits=10, frac_bits=5)
+    b = APyFixedArray.from_float([i for i in range(10)], int_bits=10, frac_bits=5)
+    assert a.is_identical(b), "Arrange for APyFixedArray failed"
+
+    a = arange(10, exp_bits=10, mantissa_bits=5)
+    b = APyFloatArray.from_float([i for i in range(10)], exp_bits=10, man_bits=5)
+    assert a.is_identical(b), "Arrange for APyFloatArray failed"
+
+    a = arange(5, 8, int_bits=10, frac_bits=5)
+    assert a.is_identical(APyFixedArray.from_float([5, 6, 7], int_bits=10, frac_bits=5))
+
+    a = arange(8, 5, -1, int_bits=10, frac_bits=5)
+    assert a.is_identical(APyFixedArray.from_float([8, 7, 6], int_bits=10, frac_bits=5))
+
+    a = arange(8, 5, -2, int_bits=10, frac_bits=5)
+    assert a.is_identical(APyFixedArray.from_float([8, 6], int_bits=10, frac_bits=5))
+
+    a = arange(8, 5, -3, int_bits=10, frac_bits=5)
+    assert a.is_identical(APyFixedArray.from_float([8], int_bits=10, frac_bits=5))
+
+    a = arange(0, 2, 2 ** (-2), int_bits=10, frac_bits=2)
+    assert a.is_identical(
+        APyFixedArray.from_float(
+            [0.0, 0.25, 0.5, 0.75, 1.0, 1.25, 1.5, 1.75], int_bits=10, frac_bits=2
+        )
+    )
+
+    a = arange(0, 2, 2 ** (-2), int_bits=10, frac_bits=2)
+    assert a.is_identical(
+        APyFixedArray.from_float(
+            [0.0, 0.25, 0.5, 0.75, 1.0, 1.25, 1.5, 1.75], int_bits=10, frac_bits=2
+        )
+    )
+
+    a = arange(0, 1, 2 ** (-4), int_bits=10, frac_bits=4)
+    assert a.is_identical(
+        APyFixedArray.from_float(
+            [
+                0,
+                0.0625,
+                0.125,
+                0.1875,
+                0.25,
+                0.3125,
+                0.375,
+                0.4375,
+                0.5,
+                0.5625,
+                0.625,
+                0.6875,
+                0.75,
+                0.8125,
+                0.875,
+                0.9375,
+            ],
+            int_bits=10,
+            frac_bits=4,
+        )
+    )
+
+    a = arange(1, 0, -(2 ** (-4)), int_bits=10, frac_bits=4)
+    assert a.is_identical(
+        APyFixedArray.from_float(
+            [
+                0.0625,
+                0.125,
+                0.1875,
+                0.25,
+                0.3125,
+                0.375,
+                0.4375,
+                0.5,
+                0.5625,
+                0.625,
+                0.6875,
+                0.75,
+                0.8125,
+                0.875,
+                0.9375,
+                1,
+            ][::-1],
+            int_bits=10,
+            frac_bits=4,
+        )
+    )