Problemas com Docplex e retorno de chamada

Estou tentando implementar um callback personalizado que resolva um subproblema durante o processo de otimização. Eu quero:

1. Obtenha os valores de certas variáveis ​​da solução do modelo principal.

2. Resolva um subproblema (um problema de otimização separado) com base nesses valores.

3. Se a solução do subproblema sugerir um corte, adicione-o ao modelo principal como uma nova restrição.

Vou compartilhar meu código aqui, eu realmente queria fazer isso em python, mas estou muito perdido agora, também tenho tentado xpress e similares, mas a documentação é inútil, qualquer ajuda seria muito apreciada.

from docplex.mp.model import Model
from cplex.callbacks import UserCutCallback
from docplex.mp.callbacks.cb_mixin import ConstraintCallbackMixin
import random

class CustomCutCallback(ConstraintCallbackMixin, UserCutCallback):
    def __init__(self, env):
        UserCutCallback.__init__(self, env)
        ConstraintCallbackMixin.__init__(self)
        self.eps = 1e-6
        self.nb_cuts = 0
        self.cts = []

    def add_cut_constraint(self, cut):
        self.register_constraint(cut)

    def __call__(self, context):
        """
        This method is invoked by CPLEX during optimization.
        It receives the 'context' to interact with the model and variables.
        """
        print("Callback called with context")
        m = context.model  # The main model from the context
        m_sol = m.solution

        x_values = {var: m_sol.get_value(var) for var in m.variables if var.name.startswith('x')}
        w_values = {var: m_sol.get_value(var) for var in m.variables if var.name.startswith('w')}
        m2_cuts = self.solve_subproblem(x_values, w_values, m, context)

    def solve_subproblem(self, x_values, w_values, m, context):
        """
        Solves the subproblem and generates cuts if necessary.
        """
        m2 = Model(name='subproblem')

        client_range = range(len(x_values))
        deposit_range = range(len(w_values))
        plant_range = range(len(w_values[0]))

        alpha = m2.continuous_var_matrix(client_range, deposit_range, name='alpha', lb=0)
        beta = m2.continuous_var_matrix(deposit_range, plant_range, name='beta', lb=0)

        m2.add_constraints(alpha[i, j] + (x_values.get((i, j), 0) * beta[j, k]) <= x_values.get((i, j), 0) * w_values.get((j, k), 0)
                           for i in client_range for j in deposit_range for k in plant_range)

        m2.maximize(m2.sum(alpha[i, j] * x_values.get((i, j), 0) for i in client_range for j in deposit_range) + 
                    m2.sum(beta[j, k] * w_values.get((j, k), 0) for j in deposit_range for k in plant_range))

        m2.solve()
        print(m2.solution)

        # Here, you perform an evaluation of the cut values
        for i in client_range:
            for j in deposit_range:
                for k in plant_range:
                    if  sum(m2.solution.get_value(alpha[i, j]) * x_values.get((i, j), 0) for i in client_range for j in deposit_range) + \
                    sum(m2.solution.get_value(beta[j, k]) * w_values.get((j, k), 0) for j in deposit_range for k in plant_range) > \
                    sum(transport_cost_deposit_client[j][i] * d[i] * x[i, j] for i in client_range for j in deposit_range) + \
                    sum(transport_cost_plant_deposit[j][k] * w[j, k] for j in deposit_range for k in plant_range):
                        m.add_constraint(sum(m2.solution.get_value(alpha[i, j]) * x[i, j] for i in client_range for j in deposit_range) + \
                                         sum(m2.solution.get_value(beta[j, k]) * w[j, k] for j in deposit_range for k in plant_range))


# Main model
def build_location_model(transport_cost_deposit_client, transport_cost_plant_deposit, p, q, d, **kwargs):
    m = Model(name='location', **kwargs)

    num_deposits = len(transport_cost_deposit_client)
    num_plants = len(transport_cost_plant_deposit[0])
    num_clients = len(d)
    
    deposit_range = range(num_deposits)
    plant_range = range(num_plants)
    client_range = range(num_clients)

    x = m.binary_var_matrix(client_range, deposit_range, name='x')
    w = m.integer_var_matrix(deposit_range, plant_range, name='w', lb=0)
    y = m.binary_var_list(deposit_range, name='y')
    h = m.binary_var_list(plant_range, name='h')

    m.add_constraints(m.sum(x[i, j] for j in deposit_range) == 1 for i in client_range)
    m.add_constraints(m.sum(w[j, k] for k in plant_range) == y[j] for j in deposit_range)
    m.add_constraints(x[i, j] <= y[j] for i in client_range for j in deposit_range)
    m.add_constraint(m.sum(h[k] for k in plant_range) == q)
    m.add_constraint(m.sum(y[j] for j in deposit_range) == p)
    m.add_constraints(m.sum(d[i] * x[i,j] for i in client_range) == m.sum(w[j, k] for k in plant_range) for j in deposit_range)
    m.add_constraints(w[j, k] <= (m.sum(d[i] for i in client_range) * h[k]) for j in deposit_range for k in plant_range)

    transport_cost = m.sum(transport_cost_deposit_client[j][i] * d[i] * x[i, j] for i in client_range for j in deposit_range) + \
                        m.sum(transport_cost_plant_deposit[j][k] * w[j, k] for j in deposit_range for k in plant_range)
    m.minimize(transport_cost)
    m.parameters.preprocessing.presolve = 0

    # Register the callback
    cut_cb = m.register_callback(CustomCutCallback)

    # Configure CPLEX parameters for cuts
    params = m.parameters
    params.mip.cuts.mircut = -1

    m.solve()

    return m



# Test function
def solve_model():
    num_deposits = 10
    num_plants = 4
    num_clients = 20

    TRANSPORT_COST_DEPOSITS_CLIENTS = [
        [random.randint(20, 100) for _ in range(num_clients)] for _ in range(num_deposits)
    ]

    TRANSPORT_COST_PLANTS_DEPOSITS = [
        [random.randint(30, 80) for _ in range(num_plants)] for _ in range(num_deposits)
    ]

    p = 5
    q = 3
    d = [random.randint(5, 20) for _ in range(num_clients)]

    m = build_location_model(TRANSPORT_COST_DEPOSITS_CLIENTS, TRANSPORT_COST_PLANTS_DEPOSITS, p, q, d)
    if m.solution is None:
        print("No valid solution found.")
        return None
    print(m.solution)
    return m


if __name__ == "__main__":
    solve_model()