<jupyter_start><jupyter_text>Simple Model of the LNG MarketThe object of this model is to minimize total LNG costs while satisfying constraints on the supplies available at each of the exporting nodes, and satisfying demand requirements at each of the importing nodes.<jupyter_code>import shutil
import sys
import os.path
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import chart_studio
import chart_studio.plotly as py
from tabulate import tabulate
plt.close('all')

from IPython.display import display
from geopy.geocoders import Nominatim
from geopy import distance
from pyomo.environ import *

geolocator = Nominatim(user_agent="Your_Name")
pd.set_option('precision', 3)
MMBtu_LNG = 21.04 # 1m³ of LNG = 21.04MMBtu
MMBtu_Gas = 0.036 # 1m³ of Gas = 0.036MMBtu
LNG_to_Gas = 584 # 1m³ of LNG = 584m³ of gas (http://www.lngplants.com/conversiontables.html)<jupyter_output><empty_output><jupyter_text>1 Data   1.1 Importing NodesTen largest importing countries according to "BP Stats Review (2020)" are taken into consideration. Each importing node has volume of imported LNG as an attribute.<jupyter_code># Data for importing and exporting countries from "BP Stats Review (2020)"

# "Import (MMBtu)" column consists of imports in m³ multiplied with conversion factor

Imp_Countries = pd.DataFrame({
    
    "Country" : ["Japan", "China", "South Korea", "India", "France", "Spain",
                 "UK", "Italy", "Turkey", "Belgium"],
    
     "Import (MMBtu)" : [105.5*10**9 * MMBtu_Gas, 84.8*10**9 * MMBtu_Gas, 55.6*10**9 * MMBtu_Gas,
                         32.9*10**9 * MMBtu_Gas, 22.9*10**9 * MMBtu_Gas, 21.9*10**9 * MMBtu_Gas,
                         18*10**9 * MMBtu_Gas, 13.5*10**9 * MMBtu_Gas, 12.9*10**9 * MMBtu_Gas,
                         7.2*10**9 * MMBtu_Gas],
    
    "Regasification Terminal" : ["Sodegaura", "Yancheng", "Incheon", "Dahej", "Dunkirk",
                     "Barcelona", "Isle of Grain", "Rovigo", "Marmara", "Zeebrugge"]
})

#Imp_Countries["Import (MMBtu)"] = (Imp_Countries["Import (MMBtu)"].astype(float)/1000000000).astype(str)


# getting the coordinates of largest regas terminals, not neccesarry anymore, since there is no suitable package
# for calculating sea distances and the API's are not free of charge:(

Imp_Countries["Location"] = Imp_Countries["Regasification Terminal"].apply(geolocator.geocode)
Imp_Countries["Point"]= Imp_Countries["Location"].apply(lambda loc: tuple(loc.point) if loc else None)
Imp_Countries[["Latitude", "Longitude", "Altitude"]] = pd.DataFrame(Imp_Countries["Point"].to_list(), index=Imp_Countries.index)

del Imp_Countries["Altitude"]
del Imp_Countries["Location"]
del Imp_Countries["Point"]


#center the columns, display() method doesn't work anymore, so not necessary for now
def pd_centered(df):
    return df.style.set_table_styles([
        {"selector": "th", "props": [("text-align", "center")]},
        {"selector": "td", "props": [("text-align", "center")]}])

Imp_Countries.style<jupyter_output><empty_output><jupyter_text>1.2 Exporting NodesTen largest exporting countries according to "BP Stats Review (2020)" are taken into consideration. Each exporting node has volume of exported LNG and break even price (feedstock price + liquefaction price) as attributes.<jupyter_code># define DF with data for exporting countries (BP 2020)

Exp_Countries = pd.DataFrame({
    
    "Country" : ["Qatar", "Australia", "USA", "Russia", "Malaysia", "Nigeria",
                 "Trinidad & Tobago", "Algeria", "Indonesia","Oman"],
    
    "Export (MMBtu)" : [107.1*10**9 * MMBtu_Gas, 104.7*10**9 * MMBtu_Gas, 47.5*10**9 * MMBtu_Gas,
                        39.4*10**9 * MMBtu_Gas, 35.1*10**9 * MMBtu_Gas, 28.8*10**9 * MMBtu_Gas,
                        17.0*10**9 * MMBtu_Gas, 16.9*10**9 * MMBtu_Gas, 16.5*10**9 * MMBtu_Gas,
                        14.1*10**9 * MMBtu_Gas],
    
    "Break Eeven Price ($/MMBtu)" : [2.5, 6.5, 6.0, 4.5, 6.0, 4.08, 5.1, 2.5, 6.0, 3.5],
    
    #"Feedstock ($/MMBtu)" : [0.5, 3.0, 3.0, 0.5, 3.0, 2.2, 3.0, 0.5, 2.0, 1.5],

    #"Liquefaction ($/MMBtu)" : [2.0, 3.5, 3.0, 4.0, 3.0, 1.88, 2.1, 2.0, 4.0, 2.0],
    
    "Liquefaction Terminal" : ["Ras Laffan", "Gladstone", "Sabine Pass", "Sabetta",
                               "Bintulu", "Bonny Island", "Point Fortin",
                               " Bethioua", "Bontang", "Qalhat"],
})

#Exp_Countries["Export (Bcm)"] = (Exp_Countries["Export (Bcm)"].astype(float)/1000000000).astype(str)


# Getting the coordinates of liquefaction terminals, not neccesarry anymore, since there is no suitable package
# for calculating sea distances and the API's are not free of charge:(

Exp_Countries["Location"] = Exp_Countries["Liquefaction Terminal"].apply(geolocator.geocode)
Exp_Countries["Point"] = Exp_Countries["Location"].apply(lambda loc: tuple(loc.point) if loc else None)
Exp_Countries[["Latitude", "Longitude", "Altitude"]] = pd.DataFrame(Exp_Countries["Point"].to_list(), index=Exp_Countries.index)

# remove unnecessary columns
del Exp_Countries["Altitude"]
del Exp_Countries["Location"]
del Exp_Countries["Point"]

#display(pd_centered(Exp_Countries))
Exp_Countries.style<jupyter_output><empty_output><jupyter_text>1.3 Sea distances among the importing and exporting nodes<jupyter_code># Define DF with distances among largest liquefaction and regas terminals in nautical miles
# In case of lack of route to the specific terminal, next largest port was taken
# sea-distances.org, searoutes.com

# if route over canal: distance rounded on 0 for suez and on 5 for panama, needed for later calculations

Distances = pd.DataFrame([(6512, 5846, 6161, 1301, 6240, 4650, 6260, 4380, 3460, 6270),
                          
                          (3861, 4134, 4313, 5918, 11650, 10060, 11670, 9790, 8870, 11680),
                          
                          (9205, 10075, 10005, 9645, 4881, 5206, 4897, 6354, 6341, 4908),
                          
                          (4929, 5596, 5671, 10160, 2494, 4252, 2511, 5546, 5481, 2477),
                          
                          (2276, 1656, 1998, 3231, 8980, 7390, 9000, 7120, 6200, 9010),
                          
                          (10752, 10088, 10406, 6988, 4287, 3824, 4309, 4973, 4961, 4321),
                          
                          (8885, 9755, 9665, 8370, 3952, 3926, 3974, 5074, 5062, 3984),
                          
                          (9590, 8930, 9240, 4720, 1519, 343, 1541, 1432, 1417, 1552),
                          
                          (2651, 2181, 2487, 3517, 9260, 7670, 9280, 7390, 6470, 9290),
                          
                          (6046, 5379, 5694, 853, 5760, 4180, 5790, 3900, 2980, 5800)

                         ],

                  index=["Qatar", "Australia", "USA", "Russia", "Malaysia", 
                         "Nigeria", "Trinidad & Tobago", "Algeria", "Indonesia"
                         ,"Oman"],

                  columns = ("Japan", "China", "South Korea", "India", "France",
                             "Spain", "UK", "Italy", "Turkey", "Belgium"))

#display(pd_centered(Distances))
Distances<jupyter_output><empty_output><jupyter_text>1.4 LNG CarrierTo keep the model simple it is assumed that there is only one type of the LNG carrier with the following characteristics: - Capacity: 160000 m³ - Spot Charter Rate: 70000 \\$/day - Average speed: 17 knots/h - Boil off: 0.1%  - Heel: 4%<jupyter_code># define DF with lng carrier properties

LNG_Carrier = pd.DataFrame({"Capacity (m³)" : [160000], # average capacity in 2019 160000 (!!! https://www.rivieramm.com/opinion/opinion/lng-shipping-by-numbers-36027)
                            "Spot Charter Rate ($/day)" : [70000], # average spot charter rate in 2019
                            "Speed (knots/h)" : [17],
                            "Bunker (mt/d)" : [110], # !!! Does boil off reduce it and what's the ratio then??? 
                            "Bunkers Price ($/mt)" : [670], 
                            "Boil off" : [1/1000], # the daily guaranteed maximum boil-off or DGMB (0,1%)
                            "Heel" : [4/100]}) # 4% of the cargo is retained as the heel (LNG storage tanks need to be kept cool)

#print(LNG_Carrier_DF.loc[0, "Capacity (m³)"])
#display(pd_centered(LNG_Carrier))
LNG_Carrier.style<jupyter_output><empty_output><jupyter_text>1.5 Additional Costs- Port costs for 3 days (outbound port, destination port and one day in case it has to wait for loading/unloading): 300000 \\$- Canal Fees (Suez or Panama Canal): 1000000 \\$ and 900000 \\$ respectively- Agents and broker fees: 2% of charter cost plus- Insurance: 2,600 $/day<jupyter_code># define DF containing additional costs

Additional_Costs = pd.DataFrame({
    "Port Costs ($ for 3 days)" : [300000], 
    "Suez Canal Costs ($/Cargo)" : [1000000],  # 0.24 $/mmBtU
    "Panama Canal Costs ($/Cargo)" : [900000], # 0.21 $/mmBtu,                     
    "Insurance ($/day)" : [2600],
    "Fees (Percentage of Charter Costs)" : [2/100]})

#display(pd_centered(Additional_Costs))
Additional_Costs.style<jupyter_output><empty_output><jupyter_text>2 Cost and voyage times calculations  2.1 (Return) Voyage Times<jupyter_code>Time = Distances / (LNG_Carrier.loc[0, "Speed (knots/h)"] * 12) # in days

# note to myself: it would be better to define a function which would compute the time!!!
#display(pd_centered(Time))
Time<jupyter_output><empty_output><jupyter_text>2.2 Charter Costs<jupyter_code>Charter_Costs = (Time + 2) * LNG_Carrier.loc[0, "Spot Charter Rate ($/day)"] # in $

# display(pd_centered(Charter_Costs))
Charter_Costs.style<jupyter_output><empty_output><jupyter_text>2.3 Fuel Costs<jupyter_code># this is arbitrary cost (it should be taken either price after relequifaction or destination hub price)
lng_cost = 5 # $/MMBtu 

Fuel_Costs = Time * LNG_Carrier.loc[0, "Boil off"] * 0.98 * LNG_Carrier.loc[0, "Capacity (m³)"] * lng_cost * MMBtu_LNG

# this is arbitrary cost (it should be taken either price after relequifaction or destination hub price)
lng_cost = 5 # $/MMBtu 

#display(pd_centered(Fuel_Costs))
Fuel_Costs<jupyter_output><empty_output><jupyter_text>2.4 Agents and Broker Fees + Insurance<jupyter_code>Fees_and_Insurance = (Time + 2) * Additional_Costs.loc[0, "Insurance ($/day)"] + Charter_Costs * Additional_Costs.loc[0, "Fees (Percentage of Charter Costs)"]

#display(pd_centered(Fees_and_Insurance))
Fees_and_Insurance<jupyter_output><empty_output><jupyter_text>2.5 Suez/Panama Canal Fee<jupyter_code># create empty data frames for canal fees and port costs
Canal_Fees = pd.DataFrame(np.zeros((10,10)),
                   
                    index = ["Qatar", "Australia", "USA", "Russia", "Malaysia", 
                             "Nigeria", "Trinidad & Tobago", "Algeria", "Indonesia",
                             "Oman"],

                    columns = ("Japan", "China", "South Korea", "India", "France",
                               "Spain", "UK", "Italy", "Turkey", "Belgium"))

for index in Distances.index:
    for columns in Distances.columns:
        if (((Distances.loc[index, columns]) % 5 == 0) and ((Distances.loc[index, columns]) %10 != 0)):
            Canal_Fees.loc[index, columns] += Additional_Costs.loc[0, "Panama Canal Costs ($/Cargo)"]
        elif ((Distances.loc[index, columns]) % 10 == 0):
            Canal_Fees.loc[index, columns] += Additional_Costs.loc[0, "Suez Canal Costs ($/Cargo)"]
            
#display(pd_centered(Canal_Fees))
Canal_Fees<jupyter_output><empty_output><jupyter_text>2.6 Port Costs<jupyter_code># empty data frame for port costs
Port_Costs = pd.DataFrame(np.zeros((10,10)),
                   
                    index = ["Qatar", "Australia", "USA", "Russia", "Malaysia", 
                             "Nigeria", "Trinidad & Tobago", "Algeria", "Indonesia",
                             "Oman"],

                    columns = ("Japan", "China", "South Korea", "India", "France",
                               "Spain", "UK", "Italy", "Turkey", "Belgium"))

Port_Costs += Additional_Costs.loc[0, "Port Costs ($ for 3 days)"]

#display(pd_centered(Port_Costs))
Port_Costs<jupyter_output><empty_output><jupyter_text>2.7 Total Costs of Transportation<jupyter_code>Transport_Costs =  Charter_Costs + Fuel_Costs + Fees_and_Insurance + Canal_Fees + Port_Costs
#Transport_Costs

Transport_Costs_MMBtu = Transport_Costs / (0.94 * LNG_Carrier.loc[0, "Capacity (m³)"] * MMBtu_LNG) # 1-heel-boil_off = 0.94

#display(pd_centered(Transport_Costs_MMBtu))
Transport_Costs_MMBtu.style #in $/MMBtu<jupyter_output><empty_output><jupyter_text>2.8 Total CostsTotal cost per unit of transported LNG = Feedstock Price + Liquefaction Cost + Total Cost of Transportation<jupyter_code># taken from EXP_Countries df

Breakeven = pd.DataFrame([(2.5, 2.5, 2.5, 2.5, 2.5, 2.5, 2.5, 2.5, 2.5, 2.5),
                          
                        (6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5, 6.5),
                          
                        (6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0),
                          
                        (4.5, 4.5, 4.5, 4.5, 4.5, 4.5, 4.5, 4.5, 4.5, 4.5),
                          
                        (6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0),
                          
                        (4.08, 4.08, 4.08, 4.08, 4.08, 4.08, 4.08, 4.08, 4.08, 4.08),
                          
                        (5.1, 5.1, 5.1, 5.1, 5.1, 5.1, 5.1, 5.1, 5.1, 5.1),
                          
                        (2.5, 2.5, 2.5, 2.5, 2.5, 2.5, 2.5, 2.5, 2.5, 2.5),
                          
                        (6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0, 6.0),
                          
                        (3.5, 3.5, 3.5, 3.5, 3.5, 3.5, 3.5, 3.5, 3.5, 3.5)

                       ],

                  index = ["Qatar", "Australia", "USA", "Russia", "Malaysia", "Nigeria",
                           "Trinidad & Tobago", "Algeria", "Indonesia","Oman"],

                  columns = ("Japan", "China", "South Korea", "India", "France",
                             "Spain", "UK", "Italy", "Turkey", "Belgium"))

Total_Costs = Breakeven + Transport_Costs_MMBtu

#display(pd_centered(Total_Costs))
Total_Costs.style # in $/MMBtu<jupyter_output><empty_output><jupyter_text>3 Optimization ModelThe aim of the model is to calculate the optimal distribution of a commodity from a set of exporters to a set of importers in such a way that the overall costs are minimal, while satisfying constraints on the supplies available at each of the exporters, and satisfying demand requirements at each of the importers.**After the shippments between each exporting and importing node are calculated the hub price in importing country can be determinated following the merrit order effect: the shippment with highest price per unit of LNG is also the hub price. (Is this correct???)**![drawing-2.png](attachment:drawing-2.png)The situation can be modeled by links connecting a set of nodes representing exporters to a set of nodes representing importers. For each link there is a parameter TC[e,i] denoting the total cost per MMBtu of the transported LNG over the link. A non-negative decision variable x[e,i] represents the volume of shipped LNG over the link.Objective function:$$minimize:     Cost = \sum_{e\in Exporters}\sum_{i\in Importers}TC[e,i]x[e,i]$$Constraints:- All shipments from an exporting node can not exceed the export volume of that exporting node:$$ \sum_{e\in Exporters}x[e,i] \leq Export[e] ,        \forall e \in Exporters$$ - Shipments to each importing node must satisfy their demand.$$ \sum_{i\in Importers}x[e,i] = Import[i]  ,          \forall i \in Importers$$  3.1 Data File<jupyter_code># Extract "Countries", "Exports" and "Imports" from "Exp_Countries" and "Imp_Countries" dfs
# and create dictionaries

Supply = Exp_Countries.set_index("Country")["Export (MMBtu)"].to_dict()
#print(Supply)

Demand = Imp_Countries.set_index("Country")["Import (MMBtu)"].to_dict()
#print(Demand)

# I would really like to get "Cost" dictionary straight from the "Total_Costs" df, since doing it 
# this way is real pain in the ...
Cost = {
    ("Japan", "Qatar"): 3.473,
    ("China", "Qatar"): 3.387,
    ("South Korea", "Qatar"): 3.428,
    ("India", "Qatar"): 2.797,
    ("France", "Qatar"): 3.725,
    ("Spain", "Qatar"): 3.519,
    ("UK", "Qatar"): 3.728,
    ("Italy", "Qatar"): 3.484,
    ("Turkey", "Qatar"): 3.364,
    ("Belgium", "Qatar"): 3.729,
    ("Japan", "Australia"): 7.129,
    ("China", "Australia"): 7.165,
    ("South Korea", "Australia"): 7.188,
    ("India", "Australia"): 7.396,
    ("France", "Australia"): 8.427,
    ("Spain", "Australia"): 8.220,
    ("UK", "Australia"): 8.429,
    ("Italy", "Australia"): 8.185,
    ("Turkey", "Australia"): 8.066,
    ("Belgium", "Australia"): 8.430,
    ("Japan", "USA"): 7.581,
    ("China", "USA"): 7.694,
    ("South Korea", "USA"): 7.685,
    ("India", "USA"): 7.638,
    ("France", "USA"): 6.762,
    ("Spain", "USA"): 6.804,
    ("UK", "USA"): 6.764,
    ("Italy", "USA"): 6.953,
    ("Turkey", "USA"): 6.951,
    ("Belgium", "USA"): 6.765,
    ("Japan", "Russia"): 5.268,
    ("China", "Russia"): 5.354,
    ("South Korea", "Russia"): 5.364,
    ("India", "Russia"): 6.233,
    ("France", "Russia"): 4.952,
    ("Spain", "Russia"): 5.180,
    ("UK", "Russia"): 4.954,
    ("Italy", "Russia"): 5.348,
    ("Turkey", "Russia"): 5.339,
    ("Belgium", "Russia"): 4.950,
    ("Japan", "Malaysia"): 6.424,
    ("China", "Malaysia"): 6.343,
    ("South Korea", "Malaysia"): 6.388,
    ("India", "Malaysia"): 6.548,
    ("France", "Malaysia"): 7.580,
    ("Spain", "Malaysia"): 7.374,
    ("UK", "Malaysia"): 7.583,
    ("Italy", "Malaysia"): 7.339,
    ("Turkey", "Malaysia"): 7.220,
    ("Belgium", "Malaysia"): 7.584,
    ("Japan", "Nigeria"): 5.603,
    ("China", "Nigeria"): 5.517,
    ("South Korea", "Nigeria"): 5.558,
    ("India", "Nigeria"): 5.115,
    ("France", "Nigeria"): 4.765,
    ("Spain", "Nigeria"): 4.705,
    ("UK", "Nigeria"): 4.767,
    ("Italy", "Nigeria"): 4.854,
    ("Turkey", "Nigeria"): 4.852,
    ("Belgium", "Nigeria"): 4.769,
    ("Japan", "Trinidad & Tobago"): 6.639,
    ("China", "Trinidad & Tobago"): 6.752,
    ("South Korea", "Trinidad & Tobago"): 6.740,
    ("India", "Trinidad & Tobago"): 6.601,
    ("France", "Trinidad & Tobago"): 5.741,
    ("Spain", "Trinidad & Tobago"): 5.738,
    ("UK", "Trinidad & Tobago"): 5.744,
    ("Italy", "Trinidad & Tobago"): 5.887,
    ("Turkey", "Trinidad & Tobago"): 5.885,
    ("Belgium", "Trinidad & Tobago"): 5.745,
    ("Japan", "Algeria"): 4.159,
    ("China", "Algeria"): 4.074,
    ("South Korea", "Algeria"): 4.114,
    ("India", "Algeria"): 3.528,
    ("France", "Algeria"): 2.826,
    ("Spain", "Algeria"): 2.673,
    ("UK", "Algeria"): 2.828,
    ("Italy", "Algeria"): 2.814,
    ("Turkey", "Algeria"): 2.812,
    ("Belgium", "Algeria"): 2.830, 
    ("Japan", "Indonesia"): 6.472,
    ("China", "Indonesia"): 6.411,
    ("South Korea", "Indonesia"): 6.451,
    ("India", "Indonesia"): 6.585,
    ("France", "Indonesia"): 7.617,
    ("Spain", "Indonesia"): 7.410,
    ("UK", "Indonesia"): 7.619,
    ("Italy", "Indonesia"): 7.374,
    ("Turkey", "Indonesia"): 7.255,
    ("Belgium","Indonesia"): 7.620,
    ("Japan", "Oman"): 4.413,
    ("China", "Oman"): 4.326,
    ("South Korea", "Oman"): 4.367,
    ("India", "Oman"): 3.739,
    ("France", "Oman"): 4.663,
    ("Spain", "Oman"): 4.458,
    ("UK", "Oman"): 4.667,
    ("Italy", "Oman"): 4.421,
    ("Turkey", "Oman"): 4.302,
    ("Belgium", "Oman"): 4.668    
}

#print(Cost)<jupyter_output><empty_output><jupyter_text>3.1 Model File<jupyter_code># Step 0: Create an instance of the model
model = ConcreteModel()
model.dual = Suffix(direction=Suffix.IMPORT)

# Step 1: Define index sets
CUS = list(Demand.keys()) # Customer
SRC = list(Supply.keys()) # Source

# Step 2: Define the decision 
model.x = Var(CUS, SRC, domain = NonNegativeReals)

# Step 3: Define Objective
model.Cost = Objective(
    expr = sum([Cost[c,s]*model.x[c,s] for c in CUS for s in SRC]),
    sense = minimize)

# Step 4: Constraints
model.src = ConstraintList()
for s in SRC:
    model.src.add(sum([model.x[c,s] for c in CUS]) <= Supply[s])
        
model.dmd = ConstraintList()
for c in CUS:
    model.dmd.add(sum([model.x[c,s] for s in SRC]) == Demand[c])
    
results = SolverFactory('glpk').solve(model)
results.write()<jupyter_output># ==========================================================
# = Solver Results                                         =
# ==========================================================
# ----------------------------------------------------------
#   Problem Information
# ----------------------------------------------------------
Problem: 
- Name: unknown
  Lower bound: 69384718800.0
  Upper bound: 69384718800.0
  Number of objectives: 1
  Number of constraints: 21
  Number of variables: 101
  Number of nonzeros: 201
  Sense: minimize
# ----------------------------------------------------------
#   Solver Information
# ----------------------------------------------------------
Solver: 
- Status: ok
  Termination condition: optimal
  Statistics: 
    Branch and bound: 
      Number of bounded subproblems: 0
      Number of created subproblems: 0
  Error rc: 0
  Time: 0.05397605895996094
# ----------------------------------------------------------
#   Solution Information
# ----------------[...]<jupyter_text>4 Solution<jupyter_code>#for c in CUS:
#    for s in SRC:
#        print(c, s, model.x[c,s]() / MMBtu / 1000000000)
if 'ok' == str(results.Solver.status):
    print("Total Costs = ",model.Cost(), "$")
    print("\nShipping Table:")
    for s in SRC:
        for c in CUS:
            if model.x[c,s]() > 0:
                print("Ship from ", s," to ", c, ":", model.x[c,s]() / MMBtu_Gas / 1000000000, "bcm")
else:
    print("No Valid Solution Found")
#print("\nShipping Table:")
#for c in CUS:
#    for s in SRC:
#        if model.x[c,s]() > 0:
#               print("Ship from ", s," to ", c, ":", model.x[c,s]() / MMBtu / 1000000000, "bcm")
#exporting transports to dataframe, not necessary anymore, but can be useful :)

#df = pd.DataFrame.from_dict(model.x.extract_values(), orient='index', columns=[str(model.x)])
#df<jupyter_output><empty_output><jupyter_text>4.1 LNG Shippments Table and Sankey Diagram<jupyter_code>sol = pd.DataFrame(np.zeros((10,10)),
                   
                    index = ["Qatar", "Australia", "USA", "Russia", "Malaysia", 
                             "Nigeria", "Trinidad & Tobago", "Algeria", "Indonesia",
                             "Oman"],

                    columns = ("Japan", "China", "South Korea", "India", "France",
                               "Spain", "UK", "Italy", "Turkey", "Belgium"))

for c in CUS:
    for s in SRC:
        sol.loc[c,s] = model.x[c,s]() / MMBtu_Gas / 1000000000


sol1 = sol.T # transpose

sol2 = sol1.tail(10) # last 10 rows

transported_LNG = sol2[["Japan", "China", "South Korea", "India", "France", "Spain", "UK", "Italy", "Turkey", "Belgium"]]
#transported_LNG
# imports
import pandas as pd
import numpy as np
#import plotly.graph_objs as go
import chart_studio.plotly as py
from plotly.offline import download_plotlyjs, init_notebook_mode, plot, iplot
init_notebook_mode(connected=True)

# Nodes & links
nodes = [['ID', 'Label', 'Color'],
        [0,'Qatar','limegreen'],
        [1,'Australia','blue'],
        [2,'USA','crimson'],
        [3,'Russia','red'],
        [4,'Malaysia','purple'],
        [5,'Nigeria','black'],
        [6,'Trinidad & Tobago','deeppink'],
        [7,'Algeria','turquoise'],
        [8,'Indonesia','gold'],
        [9,'Oman','deepskyblue'],
        [10,'Japan','chocolate'],
        [11,'China','crimson'],
        [12,'South Korea','olive'],
        [13,'India','saddlebrown'],
        [14,'France','deepskyblue'],
        [15,'Spain','darkorange'],
        [16,'UK','royalblue'],
        [17,'Italy','limegreen'],
        [18,'Turkey','grey'],
        [19,'Belgium','indigo'],
        ]

# links with your data
links = [['Source','Target','Value','Link Color'],

        # Qa
        [0,11,35.6,'lightgreen'],
        [0,12,38.6,'lightgreen'],
        [0,13,32.9,'lightgreen'],

        # Au
        [1,10,66.1,'cornflowerblue'],
        [1,12,0.5,'cornflowerblue'],

        # USA
        [2,14,8.5,'palevioletred'],
        [2,16,18,'palevioletred'],
        [2,19,7.2,'palevioletred'],
         
        #Ru
        [3,10,39.4,'salmon'],
        
        #Ma
        [4,11,35.1,'orchid'],
        
        #Ni
        [5,15,19.3,'gray'],
        [5,17,9.5,'gray'],
         
        #T&T
        [6,14,14.4,'hotpink'],
        [6,15,2.6,'hotpink'], 
        
        #Al
        [7,17,4,'aquamarine'],
        [7,18,12.9,'aquamarine'],
         
        #In
        [8,12,16.5,'khaki'],
         
        #Om
        [9,11,14.1,'skyblue'],]
         

# Retrieve headers and build dataframes
nodes_headers = nodes.pop(0)
links_headers = links.pop(0)
df_nodes = pd.DataFrame(nodes, columns = nodes_headers)
df_links = pd.DataFrame(links, columns = links_headers)

# Sankey plot setup
data_trace = dict(
    type='sankey',
    domain = dict(
      x =  [0,1],
      y =  [0,1]
    ),
    orientation = "h",
    valueformat = ".0f",
    node = dict(
      pad = 10,
     thickness = 20,
      line = dict(
        color = "black",
        width = 0
      ),
      label =  df_nodes['Label'].dropna(axis=0, how='any'),
      color = df_nodes['Color']
    ),
    link = dict(
      source = df_links['Source'].dropna(axis=0, how='any'),
      target = df_links['Target'].dropna(axis=0, how='any'),
      value = df_links['Value'].dropna(axis=0, how='any'),
      color = df_links['Link Color'].dropna(axis=0, how='any'),
  )
)

layout = dict(
        title = "LNG Shippments",
    height = 772,
    font = dict(
      size = 13),)

fig = dict(data=[data_trace], layout=layout)
iplot(fig, validate=False)<jupyter_output><empty_output><jupyter_text>4.2 Determing the Hub Price in Importing Country<jupyter_code>Merit_Order = pd.DataFrame(np.zeros((len(transported_LNG.columns),len(transported_LNG.columns))),
                   
                    index = ["Qatar", "Australia", "USA", "Russia", "Malaysia", 
                             "Nigeria", "Trinidad & Tobago", "Algeria", "Indonesia",
                             "Oman"],

                    columns = ("Japan", "China", "South Korea", "India", "France",
                               "Spain", "UK", "Italy", "Turkey", "Belgium"))

for i in transported_LNG.index:
    for j in transported_LNG.columns:
        if transported_LNG.loc[i,j] > 0:
            Merit_Order.loc[i,j] = Total_Costs.loc[i,j]

# It's not really a merit order, since columns are not sorted out :)
Merit_Order
Hub_Prices = pd.DataFrame(np.zeros((1,len(Merit_Order.columns))),
                   
                   index = ["Hub Price ($/MMBtu)"],
    
                   columns = ["Japan", "China", "South Korea", "India", "France",
                               "Spain", "UK", "Italy", "Turkey", "Belgium"])


for column in Merit_Order.columns:
    Hub_Prices.loc["Hub Price ($/MMBtu)", column] = Merit_Order[column].max() 

Hub_Prices<jupyter_output><empty_output><jupyter_text>4.3 Plotting the Hub Prices<jupyter_code>hp = Hub_Prices.T
hp.sort_values("Hub Price ($/MMBtu)", ascending=False, inplace = True)
ax = hp.plot(kind = 'bar', title = "Hub Prices", figsize = (10, 7.5), legend = True, fontsize = 14, color="cornflowerblue")
ax.set_xlabel("", fontsize=14)
ax.set_ylabel("$/MMBtu", fontsize=14)
plt.show()<jupyter_output><empty_output><jupyter_text>5 Sensitivity analysis  5.1 Analysis by exporting nodes<jupyter_code>if 'ok' == str(results.Solver.status):
    print("Exporter                   Supply           Shipped  Margin")
    print("\n")
    for m in model.src.keys():
        s = SRC[m-1]
        print("{0:17s}{1:15.1f}{2:17.1f}{3:10.4f}".format(s,Supply[s]/MMBtu_Gas/10**9,model.src[m]()/MMBtu_Gas/10**9,model.dual[model.src[m]]))
else:
    print("No Valid Solution Found")<jupyter_output>Exporter                   Supply           Shipped  Margin


Qatar                      107.1            107.1   -3.7600
Australia                  104.7             66.6    0.0000
USA                         47.5             33.7    0.0000
Russia                      39.4             39.4   -1.8610
Malaysia                    35.1             35.1   -0.8040
Nigeria                     28.8             28.8   -2.0540
Trinidad & Tobago           17.0             17.0   -1.0210
Algeria                     16.9             16.9   -4.0940
Indonesia                   16.5             16.5   -0.7370
Oman                        14.1             14.1   -2.8210<jupyter_text>Caution: The following discussion is based on the solutions of this particular (far from reality) modelThe 'marginal' values are telling us how much the total (for all importing nodes) costs will be increased for each MMBtu increase in the available supply from each exporter. The optimization calculation says that approximately only 2/3 of the "available" LNG from Australia and USA should be used for a minimum cost solution, which rules out any further cost reductions by increasing the available supply in these two exporting nodes. In fact, the supply from these two nodes could be decreased without any harm (the marginal value of LNG from these 2 countries is 0).Other exporting nodes are a different matter. First, all real export volumes are being used. Second, from the marginal value it can be seen that the total costs would be reduced by a certain value for each additional MMBtu of supply.**Question:****Since this model doesn't contain all exporting and importing countries and their respective volumes of export and import (in 2019), should the available liquefaction capacity in each exporting node be used instead of the "available" supply???**  5.2 Analysis by importing nodes<jupyter_code>if 'ok' == str(results.Solver.status):    
    print("Exporter              Demand    Shipped  Margin")
    print("\n")
    for n in model.dmd.keys():
        c = CUS[n-1]
        print("{0:12s}{1:15.1f}{2:10.1f}{3:10.4f}".format(c,Demand[c]/MMBtu_Gas/10**9,model.dmd[n]()/MMBtu_Gas/10**9,model.dual[model.dmd[n]]))
else:
    print("No Valid Solution Found")<jupyter_output>Exporter              Demand    Shipped  Margin


Japan                 105.5     105.5    7.1290
China                  84.8      84.8    7.1470
South Korea            55.6      55.6    7.1880
India                  32.9      32.9    6.5570
France                 22.9      22.9    6.7620
Spain                  21.9      21.9    6.7590
UK                     18.0      18.0    6.7640
Italy                  13.5      13.5    6.9080
Turkey                 12.9      12.9    6.9060
Belgium                 7.2       7.2    6.7650
{'fraction_non_alphanumeric': 0.1672165081829355, 'fraction_numerical': 0.06575457662203248, 'mean_word_length': 2.5534995977473853, 'pattern_counts': {'":': 6, '<': 87, '<?xml version=': 0, '>': 89, 'https://': 1, 'lorem ipsum': 0, 'www.': 2, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 6, 'word_list_counts': {'cursed_substrings.json': 13, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '22', 'dir': 'jupyter-structured-clean-dedup', 'id': '113917', 'license': 'MIT', 'n_tokens_mistral': 12105, 'n_tokens_neox': 10561, 'n_words': 2744, 'path': 'Simple_Model_of_LNG_Market_v_0.5.ipynb', 'repo_name': 'ModellingEnergyMarkets/LNG_Market', 'size': '30918'}