I get this error message: Cannot cast array data from dtype('O') to dtype('float64') according to the rule 'safe'

here is my code.

import numpy as np
from scipy.optimize import minimize
import sympy as sp
sp.init_printing()
from sympy import *
from sympy import Symbol, Matrix

rom sympy import *
def make_Aij(m, n, a='a') :
    from sympy import Symbol, Matrix # just in case they aren't already loaded
    A = zeros(m, n)
    for i in range(0, m) :
        for j in range(0, n) :
            s = a+'_'+str(i)+str(j)
            exec  (s + "= Symbol('" + s + "')") # go look up what "exec" does!
            exec ("A[i, j] = " + s)
    return A
C = make_Aij(1, 2, 'c')
C

z = C[0,0]
g = C[0,0]**2-C[0,1]**2
g
h = C[0,0]
h

#Function defined
def function(h):
    return g
g

#Jacobian working with sympy
q = diff(g,C[0,0])
q

#Jacobian final
def jacobian(h):
    return q
q

Hf = diff(q,C[0,0])
Hf

#Hessianf
def Hessianf(h):
    return Hf
Hf

from scipy.optimize import Bounds
bounds = Bounds(-1, 1)

h0 = (0*C[0,1])

res = minimize(function, h0, method='trust-constr', jac=jacobian, hess=Hessianf,
               
               options={'verbose': 1}, bounds=bounds)

ERROR MESSAGE: ---------------------------------------------------------------------------

TypeError                                 Traceback (most recent call last)
<ipython-input-41-94c28f7f32b2> in <module>
----> 1 res = minimize(function, h0, method='trust-constr', jac=jacobian, hess=Hessianf,
      2 
      3                options={'verbose': 1}, bounds=bounds)

~\anaconda3\lib\site-packages\scipy\optimize\_minimize.py in minimize(fun, x0, args, method, jac, hess, hessp, bounds, constraints, tol, callback, options)
    626                                constraints, callback=callback, **options)
    627     elif meth == 'trust-constr':
--> 628         return _minimize_trustregion_constr(fun, x0, args, jac, hess, hessp,
    629                                             bounds, constraints,
    630                                             callback=callback, **options)

~\anaconda3\lib\site-packages\scipy\optimize\_trustregion_constr\minimize_trustregion_constr.py in _minimize_trustregion_constr(fun, x0, args, grad, hess, hessp, bounds, constraints, xtol, gtol, barrier_tol, sparse_jacobian, callback, maxiter, verbose, finite_diff_rel_step, initial_constr_penalty, initial_tr_radius, initial_barrier_parameter, initial_barrier_tolerance, factorization_method, disp)
    507 
    508     elif method == 'tr_interior_point':
--> 509         _, result = tr_interior_point(
    510             objective.fun, objective.grad, lagrangian_hess,
    511             n_vars, canonical.n_ineq, canonical.n_eq,

~\anaconda3\lib\site-packages\scipy\optimize\_trustregion_constr\tr_interior_point.py in tr_interior_point(fun, grad, lagr_hess, n_vars, n_ineq, n_eq, constr, jac, x0, fun0, grad0, constr_ineq0, jac_ineq0, constr_eq0, jac_eq0, stop_criteria, enforce_feasibility, xtol, state, initial_barrier_parameter, initial_tolerance, initial_penalty, initial_trust_radius, factorization_method)
    319     while True:
    320         # Solve SQP subproblem
--> 321         z, state = equality_constrained_sqp(
    322             subprob.function_and_constraints,
    323             subprob.gradient_and_jacobian,

~\anaconda3\lib\site-packages\scipy\optimize\_trustregion_constr\equality_constrained_sqp.py in equality_constrained_sqp(fun_and_constr, grad_and_jac, lagr_hess, x0, fun0, grad0, constr0, jac0, stop_criteria, state, initial_penalty, initial_trust_radius, factorization_method, trust_lb, trust_ub, scaling)
     80     Z, LS, Y = projections(A, factorization_method)
     81     # Compute least-square lagrange multipliers
---> 82     v = -LS.dot(c)
     83     # Compute Hessian
     84     H = lagr_hess(x, v)

~\anaconda3\lib\site-packages\scipy\sparse\linalg\interface.py in dot(self, x)
    416 
    417             if x.ndim == 1 or x.ndim == 2 and x.shape[1] == 1:
--> 418                 return self.matvec(x)
    419             elif x.ndim == 2:
    420                 return self.matmat(x)

~\anaconda3\lib\site-packages\scipy\sparse\linalg\interface.py in matvec(self, x)
    230             raise ValueError('dimension mismatch')
    231 
--> 232         y = self._matvec(x)
    233 
    234         if isinstance(x, np.matrix):

~\anaconda3\lib\site-packages\scipy\sparse\linalg\interface.py in _matvec(self, x)
    528 
    529     def _matvec(self, x):
--> 530         return self.__matvec_impl(x)
    531 
    532     def _rmatvec(self, x):

~\anaconda3\lib\site-packages\scipy\optimize\_trustregion_constr\projections.py in least_squares(x)
    151         # lu_sol = [aux]
    152         #          [ z ]
--> 153         lu_sol = solve(v)
    154         # return z = inv(A A.T) A x
    155         return lu_sol[n:m+n]

TypeError: Cannot cast array data from dtype('O') to dtype('float64') according to the rule 'safe'

I simply want to calculate the minimum along the c00 axis of the 2-variables function f(c00,c01) with some bounds. Therefore, i should get the value c00 as a function of C01 (c01 being a parameter here). I think the problem comes from the sympy functions associated with scipy.optimize.

Thanks in advance.