When trying to get fit, proper nutrition is essential. Among other things, one needs to watch the calory intake based on TDEE and activity, watch the macronutrients (macro calculator) while ensuring you’re getting sufficient micronutrients. In other words, a giant PITA to plan around.

Commercial tools like EatThisMuch.com are incredibly helpful in meeting the requirements above, adapting to your tastes. Its pricing is not accurate enough for what I can get locally, and it does not have all of the recipes I like. I’m sure I add my recipes to the tool, but why use their complicated solution when I can reinvent my overcomplicated solution?

The approach

With nutrition you get a calory budget. We’ll use 2000 kcal for our example. A common macro split I’ve seen is 34% proteins (4 kcal per 1g), 33% carbohydrates (4 kcal per 1g), 33% fats (9 kcal per 1g), which for our example would result in 170g protein, 165g carbs, 73g fats targets.

I want to minimize the price of the food as well, so our objective is to get as close as we can to the target calories and macros (e.g. within 5-10%) while minimizing the cost - a good objective for a constraint solver to take on. The following implementation will skip micronutrient calculations, just because I don’t have them precalculated for my recipes, but it should be easy to extend the program to support them in the future.

A while back I used Microsoft’s Z3 Theorem Prover’s constraint solving capabilities to solve the PuzLogic puzzle game. It worked well for this use case as well.

The optimization task can be written down as a system of equations:

total_price = meal[0].quantity * meal[0].price + meal[1].quantity * meal[1].price + ... + meal[n].quantity + meal[n].price
total_protein = meal[0].quantity * meal[0].protein + meal[1].quantity * meal[1].protein + ... + meal[n].quantity + meal[n].protein
total_carbs = meal[0].quantity * meal[0].carbs + meal[1].quantity * meal[1].carbs + ... + meal[n].quantity + meal[n].carbs
total_fats = meal[0].quantity * meal[0].fats + meal[1].quantity * meal[1].fats + ... + meal[n].quantity + meal[n].fats
lower_protein_limit <= total_protein <= upper_protein_limit
lower_carbs_limit <= total_carbs <= upper_carbs_limit
lower_fats_limit <= total_fats <= upper_fats_limit
minimize(total_price)

The equations are straightforward to translate for Z3.

Implementation

First we define the targets and some sample meals:

targets = {
    'calories': { 'ideal': 1997, 'lower_limit': 0.8, 'upper_limit': 1.2 },
    'protein': { 'ideal': 170, 'lower_limit': 0.95, 'upper_limit': 1.2 },
    'carbs': { 'ideal': 165, 'lower_limit': 0.9, 'upper_limit': 1.2 },
    'fat': { 'ideal': 73, 'lower_limit': 0.9, 'upper_limit': 1.2 },
    'price': 'minimize',
}

meals = {
    'cottage_cheese': { 'price': 1.4, 'calories': 483.2, 'fat': 11.59, 'carbs': 40.276, 'protein': 55.24, },
    # ... more recipes
    'tortilla_pizza': { 'price': 0.99, 'calories': 510.25, 'fat': 32.04, 'carbs': 26.71, 'protein': 26.66, },
}

Z3 provides convenient objects for Real variables, as well as Sum operations in its Python bindings. We’ll use the optimizer for this program. In its most basic form, we can write the optimizing program as:

from z3 import *

opt = Optimize()
total_price, total_calories, total_protein, total_carbs, total_fat = Reals('total_price total_calories total_protein total_carbs total_fat')
portions_cottage_cheese, portions_tortilla_pizza = Reals('portions_cottage_cheese portions_tortilla_pizza')

# Sum up price, protein, carb and fat totals between our meals
opt.add(total_price == Sum([portions_cottage_cheese*meals['cottage_cheese']['price'], portions_tortilla_pizza*meals['tortilla_pizza']['price']]))
opt.add(total_calories == Sum([portions_cottage_cheese*meals['cottage_cheese']['calories'], portions_tortilla_pizza*meals['tortilla_pizza']['calories']]))
opt.add(total_protein == Sum([portions_cottage_cheese*meals['cottage_cheese']['protein'], portions_tortilla_pizza*meals['tortilla_pizza']['protein']]))
opt.add(total_carbs == Sum([portions_cottage_cheese*meals['cottage_cheese']['carbs'], portions_tortilla_pizza*meals['tortilla_pizza']['carbs']]))
opt.add(total_fat == Sum([portions_cottage_cheese*meals['cottage_cheese']['fat'], portions_tortilla_pizza*meals['tortilla_pizza']['fat']]))

# Add limits that our meals have to sum up to
opt.add(And([
    total_calories >= targets['calories']['ideal'] * targets['calories']['lower_limit'],
    total_calories >= targets['calories']['ideal'] * targets['calories']['upper_limit']]))
opt.add(And([
    total_protein >= targets['protein']['ideal'] * targets['protein']['lower_limit'],
    total_protein >= targets['protein']['ideal'] * targets['protein']['upper_limit']]))
opt.add(And([
    total_carbs >= targets['carbs']['ideal'] * targets['carbs']['lower_limit'],
    total_carbs >= targets['carbs']['ideal'] * targets['carbs']['upper_limit']]))
opt.add(And([
    total_fat >= targets['fat']['ideal'] * targets['fat']['lower_limit'],
    total_fat >= targets['fat']['ideal'] * targets['fat']['upper_limit']]))

# Ensure that portions are non-zero values
opt.add(And([portions_cottage_cheese >= 0, portions_tortilla_pizza >= 0]))

# Let Z3 pick the meal portions while minimizing the price
opt.minimize(total_price)

if opt.check() == sat:
    m = opt.model()
    print(m)

Which outputs:

[portions_cottage_cheese = 150154650/49043707,
 portions_tortilla_pizza = 46250910/49043707,
 total_fat = 657/10,
 total_carbs = 297/2,
 total_protein = 47637960633/245218535,
 total_calories = 192308507415/98087414,
 total_price = 2560049109/490437070]

Or in a more readable format:

[portions_cottage_cheese = 3.061649683210121,
 portions_tortilla_pizza = 0.9430549366914699,
 total_fat = 65.7,
 total_carbs = 148.5,
 total_protein = 194.26737311272169,
 total_calories = 1960.5829083739532,
 total_price = 5.219933943818725]

I had to broaden the upper acceptable range of all of the macros to 1.2x for Z3 to find an acceptable solution with such few possible recipies to choose from, but it did find a satisfactory solution. About 3 portions of cottage cheese and about a portion of tortilla pizza gets us fairly close to the target TDEE and macro requirements.

However, if you run this program multiple times, the result will always be the same. That’s because Z3 only looks for the first satisfactory solution. If we want to instruct Z3 to look for more solutions, we have to invalidate the solutions it already found. Which is fairly easy, just need to change the model displaying part:

variants_to_show = 3
for i in range(variants_to_show):
    if opt.check() != sat:
        break
    if i >= 1:
        print('=' * 80)

    m = opt.model()
    for d in m.decls():
        print(d.name(), m[d].as_decimal(3))

    opt.add(Or([
        portions_cottage_cheese != m[portions_cottage_cheese],
        portions_tortilla_pizza != m[portions_tortilla_pizza]]))

Resulting in:

portions_cottage_cheese 3.061?
portions_tortilla_pizza 0.943?
total_fat 65.7
total_carbs 148.5
total_protein 194.267?
total_calories, 1960.582?
total_price 5.219?
================================================================================
portions_cottage_cheese 2.561?
portions_tortilla_pizza 1.697?
total_fat 84.061?
total_carbs 148.5
total_protein 186.747?
total_calories, 2103.685?
total_price 5.266?
================================================================================
portions_cottage_cheese 3.275?
portions_tortilla_pizza 0.865?
total_fat 65.7
total_carbs 155.034?
total_protein 204
total_calories, 2024.325?
total_price 5.442?

Now we can get Z3 to output a bunch of meal plans for us to choose from.

There are some issues with the code above: it does not scale well with more meals added, the meal portions are calculated in fractions while we’d rather work with full portions, the code could be further cleaned up. But the idea is still there. I ended up with this:

opt = Optimize()

expressions = dict((name, []) for name in targets.keys())
variables = {}

for meal_name, meal in meals.items():
    portions = variables[f'portions_{meal_name}'] = Real(f'portions_{meal_name}')
    opt.add(Or([portions == option for option in [0, 1, 2]])) # Limit to up to 2 full portions per meal

    for target_name in targets.keys():
        expressions[target_name].append(portions * meal[target_name])

expr_to_minimize = None
for name, value in targets.items():
    total = variables[f'total_{name}'] = Real(f'total_{name}')

    if value == 'minimize':
        opt.add(total == Sum(expressions[name]))
        expr_to_minimize = total
    else:
        opt.add(total == Sum(expressions[name]))
        opt.add(total >= value['ideal'] * value['lower_limit'])
        opt.add(total <= value['ideal'] * value['upper_limit'])

optimized = opt.minimize(expr_to_minimize)
variants_to_show = 5
for i in range(variants_to_show):
    if opt.check() != sat:
        break
    if i >= 1:
        print('=' * 80)

    m = opt.model()
    for name, variable in variables.items():
        if name.startswith('portions_') and m[variable] == 0:
            continue

        if m[variable].is_real():
            print(name, m[variable].as_decimal(3))
        elif m[variable].is_int():
            print(name, m[variable].as_long())
        else:
            print(name, m[variable])

    # Ignore already found variants
    opt.add(Or([variable != m[variable] for name, variable in variables.items() if name.startswith('portions_')]))

Results

With a bigger set of recipes to choose from, Z3 came up with the following meal plans:

portions_peanuts 2
portions_protein_pancakes 2
portions_chicken_quesadilla_w_cheese 1
total_calories 1994.24
total_protein 169.76
total_carbs 164.76
total_fat 77
total_price 3
================================================================================
portions_peanuts 1
portions_protein_pancakes 2
portions_pulled_pork_quesadilla 2
total_calories 2088.18
total_protein 168.39
total_carbs 180.74
total_fat 79.38
total_price 3.11
================================================================================
portions_chicken_pilaf 1
portions_peanuts 2
portions_protein_pancakes 1
portions_chicken_quesadilla_w_cheese 1
total_calories 2051.22
total_protein 168.23
total_carbs 172.24
total_fat 78.6
total_price 3.46
================================================================================
portions_chicken_quesadilla 1
portions_protein_pancakes 2
portions_tortilla_pizza 1
total_calories 1917.45
total_protein 162.78
total_carbs 171.85
total_fat 66.5
total_price 3.61
================================================================================
portions_almonds 1
portions_peanuts 1
portions_protein_pancakes 2
portions_chicken_quesadilla_w_cheese 1
total_calories 1979.74
total_protein 166.91
total_carbs 162.91
total_fat 80
total_price 3.77

And with a bigger bank of recipes it beats me having to plan it all out. Now I can let it generate hundreds of plans for me to choose from, and they’ll all be composed of recipes I like and based on prices at my grocery store. Granted, some of the suggestions are not great: not doable for time reasons or personal tastes, lack of micronutrients (e.g. the plans should definitely include more veggies), but it’s a decent starting point. With a selection of better recipes the meal plans would be better too.

You can play around with the full implementation on Google’s Colaboratory.