ODE fitting with Symfit - ValueError: operands could not be broadcast together with remapped shapes

Question

I'm trying to solve and fit these coupled Differential equations with Symfit. I've tried following the documentation for ODE fitting, looked at this answer to a similar question and searched the error in other contexts, but without any luck. My code is supposed be used on experimental data, but as a start I've just created somewhat similar data. The relevant code is:

#Imports
import numpy as np
import matplotlib.pyplot as plt
from symfit import variables, parameters, Parameter, Fit, D, ODEModel, cos, sin, atan, asin, sqrt

#Constants
m  = 0.01      #Unit: kg    - Mass of object
mt = 0.5 * m   #Unit: kg    - Related to moment of inertia (\tilde{m] in the equations)
R  = 0.1       #Unit: m     - Radius of object
g  = 9.82      #Unit: m/s^2 - gravitational acceleration 
#Variables and parameters
x_1,x_2,x_3,x_4, t = variables("x_1,x_2,x_3,x_4, t")  #x-position, y-position, x-velocity, y-velocity, time
p         = Parameter("p", min=1, max=4)
c, d      = parameters("c,d")       #d is \gamma in the equations

#Generate proxy data
t_vals  = np.linspace(0,10,100)            #Choosen time interval
x1_data = 8*np.sin(0.15*t_vals)+7*t_vals   #x-data proxy
x2_data = 10*np.sin(0.15* t_vals)-2*t_vals #y-data proxy
x3_data = 0.15*8*np.cos(0.15*t_vals)-7     #vx-data proxy
x4_data = 0.15*10*np.cos(t_vals)-2         #vy-data proxy

#Model
model_dict = {
    D(x_1, t) : x_3,
    D(x_2, t) : x_4,
    D(x_3, t) : -d/mt * (x_3**2 + x_4**2)**((p-1)/2) * x_3 + c*(x_3**2+x_4**2)/(mt*R) * cos(atan(x_1/x_2) + asin(R/sqrt(x_1**2+x_2**2))),
    D(x_4, t) : -d/mt * (x_3**2 + x_4**2)**((p-1)/2) * x_4 + c*(x_3**2+x_4**2)/(mt*R) * sin(atan(x_1/x_2) + asin(R/sqrt(x_1**2+x_2**2)))+m/mt * g,
}

ode_model  = ODEModel(model_dict, initial={t: 0.0, x_1: 0, x_2: 3, x_3: 0, x_4: 1})

fit        = Fit(ode_model, t=tdata, x_1=x1_data, x_2=x2_data, x_3=x3_data, x_4=x4_data)
fit_result = fit.execute()

But I get the following error:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-8-ef065dc7aeb9> in <module>
     27 
     28 fit        = Fit(ode_model, t=tdata, x_1=x1_data, x_2=x2_data, x_3=x3_data, x_4=x4_data)
---> 29 fit_result = fit.execute()
     30 
     31 t_vec      = np.linspace(0,10,1000)

[...]

ValueError: operands could not be broadcast together with remapped shapes [original->remapped]: (5,) and requested shape (100,)

I'm rather inexperienced in using Symfit, but I suppose the two arrays the error is referring to, is the five variables and the 100 datapoints each of them have. I have no idea what causes this error, except maybe a lack of initial guesses and badly generated data. The problem simply was that I used the wrong time array

EDIT: As user mikuszefski pointed out, I had used tdata where t_vals should've been used. This solved my initial error. The same user also correctly predicted numerical troubles, as some fitted values are of magnitude 10^300. I will edit my post once I find a solution. I have also written my imports and defined g.