Stock level optimisation

I am trying to solve the following problem:

Given the dataframe with columns:

• location
• product
• current_stock
• demand
• price

I want to optimize the unmet demand given additional constraints:
- the number of shipments from one location cannot exceed a given threshold like 50

- the value of each package sent cannot be lower then $250
I use the code below. I think it works quite well if the number of products/locations is not that big like 50 or 100. But when I run it on a full problem with more than 100 locations and couple thousands of products the only solution it finds is 0.

I don't get it, why there is a non zero solution which improves the satisfied demand for 50 products and it cannot improve for a full set. I expected that in worst case it should find the solution for 50 products leaving all other 0.

def optimize_demand_gurobi_with_more_constraints(min_stock, current_stock, demand, price, locations, products,max_time):
    log_system_status("Opt started")
    m = Model("Unsatisfied_Demand_Optimization")
    m.setParam('OutputFlag', 1)  # Enable output for more verbose logging

    # m.params.Method = 1
    # m.params.SoftMemLimit = 170
    m.params.Presolve = 2
    m.params.PreSparsify = 1
    # m.params.NodefileStart = 0.5
    # m.params.Threads = 
    m.params.LogFile = add_timestamp_to_filename('logs/GRB_opt.txt')
    m.params.DegenMoves = 0
    m.params.TimeLimit = max_time
    log_system_status("x")
    # Defining variables: Only create shipment variables x where conditions are met
    x = m.addVars(
        [(loc_a, loc_b, p) for loc_a in tqdm(locations) for loc_b in locations for p in products
         if loc_a != loc_b and current_stock[loc_a, p] > demand[loc_a, p] and current_stock[loc_b, p] < demand[loc_b, p]],
        vtype=GRB.INTEGER, name="Ship"
    )
    log_system_status("unmet")
    unmet = m.addVars([(loc_b, p) for loc_b in tqdm(locations) for p in products],
                      vtype=GRB.INTEGER, name="Unmet")
    log_system_status("y")
    y = m.addVars([(loc_a, loc_b, p) for loc_a in tqdm(locations) for loc_b in locations for p in products if loc_a != loc_b],
                  vtype=GRB.BINARY, name="NoCrossShip")

    # log_system_status("obj")
    # Objective: Minimize the total unmet demand
    # m.setObjective(quicksum(unmet[loc_b, p] for loc_b in tqdm(locations) for p in products), GRB.MINIMIZE)

    log_system_status("Hints")
    for k, v in tqdm(x.items()):
        x[k].VarHintVal = shipments.get(k,0)
        x[k].VarHintPri = 10
    
    for k, v in tqdm(y.items()):
        if shipments.get(k,0) == 0:
            y[k].VarHintVal = 0
            y[k].VarHintPri = 5
        else:
            y[k].VarHintVal = 1
            y[k].VarHintPri = 5
    
    log_system_status("stock & demand")
    # Constraints to ensure shipments do not exceed stock available or demand
    for loc_a in tqdm(locations):
        for p in products:
            if any((loc_a, loc_b, p) in x for loc_b in locations if loc_a != loc_b):
                m.addConstr(quicksum(x[loc_a, loc_b, p] for loc_b in locations if loc_a != loc_b and (loc_a, loc_b, p) in x) <= current_stock[loc_a, p])

    for loc_b in tqdm(locations):
        for p in products:
            if any((loc_a, loc_b, p) in x for loc_a in locations if loc_a != loc_b):
                m.addConstr(quicksum(x[loc_a, loc_b, p] for loc_a in locations if loc_a != loc_b and (loc_a, loc_b, p) in x) <= demand[loc_b, p])
            m.addConstr(unmet[loc_b, p] >= price[p]*(demand[loc_b, p]-current_stock[loc_b, p]- quicksum(x[loc_a, loc_b, p] for loc_a in locations if loc_a != loc_b and (loc_a, loc_b, p) in x)))
            m.addConstr(unmet[loc_b, p] >=0 )
    
    log_system_status("Crosshipment")
    # Constraint to prevent shipment in both directions for the same product
    for loc_a in tqdm(locations):
        for loc_b in locations:
            if loc_a != loc_b:
                for p in products:
                    if (loc_a, loc_b, p) in y and (loc_b, loc_a, p) in y:
                        m.addConstr(y[loc_a, loc_b, p] + y[loc_b, loc_a, p] <= 1)
                        if (loc_a, loc_b, p) in x:
                            m.addConstr(x[loc_a, loc_b, p] <= y[loc_a, loc_b, p] * max(current_stock[loc_a, p], demand[loc_b, p]))
                        if (loc_b, loc_a, p) in x:
                            m.addConstr(x[loc_b, loc_a, p] <= y[loc_b, loc_a, p] * max(current_stock[loc_b, p], demand[loc_a, p]))

    count_parcels = m.addVars(
        locations,
        vtype=GRB.INTEGER, name="Ship"
    )
    
    total_value = m.addVars([(loc_a, loc_b) for loc_a in tqdm(locations) for loc_b in locations],
                  vtype=GRB.INTEGER, name="toal_value")
    
    shipment = m.addVars([(loc_a, loc_b) for loc_a in tqdm(locations) for loc_b in locations],
                  vtype=GRB.BINARY, name="shipment")
    
    
    
    for loc_a in tqdm(locations):
        for loc_b in locations:
            if loc_a != loc_b:
                m.addConstr(shipment[loc_a, loc_b] == max_(y[loc_a, loc_b, p] for p in products))
            else:
                m.addConstr(shipment[loc_a, loc_b] == 0)
    
    for loc_a in tqdm(locations):
        m.addConstr(count_parcels[loc_a] == quicksum(shipment[loc_a, loc_b] for loc_b in locations ))
        m.addConstr(count_parcels[loc_a] <= 50)
        
    
    for loc_a in tqdm(locations):
        for loc_b in locations:
            if loc_a != loc_b:
                m.addConstr(total_value[loc_a, loc_b] == quicksum(x.get((loc_a, loc_b, p),0)*price[p] for p in products))
                m.addConstr(total_value[loc_a, loc_b] >= 250*shipment[loc_a, loc_b])
                
                
    log_system_status("Vars and const defined")
    log_system_status("obj")
    # Objective: Minimize the total unmet demand
    # m.setObjective(-unmet.sum(), GRB.MAXIMIZE)
    log_system_status("obj")
    # Objective: Minimize the total unmet demand
    m.setObjective(
        quicksum(x[src, dest, prod] for src in locations for dest in locations if src != dest for prod in products if (src, dest, prod) in x),
        GRB.MAXIMIZE
    )

    
    log_system_status("Optimizing")
    
    m._solution_count = 0
    m.optimize()
    
    if m.status in [GRB.OPTIMAL, GRB.TIME_LIMIT, GRB.ITERATION_LIMIT]:
        log_system_status("GRB is optimal")
        shipment_plan = {k: v.X for k, v in x.items() if v.X > 0}
        m.write(add_timestamp_to_filename('model/with_value_constraint.mps'))
        m.write(add_timestamp_to_filename('model/model.sol'))
        return shipment_plan, m
    else:
        log_system_status("Optimization was not successful. Status code:"+str(m.status))
        if m.status in [GRB.INFEASIBLE, GRB.INF_OR_UNBD]:
            m.computeIIS()
            m.write("model.ilp")
        return {}