Runge Kutta 2 Python ODE fluid flow
This is an interesting problem but I'm afraid that my feedback will be somewhat limited:

Consider vectorizing your velfun to be represented with matrix multiplications
Your x/y/z should be combined into one matrix
You should add axis labels, and some form of colour indicator to make your graph more legible - here I've added an arbitrary &quot;start-to-finish progress&quot; colour
Rather than label, you likely meant to use title

The code I'm suggesting produces exactly the same results as yours. The problem is that I have no idea whether those are correct. I suggest that you re-post on https://physics.stackexchange.com/ with more details on your physics problem to get confirmation as to whether your RK is implemented correctly. You should include a picture of your graph.
Suggested
from typing import Tuple

import numpy as np
import matplotlib.pyplot as plt

A = 0.2
B = 1.0
C = 0.20

N_STEP = 10_000  # stopping time
H = 0.01         # step length


def new() -&gt; np.ndarray:
    M = np.array([
        [0, 0, A],
        [B, 0, 0],
        [0, C, 0],
    ])
    N = np.array([
        [C, 0, 0],
        [0, 0, A],
        [B, 0, 0],
    ])

    # Define function to compute velocity field
    def velfun(p: np.ndarray) -&gt; np.ndarray:
        return M@np.sin(p) + N@np.cos(p)

    X = np.empty((3, N_STEP))
    a = np.array((1, 0, 0), dtype=np.float64)
    X[:, 0] = a

    for k in range(1, N_STEP):
        # Provisional Euler step
        u = velfun(a)
        p = a + H*u

        # Compute velocity at provisional point
        up = velfun(p)

        # Compute average velocity
        ua = 0.5*(u + up)

        # Make final Euler step
        # using average velocity
        a += H*ua
        X[:, k] = a

    return X


def plot(x, y, z):
    # Plot results
    fig = plt.figure()
    ax = plt.axes(projection='3d')
    ax.set_title('x (with RK2)')
    ax.set_xlabel('x')
    ax.set_ylabel('y')
    ax.set_zlabel('z')
    colour = np.linspace(0, 1, len(x))

    ax.scatter3D(x, y, z, s=0.5, c=colour)


def main():
    # plot(*old())
    plot(*new())
    plt.show()


if __name__ == '__main__':
    main()

Be warned: I have very little understanding of what this code is actually doing, so my variable names are bad.