What's wrong with my gradient descent? (simple, one variable)
You set your model parameters using np.array([[0], [0]]), meaning that your parameter array is of shape (2, 1). Changing this to np.array([0, 0]) (with shape (2,)) decreases the loss over time as expected using a learning rate $\alpha$ of 0.01. Looking more closely at the difference of using the two, I see that your way returns way higher losses (exactly 100x). Using a smaller $\alpha$, e.g. 0.0001, with your implementation also works:

Therefore the main reason of your implementation is your choice of $\alpha$, which is too high, preventing the model to reach the optimum.
EDIT: The issues is twofold, the first is indeed the shape or your theta array. The second has to do with the fact that you are using the theta array in your for loop which always contains the original theta parameters, as you are not storing and using the updated theta1 and theta2. The following code seems to give a good result of $\theta_1 \approx -3.2$ and $\theta_2 \approx 1.12$ using 1000 iterations with $\alpha=0.01$
import numpy as np

# load and transform data
data = np.genfromtxt(&quot;data.txt&quot;, delimiter=&quot;,&quot;)
X = np.stack((np.ones(data.shape[0]), data[:,0]), axis=-1)
y = data[:,1]

m = X.shape[0]
iter = 1000
 
# calculating the cost function
theta = np.array([[0.3], [0.5]])
def h(X, theta):
    return np.matmul(X, theta)
 
def cost(X, y, theta, m):
    return 1/(2* m) * np.sum((h(X, theta)-y)**2)
 
print(cost(X, y, theta, m))
# 1-variable gradient decent
cost_vals = []
def grad_decent(X, y, theta, m, alpha):
    theta1 = theta[0, 0]
    theta2 = theta[1, 0]
    for i in range(0, iter):
        theta = np.array([theta1, theta2])
        cost_vals.append(cost(X, y, np.array([theta1, theta2]), m))
        theta1 = theta1 - alpha/m * np.sum((h(X, theta)-y) * X[:, 0])
        theta2 = theta2 - alpha/m * np.sum((h(X, theta)-y) * X[:, 1])
        # print(f&quot;Iteration {i}: &quot;, theta1, &quot;, &quot;, theta2)
        
    return np.array([[theta1], [theta2]])

grad_decent(X, y, theta, m, 0.01)
# 1740.4098044406999
# array([[-3.20070025],
#        [ 1.12320526]])