Algorithm for maximizing shipping profit with limitations on mass and cost

Question

The title isn't very helpful, because I'm not sure what I'm trying to say exactly. I'm sure an algorithm for this must exist, but I can't remember. Note: not a homework problem, I finished school a very long time ago.

So here's the problem:

• We're doing a shipping and trading job, trying to maximize profits
• We have a list of items that we can ship in a truck. Each item has:
  • A buy price (at the source)
  • A sell price (at the destination)
  • A per-unit mass
  • An upper limit on how many can be purchased
• Our truck is limited in the amount of mass it can carry
• We have an upper limit on how much we're allowed to "invest" (spend on items at the source).
• We want to maximize the profit for our job (buy at the source, transport, sell at the destination).

If there were only one limit (total mass, or total investment), it would be easy, but I'm not sure how to approach this when there are two.

The equation for calculating profit would be:

profit = ItemA['quantity'] * (ItemA['sell_price'] - ItemA['buy_price']) + ItemB['quantity'] * (ItemB['sell_price'] - ItemB['buy_price']) + ...

So I'm trying to choose which items, and the quantity of each item, that should be purchased in order to maximize the profit.

Are there any existing, known algorithms for solving this? Likely some sort of mathematical optimization problem? I'm using Python, so I'm thinking that the mystic package might be appropriate, but I'm not sure how I'd configure it.

Answer 1

You can try the framework optuna for hyperparameter tuning.

Here is an example code that you can try. Products are named product1 etc found in parameters.json file. Data values are just assumptions.

Study/optimization session are now saved in sqlite db. This will support interrupt and resume. See version log in the code.

parameters.json

{
    "study_name": "st5_tpe",
    "sampler": "tpe",
    "trials": 1000,
    "max_purchase": 7000,
    "min_weight_no_cost": 1000,
    "high_weight_additional_cost": 0.5,
    "trucks": {
        "smalltruck": {
            "maxmass": 1000,
            "cost": 75
        },
        "mediumtruck": {
            "maxmass": 2000,
            "cost": 150
        },
        "bigtruck": {
            "maxmass": 5000,
            "cost": 400
        }
    },
    "products": {
        "product1_qty": {
            "min": 20,
            "max": 100,
            "massperunit": 2,
            "buyprice": 5,
            "sellprice": 8
        },
        "product2_qty": {
            "min": 20,
            "max": 100,
            "massperunit": 4,
            "buyprice": 6,
            "sellprice": 10
        },
        "product3_qty": {
            "min": 20,
            "max": 100,
            "massperunit": 1,
            "buyprice": 4,
            "sellprice": 6
        },
        "product4_qty": {
            "min": 20,
            "max": 100,
            "massperunit": 2,
            "buyprice": 7,
            "sellprice": 10
        },
        "product5_qty": {
            "min": 20,
            "max": 100,
            "massperunit": 2,
            "buyprice": 5,
            "sellprice": 8
        },
        "product6_qty": {
            "min": 20,
            "max": 100,
            "massperunit": 1,
            "buyprice": 5,
            "sellprice": 7
        },
        "product7_qty": {
            "min": 20,
            "max": 100,
            "massperunit": 1,
            "buyprice": 8,
            "sellprice": 12
        }
    }
}

Code

"""
shipping_trading.py


version 0.7.0
    * Calculate and show ROI (return of investment) and other info.
    * Add user attribute to get other costs.
    * Raise exception when max_purchase key is missing in parameters.json file.
    * Continue the study even when trucks key is missing in parameters.json file.
    
version 0.6.0
    * Save study/optimization session in sqlite db, with this it can now supports interrupt and resume.
      When study session is interrupted it can be resumed later using data from previous session.
    * Add study_name key in parameters.json file. Sqlite db name is based on study_name. If you
      want new study/optimization session, modify the study_name. If you are re-running the
      same study_name, it will run and continue from previous session. Example:
      study_name=st8, sqlite_dbname=mydb_st8.db
      By default study_name is example_study when you remove study_name key in parameters.json file.
    * Remove printing in console on truck info.

version 0.5.0
    * Replace kg with qty in parameters.json file.
    * Add massperunit in the product.
    * Optimize qty not mass.
    * Refactor

version 0.4.0
    * Add truck size optimization. It is contrained by the cost of using truck as well as the max kg capacity.
      The optimizer may suggest a medium instead of a big truck if profit is higher as big truck is expensive.
      profit = profit - truck_cost - other_costs
    * Modify parameters.json file, trucks key is added.

version 0.3.0
    * Read sampler, and number of trials from parameters.json file.
      User inputs can now be processed from that file.

version 0.2.0
    * Read a new parameters.json format.
    * Refactor get_parameters().

version 0.1.0
    * Add additional cost if total product weight is high.
"""


__version__ = '0.7.0'


import json

import optuna


def get_parameters():
    """
    Read parameters.json file to get the parameters to optimize, etc.
    """
    fn = 'parameters.json'
    products, trucks = {}, {}

    with open(fn) as json_file:
        values = json.load(json_file)

        max_purchase = values.get('max_purchase', None)
        if max_purchase is None:
            raise Exception('Missing max_purchase, please specify max_purchase in json file, i.e "max_purchase": 1000')

        study_name = values.get('study_name', "example_study")
        sampler = values.get('sampler', "tpe")
        trials = values.get('trials', 100)
        min_weight_no_cost = values.get('min_weight_no_cost', None)
        high_weight_additional_cost = values.get('high_weight_additional_cost', None)
        products = values.get('products', None)
        trucks = values.get('trucks', None)

    return (products, trucks, sampler, trials, max_purchase, min_weight_no_cost, high_weight_additional_cost, study_name)


def objective(trial):
    """
    Maximize profit.
    """
    gp = get_parameters()
    (products, trucks, _, _, max_purchase,
        min_weight_no_cost, high_weight_additional_cost, _) = gp

    # Ask the optimizer the product qty to use try.
    new_param = {}    
    for k, v in products.items():
        suggested_value = trial.suggest_int(k, v['min'], v['max'])  # get suggested value from sampler
        new_param.update({k: {'suggested': suggested_value,
                               'massperunit': v['massperunit'],
                               'buyprice': v['buyprice'],
                               'sellprice': v['sellprice']}})

    # Ask the sampler which truck to use, small, medium ....
    truck_max_wt, truck_cost = None, None
    if trucks is not None:
        truck = trial.suggest_categorical("truck", list(trucks.keys()))

        # Define truck limits based on suggested truck size.
        truck_max_wt = trucks[truck]['maxmass']
        truck_cost = trucks[truck]['cost']

    # If total wt or total amount is exceeded, we return a 0 profit.
    total_wt, total_buy, profit = 0, 0, 0
    for k, v in new_param.items():
        total_wt += v['suggested'] * v['massperunit']
        total_buy += v['suggested'] * v['buyprice']
        profit += v['suggested'] * (v['sellprice'] - v['buyprice'])

    # (1) Truck mass limit
    if truck_max_wt is not None:
        if total_wt > truck_max_wt:
            return 0

    # (2) Purchase limit amount
    if max_purchase is not None:
        if total_buy > max_purchase:
            return 0

    # Cost for higher transport weight
    cost_high_weight = 0
    if min_weight_no_cost is not None and high_weight_additional_cost is not None:
        excess_weight = total_wt - min_weight_no_cost
        if excess_weight > 0:
            cost_high_weight += (total_wt - min_weight_no_cost) * high_weight_additional_cost

    # Cost for using a truck, can be small, medium etc.
    cost_truck_usage = 0
    if truck_cost is not None:
        cost_truck_usage += truck_cost

    # Total cost
    other_costs = cost_high_weight + cost_truck_usage
    trial.set_user_attr("other_costs", other_costs)

    # Adjust profit
    profit = profit - other_costs

    # Send this profit to optimizer so that it will consider this value
    # in its optimization algo and would suggest a better value next time we ask again.
    return profit


def return_of_investment(study, products):
    """
    Returns ROI.

    ROI = Return Of Investment
    ROI = 100 * profit/costs
    """
    product_sales, product_costs = 0, 0
    for (k, v), (k1, v1) in zip(products.items(), study.best_params.items()):
        if k == 'truck':
            continue
        assert k == k1
        product_sales += v1 * v['sellprice']
        product_costs += v1 * v['buyprice']
        
    other_costs = study.best_trial.user_attrs['other_costs']
    total_costs = product_costs + other_costs

    calculated_profit = product_sales - total_costs
    study_profit = study.best_trial.values[0]
    assert calculated_profit == study_profit
    
    return_of_investment = 100 * calculated_profit/total_costs

    return return_of_investment, product_sales, product_costs, other_costs


def main():
    # Read parameters.json file for user data input.
    gp = get_parameters()
    (products, trucks, optsampler, num_trials,
        max_purchase, _, _, study_name) = gp

    # Location of sqlite db where optimization session data are saved.
    sqlite_dbname = f'sqlite:///mydb_{study_name}.db'

    # Available samplers to use:
    # https://optuna.readthedocs.io/en/stable/reference/samplers.html
    # https://optuna.readthedocs.io/en/stable/reference/generated/optuna.integration.SkoptSampler.html
    # https://optuna.readthedocs.io/en/stable/reference/generated/optuna.integration.BoTorchSampler.html
    if optsampler.lower() == 'cmaes':
        sampler = optuna.samplers.CmaEsSampler(n_startup_trials=1, seed=100)
    elif optsampler.lower() == 'tpe':
        sampler = optuna.samplers.TPESampler(n_startup_trials=10, multivariate=False, group=False, seed=100, n_ei_candidates=24)
    else:
        print(f'Warning, {optsampler} is not supported, we will be using tpe sampler instead.')
        optsampler = 'tpe'
        sampler = optuna.samplers.TPESampler(n_startup_trials=10, multivariate=False, group=False, seed=100, n_ei_candidates=24)

    # Store optimization in storage and supports interrupt/resume.
    study = optuna.create_study(storage=sqlite_dbname, sampler=sampler, study_name=study_name, load_if_exists=True, direction='maximize')
    study.optimize(objective, n_trials=num_trials)

    # Show summary and best parameter values to maximize profit.
    print()
    print(f'study_name: {study_name}')
    print(f'sqlite dbname: {sqlite_dbname}')
    print(f'sampler: {optsampler}')
    print(f'trials: {num_trials}')
    print()

    print(f'Max Purchase Amount: {max_purchase}')
    print()

    print('Products being optimized:')
    for k, v in products.items():
        print(f'{k}: {v}')
    print()

    if trucks is not None:
        print('Trucks being optimized:')
        for k, v in trucks.items():
            print(f'{k}: {v}')
        print()

    print('Study/Optimization results:')
    objective_name = 'profit'
    print(f'best parameter value : {study.best_params}')
    print(f'best value           : {study.best_trial.values[0]}')
    print(f'best trial           : {study.best_trial.number}')
    print(f'objective            : {objective_name}')
    print()

    # Show other info like roi, etc.
    roi, product_sales, product_costs, other_costs = return_of_investment(study, products)
    print('Other info.:')    
    print(f'Return Of Investment : {roi:0.2f}%, profit/costs')
    print(f'Product Sales        : {product_sales:0.2f}')
    print(f'Product Costs        : {product_costs:0.2f}')
    print(f'Other Costs          : {other_costs:0.2f}')
    print(f'Total Costs          : {product_costs + other_costs:0.2f}')
    print(f'Profit               : {product_sales - (product_costs + other_costs):0.2f}')
    print(f'Capital              : {max_purchase:0.2f}')
    print(f'Total Spent          : {product_costs + other_costs:0.2f} ({100*(product_costs + other_costs)/max_purchase:0.2f}% of Capital)')
    print(f'Capital Balance      : {max_purchase - product_costs - other_costs:0.2f}')
    print()


if __name__ == '__main__':
    main()

Output

study_name: st5_tpe
sqlite dbname: sqlite:///mydb_st5_tpe.db
sampler: tpe
trials: 1000

Max Purchase Amount: 7000

Products being optimized:
product1_qty: {'min': 20, 'max': 100, 'massperunit': 2, 'buyprice': 5, 'sellprice': 8}
product2_qty: {'min': 20, 'max': 100, 'massperunit': 4, 'buyprice': 6, 'sellprice': 10}
product3_qty: {'min': 20, 'max': 100, 'massperunit': 1, 'buyprice': 4, 'sellprice': 6}
product4_qty: {'min': 20, 'max': 100, 'massperunit': 2, 'buyprice': 7, 'sellprice': 10}
product5_qty: {'min': 20, 'max': 100, 'massperunit': 2, 'buyprice': 5, 'sellprice': 8}
product6_qty: {'min': 20, 'max': 100, 'massperunit': 1, 'buyprice': 5, 'sellprice': 7}
product7_qty: {'min': 20, 'max': 100, 'massperunit': 1, 'buyprice': 8, 'sellprice': 12}

Trucks being optimized:
smalltruck: {'maxmass': 1000, 'cost': 75}
mediumtruck: {'maxmass': 2000, 'cost': 150}
bigtruck: {'maxmass': 5000, 'cost': 400}

Study/Optimization results:
best parameter value : {'product1_qty': 99, 'product2_qty': 96, 'product3_qty': 93, 'product4_qty': 96, 'product5_qty': 100, 'product6_qty': 100, 'product7_qty': 100, 'truck': 'mediumtruck'}
best value           : 1771.5
best trial           : 865
objective            : profit

Other info.:
Return Of Investment : 42.19%, profit/costs
Product Sales        : 5970.00
Product Costs        : 3915.00
Other Costs          : 283.50
Total Costs          : 4198.50
Profit               : 1771.50
Capital              : 7000.00
Total Spent          : 4198.50 (59.98% of Capital)
Capital Balance      : 2801.50

If you increase the number of trials the program might be able to find a more profitable parameter values.

Answer 2

Another option is by using scipy. The sample below contains 3 products, that can be scaled of course. The constrains are the purchase and max truck mass capacity.

code

"""
shipping_trading_solver.py

Ref: https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.minimize.html#scipy.optimize.minimize
"""


from scipy.optimize import minimize


# Constants
sellprice = [8, 7, 10]
buyprice = [6, 5, 6]
mass_per_unit = [1, 2, 3]

purchase_limit = 100
truck_mass_limit = 70


def objective(x):
    """
    objective, return value as negative to maximize.
    x: quantity
    """
    profit = 0
    for (v, s, b) in zip(x, sellprice, buyprice):
        profit += v * (s - b)

    return -profit


def purchase_cons(x):
    """
    Used for constrain
    x: quantity
    """
    purchases = 0
    for (v, b) in zip(x, buyprice):
        purchases += v * b
    
    return purchase_limit - purchases  # not negative


def mass_cons(x):
    """
    Used for constrain
    mass = qty * mass/qty
    x: quantity
    """
    mass = 0
    for (v, m) in zip(x, mass_per_unit):
        mass += v * m
    
    return truck_mass_limit - mass  # not negative


def profit_cons(x):
    """
    Used for constrain
    x: quantity
    """
    profit = 0
    for (v, s, b) in zip(x, sellprice, buyprice):
        profit += v * (s - b)

    return profit  # not negative


def main():
    # Define constrained. Note: ineq=non-negative, eq=zero
    cons = (
        {'type': 'ineq', 'fun': purchase_cons},
        {'type': 'ineq', 'fun': mass_cons},
        {'type': 'ineq', 'fun': profit_cons}
    )

    # Bounds of product quantity, (min,max)
    bound = ((0, 50), (0, 20), (0, 30))

    # Initial values
    init_values = (0, 0, 0)

    # Start minimizing
    # SLSQP = Sequential Least Squares Programming
    res = minimize(objective, init_values, method='SLSQP', bounds=bound, constraints=cons)

    # Show summary
    print('Results summary:')
    print(f'optimization message: {res.message}')
    print(f'success status: {res.success}')
    print(f'profit: {sum([(s-b) * int(x) for (x, s, b) in zip(res.x, sellprice, buyprice)]):0.1f}')
    print(f'best param values: {[int(v) for v in res.x]}')
    print()

    # Verify results
    print('Verify purchase and mass limits:')

    # (1) purchases
    total_purchases = 0
    for (qty, b) in zip(res.x, buyprice):
        total_purchases += int(qty) * b
    print(f'actual total_purchases: {total_purchases:0.1f}, purchase_limit: {purchase_limit}')

    # (2) mass
    total_mass = 0    
    for (qty, m) in zip(res.x, mass_per_unit):
        total_mass += int(qty) * m
    print(f'actual total_mass: {total_mass:0.1f}, truck_mass_limit: {truck_mass_limit}')


if __name__ == '__main__':
    main()

output

Results summary:
optimization message: Optimization terminated successfully
success status: True
profit: 64.0
best param values: [0, 0, 16]

Verify purchase and mass limits:
actual total_purchases: 96.0, purchase_limit: 100
actual total_mass: 48.0, truck_mass_limit: 70

Answer 3

I'm the mystic author. First off, mystic is not the best code to use with this problem... a good linear MIP solver like the one in OR-Tools would be a better choice. Mystic will solve MIP/LP problems reliably, just not as fast as OR-Tools. In terms of speed, mystic is about as fast as scipy.optimize. Mystic can slow down as the constraints become more nonlinear, complex, and tightly constrained (note that other codes generally fail in this case while mystic does not). Below, I'll use a differential evolution solver (which is slower, but more robust than SLSQP).

Note that as soon as you have one or more nonlinear constraints, you should definitely use mystic... as mystic is built for global optimization with nonlinear constraints. Or, if instead of having a fixed pricing model, you have some market fluctuations in the model, and thus uncertainty... and want to maximize the expected profit, or even better to build a model of profit that minimizes risk, then you definitely should use mystic. OR-tools and other LP/QP codes will have to, at best, approximate the problem as linear or quadratic -- which may not be practical.

Regardless. as you asked about using mystic on this problem, here's one of many ways to get to a solution with mystic:

import mystic as my
import mystic.symbolic as ms
import mystic.constraints as mc

class item(object):
    def __init__(self, id, mass, buy, net, limit):
        self.id = id
        self.mass = mass
        self.buy = buy
        self.net = net
        self.limit = limit
    def __repr__(self):
        return 'item(%s, mass=%s, buy=%s, net=%s, limit=%s)' % (self.id, self.mass, self.buy, self.net, self.limit)

# data
masses = [10, 15, 20, 18, 34, 75, 11, 49, 68, 55]
buys = [123, 104, 149, 175, 199, 120, 164, 136, 194, 111]
nets = [13, 24, 10, 29, 29, 39, 28, 35, 33, 39]
limits = [300, 500, 200, 300, 200, 350, 100, 600, 1000, 50]
ids = range(len(limits))

# maxima
_load = 75000  # max limit on mass can carry
_spend = 350000  # max limit to spend at source

# items
items = [item(*i) for i in zip(ids, masses, buys, nets, limits)]

 # profit
def fixnet(net):
    def profit(x):
        return sum(xi*pi for xi,pi in zip(x,net))
    return profit

profit = fixnet([i.net for i in items])

# item constraints
load = [i.mass for i in items]
invest = [i.buy for i in items]
constraints = ms.linear_symbolic(G=[load, invest], h=[_load, _spend])

# bounds (on x)
bounds = [(0, i.limit) for i in items]

# bounds constraints
lo = 'x%s >= %s'
lo = '\n'.join(lo % (i,str(float(j[0])).lstrip('0')) for (i,j) in enumerate(bounds))
hi = 'x%s <= %s'
hi = '\n'.join(hi % (i,str(float(j[1])).lstrip('0')) for (i,j) in enumerate(bounds))
constraints = '\n'.join([lo, hi]).strip() + '\n' + constraints
pf = ms.generate_penalty(ms.generate_conditions(ms.simplify(constraints)))

# integer constraints
cf = mc.integers(float)(lambda x:x)

# solve
mon = my.monitors.VerboseMonitor(10)
results = my.solvers.diffev2(lambda x: -profit(x), bounds, npop=400, bounds=bounds, ftol=1e-6, gtol=100, itermon=mon, disp=True, full_output=True, constraints=cf, penalty=pf)

print ('\nmax profit: %s' % -results[1])
print("load: %s <= %s" % (sum(i*j for i,j in zip(results[0], load)), _load))
print("spend: %s <= %s" % (sum(i*j for i,j in zip(results[0], invest)), _spend))
print('')

for item,quantity in enumerate(results[0]):
  print("item %d: %s" % (item, quantity))

with results:

dude@borel>$ python knapsack.py 
Generation 0 has ChiSquare: -57840.000000
Generation 10 has ChiSquare: -60216.000000
Generation 20 has ChiSquare: -60216.000000
Generation 30 has ChiSquare: -60216.000000
Generation 40 has ChiSquare: -60683.000000
Generation 50 has ChiSquare: -61029.000000
Generation 60 has ChiSquare: -64119.000000
Generation 70 has ChiSquare: -64119.000000
Generation 80 has ChiSquare: -64119.000000
Generation 90 has ChiSquare: -64119.000000
Generation 100 has ChiSquare: -64119.000000
Generation 110 has ChiSquare: -64119.000000
Generation 120 has ChiSquare: -64119.000000
Generation 130 has ChiSquare: -64119.000000
Generation 140 has ChiSquare: -64825.000000
Generation 150 has ChiSquare: -64825.000000
Generation 160 has ChiSquare: -64836.000000
Generation 170 has ChiSquare: -64836.000000
Generation 180 has ChiSquare: -64836.000000
Generation 190 has ChiSquare: -65072.000000
Generation 200 has ChiSquare: -65147.000000
Generation 210 has ChiSquare: -65147.000000
Generation 220 has ChiSquare: -65218.000000
Generation 230 has ChiSquare: -65218.000000
Generation 240 has ChiSquare: -65548.000000
Generation 250 has ChiSquare: -65548.000000
Generation 260 has ChiSquare: -65548.000000
Generation 270 has ChiSquare: -65548.000000
Generation 280 has ChiSquare: -65548.000000
Generation 290 has ChiSquare: -65677.000000
Generation 300 has ChiSquare: -65677.000000
Generation 310 has ChiSquare: -65703.000000
Generation 320 has ChiSquare: -65742.000000
Generation 330 has ChiSquare: -65742.000000
Generation 340 has ChiSquare: -65777.000000
Generation 350 has ChiSquare: -65856.000000
Generation 360 has ChiSquare: -65856.000000
Generation 370 has ChiSquare: -65857.000000
Generation 380 has ChiSquare: -65857.000000
Generation 390 has ChiSquare: -65877.000000
Generation 400 has ChiSquare: -65878.000000
Generation 410 has ChiSquare: -65878.000000
Generation 420 has ChiSquare: -65894.000000
Generation 430 has ChiSquare: -65897.000000
Generation 440 has ChiSquare: -65916.000000
Generation 450 has ChiSquare: -65929.000000
Generation 460 has ChiSquare: -65936.000000
Generation 470 has ChiSquare: -65979.000000
Generation 480 has ChiSquare: -65979.000000
Generation 490 has ChiSquare: -65979.000000
Generation 500 has ChiSquare: -65979.000000
Generation 510 has ChiSquare: -66018.000000
Generation 520 has ChiSquare: -66018.000000
Generation 530 has ChiSquare: -66020.000000
Generation 540 has ChiSquare: -66020.000000
Generation 550 has ChiSquare: -66025.000000
Generation 560 has ChiSquare: -66025.000000
Generation 570 has ChiSquare: -66025.000000
Generation 580 has ChiSquare: -66025.000000
Generation 590 has ChiSquare: -66025.000000
Generation 600 has ChiSquare: -66025.000000
Generation 610 has ChiSquare: -66027.000000
Generation 620 has ChiSquare: -66035.000000
Generation 630 has ChiSquare: -66035.000000
Generation 640 has ChiSquare: -66041.000000
Generation 650 has ChiSquare: -66051.000000
Generation 660 has ChiSquare: -66051.000000
Generation 670 has ChiSquare: -66056.000000
Generation 680 has ChiSquare: -66056.000000
Generation 690 has ChiSquare: -66056.000000
Generation 700 has ChiSquare: -66056.000000
Generation 710 has ChiSquare: -66056.000000
Generation 720 has ChiSquare: -66056.000000
Generation 730 has ChiSquare: -66056.000000
Generation 740 has ChiSquare: -66056.000000
Generation 750 has ChiSquare: -66056.000000
Generation 760 has ChiSquare: -66056.000000
STOP("ChangeOverGeneration with {'tolerance': 1e-06, 'generations': 100}")
Optimization terminated successfully.
         Current function value: -66056.000000
         Iterations: 763
         Function evaluations: 305600

max profit: 66056.0
load: 75000.0 <= 75000
spend: 315230.0 <= 350000

item 0: 300.0
item 1: 500.0
item 2: 0.0
item 3: 300.0
item 4: 200.0
item 5: 254.0
item 6: 100.0
item 7: 600.0
item 8: 0.0
item 9: 50.0

I'm going to assume that the solver is not tuned, so there's probably still some small wiggle room for improvement as the convergence at the end is still not super smooth -- however, I assume that the solution is close to optimal (based on the constraints check). I'd play with the settings and how the constraints/penalties are imposed a bit to see if the solution can improve a little bit more.