Nonlinear constraints with scipy

The problem at hand is optimization of multivariate function with nonlinear constraints. There is a differential equation (in its oversimplified form)

dy/dx = y(x)*t(x) + g(x)

I need to minimize the solution of the DE y(x), but by varying the t(x). Since it is physics under the hood, there are constraints on t(x). I successfully implemented all of them except one:

0 < t(x) < 1 for any x in range [a,b]

For certainty, the t(x) is a general polynomial:

t(x) = a0 + a1*x + a2*x**2 + a3*x**3 + a4*x**4 + a5*x**5

The x is fixed numpy.ndarray of floats and the optimization goes for coefficients a. I use scipy.optimize with trust-constr.
What I have tried so far:

1. Root finding at each step and determining the minimal/maximal value of the function using optimize.root and checking for sign changes. Return 0.5 if constraints are satisfied and numpy.inf or -1 or whatever not in [0;1] range if constraints are not satisfied. The optimizer stops soon and the function is not minimized properly.
2. Since x is fixed-length and known, I tried to define a constraint for each point, so I got N constraints where N = len(x). This works (at least look like) but takes forever for not-so large N. Also, since x is discrete and non-uniform, I can't be sure that there are no violated constraints for any x in [a,b].

EDIT #1: the minimal reproducible example

import scipy.optimize as optimize
from scipy.optimize import Bounds
import numpy as np

# some function y(x)
x = np.linspace(-np.pi,np.pi,100)
y = np.sin(x)

# polynomial t(z)
def t(a,z):
    v = 0.0;
    for ii in range(len(a)):
        v += a[ii]*z**ii
    return v

# let's minimize the sum
def targetFn(a):
    return np.sum(y*t(a,x))

# polynomial order
polyord = 3

# simple bounds to have reliable results, 
# otherwise the solution will grow toward  +-infinity
bnd = 10.0
bounds = Bounds([-bnd for i in range(polyord+1)],
                [bnd for i in range(polyord+1)])

res = optimize.minimize(targetFn, [1.0 for i in range(polyord+1)], 
                        bounds = bounds)

if np.max(t(res.x,x))>200:
    print('max constraint violated!')
if np.min(t(res.x,x))<-100:
    print('min constraint violated!')

In the reproducible example given above, let the constraints to be that the value of the polynomial t(a,x) is in range [-100;200] for the given x.

So the question is: how does one properly define a constraint to tell the optimizer that the function's values must be constrained for the given range of arguments?