PULPO Optimization Framework: Complete Demonstration#

This notebook provides a comprehensive demonstration of PULPO’s capabilities for sustainable technology optimization, covering basic single-objective optimization through advanced uncertainty analysis and multi-objective trade-off exploration.

1. System Description#

The demonstration system models interconnected background and foreground processes for low-carbon methanol production:

Background System supplies conventional commodities:

  • Natural gas extraction and electricity generation

  • Wind electricity generation

  • Water supply

  • Hydrogen production via steam methane reforming (SMR)

  • Oxygen production via air separation unit (ASU)

Foreground System represents novel technologies:

  • Hydrogen electrolysis (producing H₂ and O₂ as byproduct)

  • Direct air capture (DAC) for CO₂ supply

  • Methanol synthesis using captured CO₂ and H₂

  • Ozone production consuming oxygen

Multifunctional Links: The O₂-byproduct and O₂-market activities balance oxygen flows between electrolysis byproduct recovery and conventional supply, ensuring proper system allocation.

Environmental Flows: The system tracks emissions of CO₂, CH₄, particulate matter (PM), and water use to the biosphere.

This represents an integrated Power-to-Methanol pathway with linked energy and material exchanges between renewable and fossil-based subsystems.

Flowchart of the system

1.1 Technology Choices and System Interdependencies#

The system incorporates several critical technological choices that determine environmental performance:

Energy Supply Alternatives:

  • Electricity generation: Natural gas combined-cycle (NGCC) vs. wind power

  • Represents background system choices reflecting alternative energy mixes

Hydrogen Production Pathways:

  • Steam methane reforming (SMR) in background system

  • Water electrolysis in foreground system

  • Demonstrates foreground-background substitution dynamics

Oxygen Supply Options:

  • Conventional production via air separation unit (ASU) in background

  • Recovery as electrolysis byproduct in foreground

  • Market-based supply balancing both sources

The O₂-market activity ensures that multiple oxygen sources can coexist or compete depending on system boundaries and optimization objectives.

1.2 Technosphere Matrix (A)#

The technosphere matrix represents product flows between processes. Positive values indicate outputs, negative values indicate inputs.

Product ↓ \\ Process → NG extraction NG electricity Wind electricity Water supply H₂ SMR O₂ ASU H₂ electrolysis O₂-byproduct O₂-market DAC plant Methanol synthesis Ozone production
Natural gas [kg] +1 −0.19 0 0 −3.2 0 0 0 0 0 0 0
Electricity, NG [kWh] 0 +1 0 −0.49 −0.10 −0.22 0 0 0 0 −1.2 0
Electricity, wind [kWh] 0 0 +1 0 0 0 −50.0 0 0 −1.5 0 −0.5
Water [m³] 0 0 0 +1 −0.003 0 −0.009 0 0 0 0 0
Hydrogen (SMR) [kg] 0 0 0 0 +1 0 0 0 0 0 0 0
Oxygen (ASU) [kg] 0 0 0 0 0 +1 0 0 0 0 0 −1.0
Hydrogen (electrolysis) [kg] 0 0 0 0 0 0 +1 0 0 0 −0.19 0
Oxygen (byproduct) [kg] 0 0 0 0 0 0 +8.0 -1 −1 0 0 0
Oxygen (market) [kg] 0 0 0 0 0 0 0 0 +1 0 0 0
CO₂, captured [kg] 0 0 0 0 0 0 0 0 0 +1 −1.375 0
Methanol [kg] 0 0 0 0 0 0 0 0 0 0 +1 0
Ozone [kg] 0 0 0 0 0 0 0 0 0 0 0 +1

1.3 Biosphere Matrix (B)#

The biosphere matrix represents environmental flows from processes to the environment. All values represent emissions or resource consumption.

Flow ↓ ---- \\ ---- Process → NG extraction NG electricity Wind electricity Water supply H₂ SMR O₂ ASU H₂ electrolysis O₂-byproduct O₂-market DAC plant Methanol synthesis Ozone production
CO₂ emission [kg] +0.0 +0.44 +0.012 0 +9.0 0 +0.1 0 0 −1.0 +0.2 0
CH₄ emission [kg] +0.0074 +0.0000028 0 0 0 0 0 0 0 0 0 0
PM emission [g] 0 +0.015 0 0 0 0 0 0 0 0 +0.001 +0.002
Water use [m³] 0 +0.0018 0 +1.0 0 0 0 0 0 0 0 0

2. Implementation: Basic Setup#

2.1 Database Initialization#

Initialize the sample database containing all processes and flows described in the system matrices. This creates both background and foreground databases with the complete methanol production system.

from pulpo import pulpo
pulpo.setup_sample_db()
biosphere3 created with 4 flows.
technosphere database created with 5 activities.
background_db created with 6 activities and uncertainty added to exchanges.
foreground_db created with 6 activities and uncertainty added to exchanges.
Registered 3 LCIA methods.

2.2 Library Import#

Import required libraries for optimization analysis and data manipulation.

import sys
import numpy as np
import pandas as pd
import copy

2.3 PULPO Configuration#

Configure the PULPO optimizer with:

  • Project name: Identifier for optimization project

  • Databases: Background (conventional) and foreground (novel) systems

  • Methods: Impact assessment categories with relative weights

    • Climate change: weight = 1 (primary objective)

    • Resources: weight = 0 (inactive initially)

from pulpo.utils.utils import is_bw25

project = "sample_project_bw25" if is_bw25() else "sample_project"
database = ["background_db", "foreground_db"]
method = {"('my project', 'climate change')": 1,
          "('my project', 'resources')": 0,}

2.4 Optimizer Initialization#

Create the PULPO optimizer instance.

pulpo_worker = pulpo.PulpoOptimizer(project, database, method)

⚠️ Attention:
PULPO assumes that the biosphere in your project is named biosphere3.
If your project uses a different biosphere name, make sure to explicitly provide it to your PULPO worker.

For example:

pulpo_worker.intervention_matrix = "ecoinvent-3.10-biosphere"

2.5 Life Cycle Inventory Data Loading#

Load technosphere matrix (A), biosphere matrix (B), and process/intervention mappings. This creates the mathematical foundation for optimization by establishing the linear system: Ax = f and environmental impacts = Bx.

pulpo_worker.get_lci_data()

3. Single-Objective Optimization Examples#

3.1 Functional Unit Definition#

The functional unit defines the reference output for environmental impact assessment. We begin with 1 kg methanol as the primary product demand.

methanol_process = pulpo_worker.retrieve_processes(reference_products='methanol')

demand = {methanol_process[0]: 1}

3.2 Technology Choice Definition#

Define alternative technologies for key system components:

Electricity Generation:

  • Wind electricity (renewable)

  • Natural gas electricity (fossil-based)

Hydrogen Production:

  • Steam methane reforming (SMR) - conventional background process

  • Water electrolysis - novel foreground process

Oxygen Supply:

  • Air separation unit (ASU) - conventional production

  • Market supply - generic oxygen source

  • Electrolysis byproduct - recovered from hydrogen production

These alternatives represent the core technological decisions determining system environmental performance.

## Electricity Choice
electricity_technologies = ["wind electricity",
                            "natural gas electricity",]
electricity_processes = pulpo_worker.retrieve_processes(processes=electricity_technologies)

## Hydrogen Choice
hydrogen_technologies = ["hydrogen SMR",
                         "hydrogen electrolysis",]
hydrogen_processes = pulpo_worker.retrieve_processes(processes=hydrogen_technologies)

## Oxygen Choice
oxygen_technologies = ["O2-market",
                       "O2 ASU",]
oxygen_processes = pulpo_worker.retrieve_processes(processes=oxygen_technologies)

### Oxygen byproduct
oxygen_aux = "O2-byproduct"
oxygen_byproduct = pulpo_worker.retrieve_processes(processes=[oxygen_aux])

3.3 Choices Dictionary Configuration#

Organize technology alternatives into an optimization-ready structure. Each choice category contains multiple alternatives with capacity constraints.

⚠️ Attention:
A capacity of float('inf') means no limit: the solver is not restricted by it.
Give an alternative a finite capacity only where it has a real limit. A huge finite number such as 1e10 adds
nothing and can degrade the solver’s numerics; leave it infinite instead.

## Choices
choices = {
    "Electricity": {electricity_processes[0]: float('inf'),
                    electricity_processes[1]: float('inf')},
    "Hydrogen": {hydrogen_processes[0]: float('inf'),
                 hydrogen_processes[1]: float('inf')},
    "Oxygen": {oxygen_processes[0]: float('inf'), 
               oxygen_processes[1]: float('inf'),},
}

3.4 Lower Bound Constraints#

Configure process flow constraints. The O₂-byproduct process requires non-negative flow (≥ 0) to ensure physically meaningful byproduct recovery from electrolysis.

# Lower bound
lower_bound = {oxygen_byproduct[0]: 0}

⚠️ Attention:
With this configuration, the O₂ by-product is allowed to be vented (e.g., to a market or directly to the environment). If you want to force the system to internally utilize the by-product, you must also specify an upper bound of 0, which prevents any venting and requires full integration within the production system.

3.5 First Optimization: Single Product#

Solve optimization with 1 kg methanol demand. PULPO identifies the technology combination minimizing climate change impact while satisfying demand constraints.

pulpo_worker.instantiate(choices=choices, demand=demand, lower_limit=lower_bound)
results = pulpo_worker.solve()
Creating Instance
Instance created
optimal solution found:  1.7471600582975855
Optimization problem solved using Highspy

Results Analysis: Technology Selection#

The optimization results show selected technologies and resource utilization. Observe which electricity and hydrogen sources were chosen to minimize climate impact.

pulpo_worker.summarize_results(zeroes=True)

Total Impact(s)

Weight Value
Method
('my project', 'climate change') 1 1.74716
('my project', 'resources') 0 0.00171

Choices Made

Electricity

Value Capacity
Metadata
wind electricity | electricity | GLO 12.763338 inf

Hydrogen

Value Capacity
Metadata
hydrogen electrolysis | hydrogen | GLO 0.19 inf

Oxygen

Value Capacity
Metadata

Constraints

Constraints Lower

Key Metadata Value Limit
ID
7 (foreground_db, O2-byproduct) O2-byproduct | oxygen | GLO 1.52 0

3.6 Multi-Product Demand#

Extend the system to include ozone production (2 kg) alongside methanol. This demonstrates PULPO’s handling of multi-functional systems and oxygen flow optimization between different end uses.

ozone_process = pulpo_worker.retrieve_processes(reference_products='ozone')

demand = {methanol_process[0]: 1,
          ozone_process[0]: 2,}

Optimization: Multi-Product System#

Solve with combined methanol and ozone demands. Oxygen byproduct from hydrogen electrolysis can now be utilized for ozone production, enabling potential system synergies.

pulpo_worker.instantiate(choices=choices, demand=demand, lower_limit=lower_bound)
results = pulpo_worker.solve()
Creating Instance
Instance created
optimal solution found:  1.760427258456396
Optimization problem solved using Highspy

Results Analysis: Multi-Product System#

Compare results with single-product case. Examine how oxygen flows changed and whether optimal technology choices shifted due to ozone demand integration.

pulpo_worker.summarize_results(zeroes=True)

Total Impact(s)

Weight Value
Method
('my project', 'climate change') 1 1.760427
('my project', 'resources') 0 0.001710

Choices Made

Electricity

Value Capacity
Metadata
wind electricity | electricity | GLO 13.868938 inf

Hydrogen

Value Capacity
Metadata
hydrogen electrolysis | hydrogen | GLO 0.19 inf

Oxygen

Value Capacity
Metadata
O2-market | oxygen | GLO 1.52 inf
O2 ASU | oxygen | GLO 0.48 inf

Constraints

'No constraint data to display.'

4. Constrained Optimization#

4.1 Technology Capacity Constraints#

Real-world technologies operate within finite capacity limits. Apply realistic constraints:

  • Electricity sources: 8 units maximum each

  • Hydrogen production: 0.1 units maximum each

  • Oxygen supply: 1.2 units maximum each

These constraints force optimization within technological feasibility bounds.

## Choices
choices = {
    "Electricity": {electricity_processes[0]: 8, electricity_processes[1]: 8},
    "Hydrogen": {hydrogen_processes[0]: 0.1, hydrogen_processes[1]: 0.1},
    "Oxygen": {oxygen_processes[0]: 1.2, oxygen_processes[1]:1.2,},
}

Optimization: Capacity-Constrained System#

Solve with active capacity constraints. This may require technology portfolio diversification rather than single-technology selection for environmental optimality.

pulpo_worker.instantiate(choices=choices, demand=demand, lower_limit=lower_bound)
results = pulpo_worker.solve()
Creating Instance
Instance created
optimal solution found:  3.2944402907815795
Optimization problem solved using Highspy

Results Analysis: Impact of Capacity Constraints#

Examine how technology mix and total environmental impact change under capacity limitations. The optimizer balances environmental optimality against technological feasibility.

pulpo_worker.summarize_results(zeroes=True)

Total Impact(s)

Weight Value
Method
('my project', 'climate change') 1 3.294440
('my project', 'resources') 0 0.003935

Choices Made

Electricity

Value Capacity
Metadata
wind electricity | electricity | GLO 8.000000 8
natural gas electricity | electricity | GLO 1.536073 8

Hydrogen

Value Capacity
Metadata
hydrogen electrolysis | hydrogen | GLO 0.10 0.1
hydrogen SMR | hydrogen | GLO 0.09 0.1

Oxygen

Value Capacity
Metadata
O2 ASU | oxygen | GLO 1.2 1.2
O2-market | oxygen | GLO 0.8 1.2

Constraints

'No constraint data to display.'

4.2 Resource Availability Constraints#

Beyond technology capacity, resource availability represents another critical limitation. Constrain water supply to 1 m³ to simulate water scarcity scenarios.

water_process = pulpo_worker.retrieve_processes(reference_products='water')
upper_bound = {water_process[0]: 1}

Optimization: Combined Constraints#

Apply both technology capacity limits and water availability constraints simultaneously. This represents realistic scenarios with multiple operational limitations.

pulpo_worker.instantiate(choices=choices, demand=demand, lower_limit=lower_bound, upper_limit=upper_bound)
results = pulpo_worker.solve()
Creating Instance
Instance created
optimal solution found:  3.2944402907815795
Optimization problem solved using Highspy

Results Analysis: Resource Limitation Effects#

Water limitation forces the system toward less water-intensive processes. Compare water consumption across scenarios to quantify constraint effects on technology selection.

pulpo_worker.summarize_results(zeroes=True)

Total Impact(s)

Weight Value
Method
('my project', 'climate change') 1 3.294440
('my project', 'resources') 0 0.003935

Choices Made

Electricity

Value Capacity
Metadata
wind electricity | electricity | GLO 8.000000 8
natural gas electricity | electricity | GLO 1.536073 8

Hydrogen

Value Capacity
Metadata
hydrogen electrolysis | hydrogen | GLO 0.10 0.1
hydrogen SMR | hydrogen | GLO 0.09 0.1

Oxygen

Value Capacity
Metadata
O2 ASU | oxygen | GLO 1.2 1.2
O2-market | oxygen | GLO 0.8 1.2

Constraints

Constraints Upper

Key Metadata Value Limit
ID
3 (background_db, water supply) water supply | water | GLO 0.00117 1

4.3 Impact Category Switching#

Demonstrate objective function flexibility by switching from climate change to resource use optimization. Different environmental priorities lead to distinct technology choices.

pulpo_worker.method = {"('my project', 'climate change')": 0,
          "('my project', 'resources')": 1,}

Optimization: Resource Efficiency Focus#

Reset constraints to unconstrained levels and optimize for resource efficiency. This will likely select different technologies than climate-optimized solutions.

## Choices
choices = {
    "Electricity": {electricity_processes[0]: float('inf'), electricity_processes[1]: float('inf')},
    "Hydrogen": {hydrogen_processes[0]: float('inf'), hydrogen_processes[1]: float('inf')},
    "Oxygen": {oxygen_processes[0]: float('inf'), oxygen_processes[1]: float('inf'),},
}
pulpo_worker.instantiate(choices=choices, demand=demand, lower_limit=lower_bound)
results = pulpo_worker.solve()
Creating Instance
Instance created
optimal solution found:  0.0005699999978020787
Optimization problem solved using Highspy

Results Analysis#

Compare resource-optimized results with climate change optimization. Different environmental priorities (climate vs. resources) yield different optimal technology portfolios, demonstrating the importance of multi-criteria decision-making in sustainable technology design.

pulpo_worker.summarize_results(zeroes=True)

Total Impact(s)

Weight Value
Method
('my project', 'resources') 1 0.000570
('my project', 'climate change') 0 3.475288

Choices Made

Electricity

Value Capacity
Metadata
wind electricity | electricity | GLO 4.721779 inf

Hydrogen

Value Capacity
Metadata
hydrogen SMR | hydrogen | GLO 0.19 inf

Oxygen

Value Capacity
Metadata
O2 ASU | oxygen | GLO 2.0 inf

Constraints

'No constraint data to display.'

5. Summary: Core PULPO Capabilities#

This demonstration covered PULPO’s fundamental optimization features:

System Integration:

  • Multi-database integration (background + foreground systems)

  • Technology choice optimization across energy, hydrogen, and oxygen sources

  • Multi-product system optimization with material flow synergies

Constraint Handling:

  • Technology capacity constraints (realistic deployment limits)

  • Resource availability constraints (scarcity scenarios)

  • Process flow constraints (physical feasibility)

Multi-Criteria Analysis:

  • Objective function flexibility (climate change vs. resource use)

  • Technology trade-off analysis under different environmental priorities

  • Impact assessment across multiple environmental categories

6. Multi-Objective Optimization: Pareto Front Analysis#

Explore trade-offs between competing objectives through Pareto front generation. This analysis reveals how environmental performance varies when balancing climate change and resource use priorities using the epsilon-constraint method.

6.1 Multi-Objective Setup#

Create dual-objective PULPO worker with equal weighting for both climate change and resource use objectives to establish baseline performance ranges.

# Create multi-objective method configuration
methods_multi = {
    "('my project', 'climate change')": 1,
    "('my project', 'resources')": 1
}

# Initialize multi-objective PULPO worker
pulpo_worker_multi = pulpo.PulpoOptimizer(project, database, methods_multi)
pulpo_worker_multi.get_lci_data()

print("Multi-objective PULPO worker created with both climate change and resource use objectives")
Multi-objective PULPO worker created with both climate change and resource use objectives

6.2 Extreme Point Identification#

Determine boundary solutions for Pareto front generation:

  1. Climate-optimal solution: Minimize climate change (resource use unconstrained)

  2. Resource-optimal solution: Minimize resource use (climate change unconstrained)

These extreme points define the feasible range for epsilon constraints.

# 1. Find climate-optimal solution (minimize climate change)
methods_climate = {
    "('my project', 'climate change')": 1,
    "('my project', 'resources')": 0
}

pulpo_worker_climate = pulpo.PulpoOptimizer(project, database, methods_climate)
pulpo_worker_climate.get_lci_data()

# Use unconstrained choices for extreme points
choices_unconstrained = {
    "Electricity": {electricity_processes[0]: float('inf'), electricity_processes[1]: float('inf')},
    "Hydrogen": {hydrogen_processes[0]: float('inf'), hydrogen_processes[1]: float('inf')},
    "Oxygen": {oxygen_processes[0]: float('inf'), oxygen_processes[1]: float('inf')},
}

pulpo_worker_climate.instantiate(choices=choices_unconstrained, demand=demand, lower_limit=lower_bound)
results_climate = pulpo_worker_climate.solve()

# Extract optimal impacts
climate_optimal_cc = pulpo_worker_climate.instance.impacts_calculated["('my project', 'climate change')"].value
climate_optimal_res = pulpo_worker_climate.instance.impacts_calculated["('my project', 'resources')"].value

print(f"Climate-optimal solution:")
print(f"  Climate change: {climate_optimal_cc:.4f}")
print(f"  Resource use: {climate_optimal_res:.4f}")
Creating Instance
Instance created
optimal solution found:  1.760427258456396
Optimization problem solved using Highspy
Climate-optimal solution:
  Climate change: 1.7604
  Resource use: 0.0017
# 2. Find resource-optimal solution (minimize resource use)
methods_resource = {
    "('my project', 'climate change')": 0,
    "('my project', 'resources')": 1
}

pulpo_worker_resource = pulpo.PulpoOptimizer(project, database, methods_resource)
pulpo_worker_resource.get_lci_data()

pulpo_worker_resource.instantiate(choices=choices_unconstrained, demand=demand, lower_limit=lower_bound)
results_resource = pulpo_worker_resource.solve()

# Extract optimal impacts
resource_optimal_cc = pulpo_worker_resource.instance.impacts_calculated["('my project', 'climate change')"].value
resource_optimal_res = pulpo_worker_resource.instance.impacts_calculated["('my project', 'resources')"].value

print(f"Resource-optimal solution:")
print(f"  Climate change: {resource_optimal_cc:.4f}")
print(f"  Resource use: {resource_optimal_res:.4f}")

print(f"\nObjective ranges:")
print(f"  Climate change: {climate_optimal_cc:.4f} to {resource_optimal_cc:.4f}")
print(f"  Resource use: {resource_optimal_res:.4f} to {climate_optimal_res:.4f}")
Creating Instance
Instance created
optimal solution found:  0.0005699999978020787
Optimization problem solved using Highspy
Resource-optimal solution:
  Climate change: 3.4753
  Resource use: 0.0006

Objective ranges:
  Climate change: 1.7604 to 3.4753
  Resource use: 0.0006 to 0.0017

6.3 Epsilon-Constraint Implementation#

Generate Pareto front using systematic constraint variation:

  • Primary objective: Minimize resource use

  • Secondary objective: Climate change (constrained via upper impact limits)

  • Method: Vary climate change constraint from climate-optimal to resource-optimal values

# Generate epsilon constraints for climate change
n_points = 8  # Number of Pareto points
epsilon_list = np.linspace(climate_optimal_cc, resource_optimal_cc, num=n_points)

# Initialize empty list to store Pareto front data
pareto_front = []

# Category names for easier reference
cc_category = "('my project', 'climate change')"
res_category = "('my project', 'resources')"

print(f"Generating Pareto front with {n_points} points...")
print(f"Climate change constraint range: {climate_optimal_cc:.4f} to {resource_optimal_cc:.4f}")
print("=" * 60)

# Configure worker for epsilon-constraint (minimize resources, constrain climate)
methods_epsilon = {
    "('my project', 'climate change')": 0,
    "('my project', 'resources')": 1
}

pulpo_worker_pareto = pulpo.PulpoOptimizer(project, database, methods_epsilon)
pulpo_worker_pareto.get_lci_data()

for i, epsilon in enumerate(epsilon_list):
    print(f"Point {i+1}/{n_points}: Climate constraint = {epsilon:.4f}")
    
    # Set up epsilon constraint
    epsilon_constraint = {cc_category: epsilon}
    
    # Solve with climate change constraint
    pulpo_worker_pareto.instantiate(
        choices=choices_unconstrained, 
        demand=demand, 
        lower_limit=lower_bound,
        upper_imp_limit=epsilon_constraint
    )
    results = pulpo_worker_pareto.solve()
    

    impacts = {x: pulpo_worker_pareto.instance.impacts_calculated[x].value for x in pulpo_worker_pareto.instance.impacts_calculated}
    
    # Store Pareto point (resource use, climate change)
    pareto_point = (impacts[res_category], impacts[cc_category])
    pareto_front.append(pareto_point)
    
    print(f"  → Resource use: {impacts[res_category]:.4f}")
    print(f"  → Climate change: {impacts[cc_category]:.4f}")
    print()

print(f"Pareto front generation complete: {len(pareto_front)} feasible points")
Generating Pareto front with 8 points...
Climate change constraint range: 1.7604 to 3.4753
============================================================
Point 1/8: Climate constraint = 1.7604
Creating Instance
Instance created
optimal solution found:  0.0017099999049305925
Optimization problem solved using Highspy
  → Resource use: 0.0017
  → Climate change: 1.7604

Point 2/8: Climate constraint = 2.0054
Creating Instance
Instance created
optimal solution found:  0.001547142775340805
Optimization problem solved using Highspy
  → Resource use: 0.0015
  → Climate change: 2.0054

Point 3/8: Climate constraint = 2.2504
Creating Instance
Instance created
optimal solution found:  0.0013842856457510172
Optimization problem solved using Highspy
  → Resource use: 0.0014
  → Climate change: 2.2504

Point 4/8: Climate constraint = 2.4954
Creating Instance
Instance created
optimal solution found:  0.0012214285161612295
Optimization problem solved using Highspy
  → Resource use: 0.0012
  → Climate change: 2.4954

Point 5/8: Climate constraint = 2.7403
Creating Instance
Instance created
optimal solution found:  0.0010585713865714416
Optimization problem solved using Highspy
  → Resource use: 0.0011
  → Climate change: 2.7403

Point 6/8: Climate constraint = 2.9853
Creating Instance
Instance created
optimal solution found:  0.0008957142569816538
Optimization problem solved using Highspy
  → Resource use: 0.0009
  → Climate change: 2.9853

Point 7/8: Climate constraint = 3.2303
Creating Instance
Instance created
optimal solution found:  0.0007328571273918662
Optimization problem solved using Highspy
  → Resource use: 0.0007
  → Climate change: 3.2303

Point 8/8: Climate constraint = 3.4753
Creating Instance
Instance created
optimal solution found:  0.0005699999978020787
Optimization problem solved using Highspy
  → Resource use: 0.0006
  → Climate change: 3.4753

Pareto front generation complete: 8 feasible points

6.4 Pareto Front Visualization#

Generate comprehensive visualization showing trade-offs between climate change and resource use across the complete feasible solution space.

import matplotlib.pyplot as plt

# Extract data for plotting
if pareto_front:
    resource_values = [point[0] for point in pareto_front]
    climate_values = [point[1] for point in pareto_front]
    
    # Create the plot
    plt.figure(figsize=(10, 6))
    
    # Plot Pareto front
    plt.plot(climate_values, resource_values, 'bo-', linewidth=2, markersize=8, 
             label='Pareto Front', alpha=0.7)
    
    # Highlight extreme points
    plt.plot(climate_optimal_cc, climate_optimal_res, 'ro', markersize=12, 
             label='Climate Optimal', alpha=0.8)
    plt.plot(resource_optimal_cc, resource_optimal_res, 'go', markersize=12, 
             label='Resource Optimal', alpha=0.8)
    
    # Formatting
    plt.xlabel('Climate Change Impact', fontsize=12)
    plt.ylabel('Resource Use Impact', fontsize=12)
    plt.title('Pareto Front: Climate Change vs. Resource Use\n(Multi-Product Methanol + Ozone System)', 
              fontsize=14, fontweight='bold')
    plt.grid(True, alpha=0.3)
    plt.legend()
    
    plt.tight_layout()
    plt.show()
    
    # Print summary statistics
    print("📈 PARETO FRONT ANALYSIS")
    print("=" * 40)
    print(f"Number of Pareto points: {len(pareto_front)}")
    print(f"\nClimate Change Range:")
    print(f"  Best: {min(climate_values):.4f}")
    print(f"  Worst: {max(climate_values):.4f}")
    print(f"  Range: {max(climate_values) - min(climate_values):.4f}")
    print(f"\nResource Use Range:")
    print(f"  Best: {min(resource_values):.4f}")
    print(f"  Worst: {max(resource_values):.4f}")
    print(f"  Range: {max(resource_values) - min(resource_values):.4f}")
    
else:
    print("No feasible Pareto points found")
../../_images/dd5efaca9ddcf9ac4641e603a459466f899a71f3ea09ab301102644dd1aad014.png
📈 PARETO FRONT ANALYSIS
========================================
Number of Pareto points: 8

Climate Change Range:
  Best: 1.7604
  Worst: 3.4753
  Range: 1.7149

Resource Use Range:
  Best: 0.0006
  Worst: 0.0017
  Range: 0.0011

7. Uncertainty Analysis: Monte Carlo Simulation#

Real-world systems contain parameter uncertainties that propagate through optimization results. PULPO’s Monte Carlo capabilities enable robustness assessment and uncertainty quantification for environmental impact optimization.

Analysis Objectives:

  • Quantify parameter uncertainty effects on optimal solutions

  • Assess result robustness across uncertainty ranges

  • Analyze correlation between environmental objectives under uncertainty

  • Give 90% intervals for the optimal impacts

7.1 Uncertainty Worker Configuration#

Configure specialized PULPO worker for uncertainty analysis with climate change minimization as primary objective while tracking both impact categories.

# Create uncertainty analysis worker focusing on climate change
methods_uncertainty = {
    "('my project', 'climate change')": 1,
    "('my project', 'resources')": 0,  # weight 0: not optimized, but reported in the results
}

# Initialize uncertainty PULPO worker
pulpo_worker_unc = pulpo.PulpoOptimizer(project, database, methods_uncertainty)
pulpo_worker_unc.get_lci_data()

print("Uncertainty analysis PULPO worker created successfully")
print("Objective: Minimize climate change impact under parameter uncertainty")
Uncertainty analysis PULPO worker created successfully
Objective: Minimize climate change impact under parameter uncertainty

7.2 Monte Carlo Parameter Configuration#

Configure system parameters for uncertainty analysis:

  • Demand specification: Multi-product system (1 kg methanol + 2 kg ozone)

  • Technology options: All electricity, hydrogen, and oxygen alternatives available

  • Constraints: Unconstrained capacities for computational clarity

  • Uncertainty source: Parameter sampling from database uncertainty distributions

# Configure same system as before for uncertainty analysis
# Multi-product demand: methanol + ozone
demand_unc = {methanol_process[0]: 1, ozone_process[0]: 2}

# Technology choices (unconstrained for clarity)
choices_unc = {
    "Electricity": {electricity_processes[0]: float('inf'), electricity_processes[1]: float('inf')},
    "Hydrogen": {hydrogen_processes[0]: float('inf'), hydrogen_processes[1]: float('inf')},
    "Oxygen": {oxygen_processes[0]: float('inf'), oxygen_processes[1]: float('inf')},
}

# Lower bounds
lower_bound_unc = {oxygen_byproduct[0]: 0}

print("System configuration for uncertainty analysis:")
print(f"  Demand: {len(demand_unc)} products")
print(f"  Technology choices: {len(choices_unc)} categories")
print(f"  Lower bounds: {len(lower_bound_unc)} constraints")
System configuration for uncertainty analysis:
  Demand: 2 products
  Technology choices: 3 categories
  Lower bounds: 1 constraints

7.3 Monte Carlo Execution#

Execute Monte Carlo simulation using solve_MC function:

  1. Problem instantiation with specified configuration

  2. Parameter sampling across multiple iterations from uncertainty distributions

  3. Optimization solving for each parameter realization

  4. Result collection across all Monte Carlo iterations

The function automatically handles parameter uncertainty sampling from database distributions.

# Run Monte Carlo simulation
print("🎲 Starting Monte Carlo simulation...")
print("=" * 50)

# Number of Monte Carlo iterations
n_iterations = 500

# Run Monte Carlo analysis
pulpo_worker_unc.instantiate(choices=choices_unc, demand=demand_unc, lower_limit=lower_bound_unc)
mc_results = pulpo_worker_unc.solve_MC(
    n_it=n_iterations,
    seed=1,  # reproducible samples
)

print(f"✅ Monte Carlo simulation completed successfully!")
print(f"   Total iterations: {n_iterations}")
print(f"   Results structure: dict with {len(mc_results)} entries")
print(f"   Available result keys per iteration: {list(next(iter(mc_results.values())).keys())}")
🎲 Starting Monte Carlo simulation...
==================================================
Creating Instance
Instance created
Pre-sampling 500 LCI matrix sets...
Solving 500 Monte Carlo optimizations in parallel...
Running 500 Monte Carlo optimizations in parallel (n_jobs=-1)...
✅ Monte Carlo simulation completed successfully!
   Total iterations: 500
   Results structure: dict with 500 entries
   Available result keys per iteration: ['Scaling Vector', 'Intervention Vector', 'Slack', 'Impacts', 'Transgressions', 'Demand', 'Choices', 'Constraints Upper', 'Constraints Lower', 'Constraints Upper Elem']

7.4 Statistical Analysis#

Collect the optimal impacts of every sample and summarize them: mean, standard deviation, 90% interval (5th to 95th percentile) and the correlation between the two impact categories.

mc_results[0]['Impacts']
Weight Value
Method
('my project', 'climate change') 1 1.633252
('my project', 'resources') 0 0.001787
# Extract both climate change and resource impacts from Monte Carlo results

climate_key = "('my project', 'climate change')"
resource_key = "('my project', 'resources')"

if isinstance(mc_results, dict) and len(mc_results) > 0:
    print(f"Extracting both climate change and resource impacts from {len(mc_results)} Monte Carlo iterations...")
    
    climate_impacts = []
    resource_impacts = []
    successful_iterations = 0
    
    for i, iteration_result in mc_results.items():
        try:
            if "error" in iteration_result:
                continue
            # Access impact values from the Impacts DataFrame
            impacts_df = iteration_result["Impacts"]
            climate_value = impacts_df.loc[climate_key, "Value"]
            resource_value = impacts_df.loc[resource_key, "Value"]
            
            if climate_value is not None and resource_value is not None:
                climate_impacts.append(climate_value)
                resource_impacts.append(resource_value)
                successful_iterations += 1
        except (KeyError, AttributeError) as e:
            print(f"Warning: Could not extract impacts from iteration {i}: {e}")
            continue
    
    print(f"Successfully extracted impacts from {successful_iterations}/{len(mc_results)} iterations")
    
    if len(climate_impacts) > 0 and len(resource_impacts) > 0:
        # Convert to numpy arrays for analysis
        climate_impacts = np.array(climate_impacts)
        resource_impacts = np.array(resource_impacts)
        
        # Calculate statistical measures for climate change
        climate_mean = np.mean(climate_impacts)
        
        # Calculate statistical measures for resource use
        resource_mean = np.mean(resource_impacts)

        climate_std, resource_std = np.std(climate_impacts, ddof=1), np.std(resource_impacts, ddof=1)
        climate_90 = np.percentile(climate_impacts, [5, 95])
        resource_90 = np.percentile(resource_impacts, [5, 95])
        correlation = np.corrcoef(climate_impacts, resource_impacts)[0, 1]
        print(f"Climate change: mean {climate_mean:.4f}, std {climate_std:.4f}, "
              f"90% interval [{climate_90[0]:.4f}, {climate_90[1]:.4f}]")
        print(f"Resource use:   mean {resource_mean:.4f}, std {resource_std:.4f}, "
              f"90% interval [{resource_90[0]:.4f}, {resource_90[1]:.4f}]")
        print(f"Correlation between the two impacts: {correlation:.2f}")
Extracting both climate change and resource impacts from 500 Monte Carlo iterations...
Successfully extracted impacts from 500/500 iterations
Climate change: mean 1.7676, std 0.2798, 90% interval [1.3284, 2.2473]
Resource use:   mean 0.0017, std 0.0003, 90% interval [0.0012, 0.0023]
Correlation between the two impacts: -0.04

7.5 Uncertainty Visualization#

The joint density of the two optimal impacts over all samples, with their means.

if climate_impacts is not None and resource_impacts is not None and len(climate_impacts) > 0:
    import seaborn as sns
    
    # Create elegant dual-objective uncertainty visualization
    fig, ax1 = plt.subplots(1, 1, figsize=(6, 4))
    
    # 1. Seaborn KDE plot for density distribution
    ax1.set_title('Climate vs. Resource Density Distribution', fontsize=14, fontweight='bold')
    sns.kdeplot(x=climate_impacts, y=resource_impacts, cmap='viridis', fill=True, 
                thresh=0.02, ax=ax1, alpha=0.7)
    
    # Add mean lines
    ax1.axvline(climate_mean, color='red', linestyle='--', linewidth=2, alpha=0.8, 
                label=f'Climate Mean: {climate_mean:.4f}')
    ax1.axhline(resource_mean, color='red', linestyle='--', linewidth=2, alpha=0.8, 
                label=f'Resource Mean: {resource_mean:.4f}')
    
    ax1.set_xlabel('Climate Change Impact')
    ax1.set_ylabel('Resource Use Impact')
    ax1.legend()
    ax1.grid(True, alpha=0.3)

    plt.tight_layout()
    plt.show()

else:
    print("❌ Cannot create visualization - missing valid impact data for both objectives")
../../_images/aa14f476aa8a53173cedfab4dbf67675eec4c65c1d94a377b4ce900f02241494.png

See also: pulpo.utils.uncertainty for the uncertainty of an impact from its declared distributions: closed-form moments and chance-constrained optimization in reduced space (uncertainty toy notebook).

8. Infeasible Problems#

Limits that cannot all be met make the optimization problem infeasible. PULPO then raises a SolveError (from pulpo.utils.optimizer) and loads nothing: the worker keeps its previous solution. The same happens when a solve stops early, for example at a time limit.

Limits belong on the processes, flows and impacts that really have one (lower_limit, upper_limit, upper_elem_limit, upper_imp_limit, …), as in section 4. Bounds on every activity through default_limits have no physical meaning; they raise a FutureWarning, and default_limits may be deprecated in a near-future release.

8.1 The Lowest Achievable Impact#

Solve the unconstrained problem to find the lowest climate change impact the system can reach.

from pulpo.utils.optimizer import SolveError

climate = "('my project', 'climate change')"
methods_default = {climate: 1, "('my project', 'resources')": 0}
pulpo_worker_default = pulpo.PulpoOptimizer(project, database, methods_default)
pulpo_worker_default.get_lci_data()

pulpo_worker_default.instantiate(choices=choices_unconstrained, demand=demand, lower_limit=lower_bound)
pulpo_worker_default.solve()
lowest_climate = pulpo_worker_default.instance.impacts_calculated[climate].value
print(f"Lowest achievable climate change impact: {lowest_climate:.6f}")
Creating Instance
Instance created
optimal solution found:  1.760427258456396
Optimization problem solved using Highspy
Lowest achievable climate change impact: 1.760427

8.2 A Limit Below It#

Ask for half of that impact. No combination of technologies reaches it, so the solve raises a SolveError:

pulpo_worker_default.instantiate(choices=choices_unconstrained, demand=demand, lower_limit=lower_bound,
                                 upper_imp_limit={climate: lowest_climate / 2})
try:
    pulpo_worker_default.solve()
except SolveError as error:
    print(f"No solution: {error}")
Creating Instance
Instance created
No solution: HiGHS did not solve the model to optimality (infeasible); the instance keeps its previous values. The solver's results are on the exception's .results.

9. Dependent Constraints: Technology Coupling#

Real-world systems often require interdependent technology operation rules. PULPO’s dependent constraints enable specification of relationships between different process scaling levels, allowing for operational constraints, market share limitations, or technology coupling requirements.

Use Cases:

  • Maximum renewable energy penetration limits

  • Technology deployment ratios

  • Supply chain coordination constraints

  • Market share limitations

  • Infrastructure coupling requirements

9.1 Market Penetration Constraint#

Implement a realistic constraint where wind electricity can constitute a maximum of 80% of total electricity production. This represents grid stability or market penetration limitations for renewable energy sources.

Constraint Formulation:

  • Wind electricity ≤ 0.8 × (Wind electricity + Gas electricity)

  • Mathematical rearrangement: Wind ≤ 4 × Gas electricity

  • Ensures renewable energy integration within practical bounds

# Configure PULPO worker for dependent constraints demonstration
methods_dependent = {"('my project', 'climate change')": 1, "('my project', 'resources')": 0}
pulpo_worker_dep = pulpo.PulpoOptimizer(project, database, methods_dependent)
pulpo_worker_dep.get_lci_data()

# Retrieve electricity technologies specifically by name for clarity
wind_electricity = pulpo_worker_dep.retrieve_processes(processes=['wind electricity'])
gas_electricity = pulpo_worker_dep.retrieve_processes(processes=['natural gas electricity'])

print("Retrieved electricity technologies:")
print(f"  Wind electricity: {wind_electricity[0]['name']}")
print(f"  Gas electricity: {gas_electricity[0]['name']}")

# Multi-product demand (1 kg methanol + 2 kg ozone)
demand_dep = {methanol_process[0]: 1, ozone_process[0]: 2}
print(f"\nDemand specification: {len(demand_dep)} products")
Retrieved electricity technologies:
  Wind electricity: wind electricity
  Gas electricity: natural gas electricity

Demand specification: 2 products
# Define dependent constraint: Wind ≤ 80% of total electricity
# Mathematical form: wind ≤ 4 × gas_electricity (derived from wind ≤ 0.8 × (wind + gas))
dependent_constraints = {
    'wind_penetration_limit': {
        'left': {wind_electricity[0]: 1},        # Wind electricity scaling
        'right': {gas_electricity[0]: 4}         # 4 × Gas electricity scaling
    }
}

# Technology choices (allow both options)
choices_dep = {
    "Electricity": {wind_electricity[0]: float('inf'), gas_electricity[0]: float('inf')},
    "Hydrogen": {hydrogen_processes[0]: float('inf'), hydrogen_processes[1]: float('inf')},
    "Oxygen": {oxygen_processes[0]: float('inf'), oxygen_processes[1]: float('inf')},
}

# Lower bounds for physical feasibility
lower_bound_dep = {oxygen_byproduct[0]: 0}

print("Dependent constraint configuration:")
print("  Constraint: Wind electricity ≤ 4 × Gas electricity")
print("  Interpretation: Wind ≤ 80% of total electricity production")
print(f"  Technology choices: {len(choices_dep)} categories available")
Dependent constraint configuration:
  Constraint: Wind electricity ≤ 4 × Gas electricity
  Interpretation: Wind ≤ 80% of total electricity production
  Technology choices: 3 categories available
# Solve optimization with dependent constraints
print("🔄 Solving optimization with wind penetration constraint...")

pulpo_worker_dep.instantiate(
    choices=choices_dep, 
    demand=demand_dep, 
    lower_limit=lower_bound_dep,
    dependent_constraints=dependent_constraints
)

results_dep = pulpo_worker_dep.solve()

# Extract scaling values for analysis
scaling_results = {}
for activity_name, process_id in pulpo_worker_dep.lci_data['process_map'].items():
    if pulpo_worker_dep.instance.scaling_vector[process_id].value > 1e-6:  # Only show active processes
        scaling_results[activity_name[1]] = pulpo_worker_dep.instance.scaling_vector[process_id].value

print("✅ Optimization completed successfully!")
print("\nActive process scaling values:")
for process_name, scaling_value in scaling_results.items():
    print(f"  {process_name}: {scaling_value:.4f}")
🔄 Solving optimization with wind penetration constraint...
Creating Instance
Instance created
optimal solution found:  3.0636673756672597
Optimization problem solved using Highspy
✅ Optimization completed successfully!

Active process scaling values:
  natural gas extraction: 0.5270
  natural gas electricity: 2.7738
  wind electricity: 11.0952
  water supply: 0.0017
  O2 ASU: 0.4800
  hydrogen electrolysis: 0.1900
  O2-market: 1.5200
  direct air capture: 1.3750
  methanol synthesis: 1.0000
  ozone production: 2.0000
# Verify constraint satisfaction and analyze electricity mix
wind_scaling = scaling_results.get('wind electricity', 0)
gas_scaling = scaling_results.get('natural gas electricity', 0)
total_electricity = wind_scaling + gas_scaling

if total_electricity > 0:
    wind_percentage = (wind_scaling / total_electricity) * 100
    constraint_satisfied = wind_scaling <= (4 * gas_scaling + 1e-6)
    climate_impact = pulpo_worker_dep.instance.impacts_calculated["('my project', 'climate change')"].value
    
    print(f"\n📊 CONSTRAINT VERIFICATION")
    print(f"Wind: {wind_scaling:.4f} | Gas: {gas_scaling:.4f} | Wind %: {wind_percentage:.1f}%")
    print(f"Constraint satisfied: {constraint_satisfied} (wind ≤ 4×gas: {wind_scaling:.4f} ≤ {4*gas_scaling:.4f})")
    print(f"Climate impact: {climate_impact:.4f}")
else:
    print("⚠️ No electricity production detected")
📊 CONSTRAINT VERIFICATION
Wind: 11.0952 | Gas: 2.7738 | Wind %: 80.0%
Constraint satisfied: True (wind ≤ 4×gas: 11.0952 ≤ 11.0952)
Climate impact: 3.0637

10. Supply Specification: Fixing a Process Output#

All previous examples are demand-driven: a functional unit is specified and the system must deliver it. PULPO also supports supply-driven analysis. When a process is given identical lower and upper limits, its scaling is fixed to exactly that value, and PULPO automatically activates a slack variable for the corresponding reference product. The slack relaxes the product balance, so the fixed output does not need to be absorbed by a final demand — the question becomes:

“What is the optimal system configuration (and what are its impacts) when this process runs at exactly this level?”

Use Cases:

  • Assessing an existing plant operating at a known production level

  • Impact accounting for a fixed production quota rather than a functional unit

  • Analyzing systems where the output is determined by upstream availability (e.g., waste treatment capacity) instead of downstream demand

10.1 Supply-Driven Optimization: Fixed Methanol Output#

Fix the methanol synthesis process to a scaling of exactly 1 by passing the same value as lower and upper limit. No demand is specified — the fixed supply alone drives the system. The technology choices remain open, so the optimizer still selects the impact-minimal electricity, hydrogen, and oxygen sources needed to run the fixed process.

# Configure a fresh worker for the supply-driven demonstration
methods_supply = {"('my project', 'climate change')": 1, "('my project', 'resources')": 0}
pulpo_worker_supply = pulpo.PulpoOptimizer(project, database, methods_supply)
pulpo_worker_supply.get_lci_data()

# Identical lower and upper limit -> scaling is fixed and the supply slack is activated
methanol_supply = {methanol_process[0]: 1}

# Technology choices remain open for the optimizer
choices_supply = {
    "Electricity": {electricity_processes[0]: float('inf'), electricity_processes[1]: float('inf')},
    "Hydrogen": {hydrogen_processes[0]: float('inf'), hydrogen_processes[1]: float('inf')},
    "Oxygen": {oxygen_processes[0]: float('inf'), oxygen_processes[1]: float('inf')},
}

# Note: no demand is passed!
pulpo_worker_supply.instantiate(
    choices=choices_supply,
    lower_limit={**methanol_supply, oxygen_byproduct[0]: 0},
    upper_limit=methanol_supply,
)
results_supply = pulpo_worker_supply.solve()

supply_climate = pulpo_worker_supply.instance.impacts_calculated["('my project', 'climate change')"].value
print(f"Supply-driven climate impact: {supply_climate:.6f}")
Creating Instance
Instance created
optimal solution found:  1.7471600582975855
Optimization problem solved using Highspy
Supply-driven climate impact: 1.747160

Results Analysis: Slack Activation#

The slack variable exists only for the methanol product (the model carries no slack for any other product). Its value equals the fixed methanol output that is not consumed by any demand — the balance production = demand + slack closes with demand = 0, so the slack absorbs the full supply of 1.

# The slack absorbs the fixed methanol output that no demand consumes
slack = pulpo_worker_supply.extract_results()['Slack']
print("Slack (exists only for the supply product):")
print(slack)

methanol_id = pulpo_worker_supply.lci_data['process_map'][methanol_process[0].key]
print(f"Methanol synthesis scaling (fixed): {pulpo_worker_supply.instance.scaling_vector[methanol_id].value:.4f}")
Slack (exists only for the supply product):
    Value
10    1.0
Methanol synthesis scaling (fixed): 1.0000

10.2 Comparison with the Demand-Driven Formulation#

For this system the supply-driven and demand-driven formulations coincide: methanol is not consumed by any other process, so fixing the gross output of 1 is equivalent to demanding a net output of 1. In systems where the fixed product is also consumed internally, the two formulations differ — the supply specification fixes how much the process runs, while a demand fixes how much must be left over for the functional unit.

# Same 1 kg methanol as a conventional demand for comparison
pulpo_worker_supply.instantiate(
    choices=choices_supply,
    demand={methanol_process[0]: 1},
    lower_limit={oxygen_byproduct[0]: 0},
)
results_demand = pulpo_worker_supply.solve()
demand_climate = pulpo_worker_supply.instance.impacts_calculated["('my project', 'climate change')"].value

print(f"Supply-driven (fixed output):    {supply_climate:.6f}")
print(f"Demand-driven (functional unit): {demand_climate:.6f}")
print(f"Difference: {abs(supply_climate - demand_climate):.2e}")
Creating Instance
Instance created
optimal solution found:  1.7471600582975855
Optimization problem solved using Highspy
Supply-driven (fixed output):    1.747160
Demand-driven (functional unit): 1.747160
Difference: 0.00e+00

11. Goal Programming: Average Transgression Level#

All previous examples minimize a (weighted) sum of impacts. In many policy-oriented settings — e.g. Planetary-Boundary assessments expressed through the EF 3.1 categories — the question is different: each impact category has a budget that should not be transgressed, but some transgression may be unavoidable. PULPO supports this as a goal-programming objective: for each category \(h\) with a soft limit (goal) \(L_h\), the transgression level is \(TL_h = \text{Impact}_h / L_h\), and the optimizer minimizes the average transgression

\[\min \; \frac{1}{K} \sum_{h=1}^{K} \max\left(0, \; TL_h - 1\right)\]

Categories within their budget contribute nothing; transgressed categories contribute proportionally to how far they exceed it. Unlike the hard upper_imp_limit constraints from Section 4, goals can be exceeded — the problem stays feasible and the solution reports how much each budget is blown.

11.1 An (Invented) Regulatory Budget Scenario#

Suppose a regulator allocates the methanol producer a per-kg budget in all three impact categories (the numbers are made up for demonstration):

Category

Budget \(L_h\)

Attainable minimum

climate change

2.0

1.747 (climate-optimal solution)

air quality

0.05

0.0

resources

0.001

0.00057 (resource-optimal solution)

The budgets are deliberately conflicting: the climate-optimal technology mix (wind + electrolysis) overshoots the resource budget by 71%, while the resource-optimal mix (gas + SMR) overshoots the climate budget by 73% and the air-quality budget as well. No corner solution satisfies everything — the goal objective must find the best compromise.

# Fresh worker with ALL three impact categories (weights are ignored in goal mode)
methods_goal = {
    "('my project', 'climate change')": 1,
    "('my project', 'air quality')": 1,
    "('my project', 'resources')": 1,
}
pulpo_worker_goal = pulpo.PulpoOptimizer(project, database, methods_goal)
pulpo_worker_goal.get_lci_data()

# The invented per-kg-methanol budgets (soft limits)
imp_goals = {
    "('my project', 'climate change')": 2.0,
    "('my project', 'air quality')": 0.05,
    "('my project', 'resources')": 0.001,
}

# Same system as Section 3: 1 kg methanol, free technology choices
choices_goal = {
    "Electricity": {electricity_processes[0]: float('inf'), electricity_processes[1]: float('inf')},
    "Hydrogen": {hydrogen_processes[0]: float('inf'), hydrogen_processes[1]: float('inf')},
    "Oxygen": {oxygen_processes[0]: float('inf'), oxygen_processes[1]: float('inf')},
}

pulpo_worker_goal.instantiate(
    choices=choices_goal,
    demand={methanol_process[0]: 1},
    lower_limit={oxygen_byproduct[0]: 0},
    imp_goals=imp_goals,
    objective='goal',
)
results_goal = pulpo_worker_goal.solve()
print(f"Average transgression level: {pulpo_worker_goal.instance.OBJ():.4f}")
Creating Instance
Instance created
optimal solution found:  0.15903260059682106
Optimization problem solved using Highspy
Average transgression level: 0.1590

11.2 Results Analysis: Budget Status and Compromise Mix#

extract_results() gains a Transgressions table (also shown by summarize_results() and exported by save_results()): per goal category the impact, the budget, the transgression level \(TL\), and the transgression \(\max(0, TL-1)\).

transgressions = pulpo_worker_goal.extract_results()["Transgressions"]
transgressions
Impact Goal TL Transgression
Method
('my project', 'resources') 0.001283 0.001 1.283437 0.283437
('my project', 'climate change') 2.387322 2.000 1.193661 0.193661
('my project', 'air quality') 0.050000 0.050 1.000000 0.000000
pulpo_worker_goal.summarize_results(zeroes=True)

Total Impact(s)

Weight Value
Method
('my project', 'air quality') 1 0.050000
('my project', 'climate change') 1 2.387322
('my project', 'resources') 1 0.001283

Goal Transgressions

Impact Goal TL Transgression
Method
('my project', 'resources') 0.001283 0.001 1.283437 0.283437
('my project', 'climate change') 2.387322 2.000 1.193661 0.193661
('my project', 'air quality') 0.050000 0.050 1.000000 0.000000

Choices Made

Electricity

Value Capacity
Metadata
wind electricity | electricity | GLO 9.215547 inf

Hydrogen

Value Capacity
Metadata
hydrogen electrolysis | hydrogen | GLO 0.118906 inf
hydrogen SMR | hydrogen | GLO 0.071094 inf

Oxygen

Value Capacity
Metadata

Constraints

Constraints Lower

Key Metadata Value Limit
ID
7 (foreground_db, O2-byproduct) O2-byproduct | oxygen | GLO 0.951249 0

The compromise is qualitatively different from every weighted-sum solution seen so far:

  • Hydrogen is blended: the optimizer runs electrolysis and SMR side by side — a split that no single-objective run produces — because fully committing to either technology would blow one of the budgets much further.

  • Air quality sits exactly at its budget (\(TL = 1.0\), transgression 0): the SMR share is expanded precisely up to the point where the air-quality budget becomes binding at the kink of \(\max(0, TL-1)\).

  • Climate change and resources are both mildly transgressed (≈ 19% and ≈ 28%) instead of one of them being blown by > 70%: the average transgression drops to ≈ 0.159, compared to ≈ 0.237 if the climate-optimal mix from Section 3 were assessed against the same budgets.

Tightening or relaxing individual budgets shifts this balance — e.g. lowering the resource budget pushes the mix back towards SMR until the climate budget takes over as the binding concern.

11.3 Yearly Budgets in the Time-Indexed Model#

The same objective is available in the time-indexed formulation (pulpo_time.PulpoOptimizerTime, see the elec_time_toy notebook). There the goals apply to the impacts aggregated over all timesteps — the natural fit for e.g. an hourly-resolved dispatch problem (8760 steps) with a total yearly budget per category:

worker = pulpo_time.PulpoOptimizerTime(project, database, methods)
worker.get_lci_data()
worker.instantiate(
    choices=choices, demand=demand,
    time_steps=list(range(8760)),                    # hourly resolution
    imp_goals={gwp_method: yearly_co2_budget},        # budget on sum over all hours
    objective='goal',
)
worker.solve()

The transgression of each category is then evaluated on \(\sum_t \text{Impact}_{t,h}\), i.e. the optimizer is free to exceed the hourly-average budget in some hours as long as the yearly total works out.