LCA results are not matching the manually calculated results

[Uzair1shah]

Hello,

I hope you are doing well. I have been trying to validate my LCA results by making the similar system in excel and cross checking my LCA emission results. However, I am unable to get the similar results. So, I wanted to ask if there is a way to print/get emission values at each step or in each unit operation? Also, when I insert the emissions by electricity and initial feedstock to zero, it makes the emissions to be zero as well.

S101.set_CF(GWP, 0, basis='kg', units='kg*CO2e')
bst.settings.set_electricity_CF(GWP, 0, basis='kWhr', units='kg*CO2e')

It means it is not accounting for the natural gas emissions in rotary dryer. So, I wanted to see if there is a way I can assign GWP to natural gas well like we do in case of electricity and feedstock. The code is as follows for your reference:

import biosteam as bst
import thermosteam as tmo
import graphviz
from biosteam import main_flowsheet as F, units, Stream, settings
import math as math
import numpy as np
from biosteam.units import VibratingScreen, HammerMill, Mixer, Splitter

import functools

species = {}		#create an empty collection to hold your chemical species
  
#Biomass cellulose
species['cellulose'] = bst.Chemical('cellulose', phase='s', search_db=False, Cp=1.364, rho=1540, formula="C6H10O5", Hf=-1019000, MW =162.1424)            			#Cellulose, verified
species['cellA'] = bst.Chemical('cellA', search_db=False, default=True, phase='s', formula ="C6H10O5", Hf=-1019000, MW=162.1424)                                   	#Activated cellulose
  
# Biomass hemicellulose
species['gmsw'] = bst.Chemical('gmsw', search_db=False, default=True, phase='s', formula="C5H8O4",MW = 132.11612, Hf=-759200)                                      	#Verified
species['xyhw'] = bst.Chemical('xyhw', search_db=False, default=True, phase='s', formula="C5H8O4",MW = 132.11612, Hf=-759200)                                      	#Verified
species['hce1'] = bst.Chemical('hce1', search_db=False, phase='s', formula="C5H8O4",MW = 132.11612, Hf=-759200)                                      				#Verified
species['hce2'] = bst.Chemical('hce2', search_db=False, phase='s', formula="C5H8O4",MW = 132.11612, Hf=-759200)                                      				#Verified

# Biomass lignin
species['ligc'] = bst.Chemical('ligc', search_db=False, phase='s', formula="C15H14O4",MW=258.27376, Hf=-759390)                                    					#Verified
species['ligh'] = bst.Chemical('ligh', search_db=False, phase='s', formula="C22H28O9",MW=436.45892, Hf=-1722700)                                   					#Verified
species['ligO'] = bst.Chemical('ligO', search_db=False, phase='s', formula="C20H22O10",MW=422.38868, Hf=-1847500)                                  					#Verified
species['lig'] = bst.Chemical('lig', search_db=False, phase='s', formula="C11H12O4",MW=208.21388, Hf=-729310)                                     					#Secondary lignin intermediate, verified
species['ligcc'] = bst.Chemical('ligcc', search_db=False, phase='s', formula="C15H14O4",MW=258.27376, Hf=-759390)                                   				#Verified
species['ligoh'] = bst.Chemical('ligoh', search_db=False, phase='s', formula="C19H22O8",MW=378.37888, Hf=-1429200)                                  				#Verified

# Biomass extractives
species['tannin'] = bst.Chemical('tannin', search_db=False, phase='s', formula="C15H12O7", MW=304.25608,Hf=-1079700)                                 				#Tannic acid formula, verified
species['itann'] = bst.Chemical('itann', search_db=False, phase='s', formula="C8H4O4",MW=164.11736, Hf=-616980)                                     				#Intermediate tannin, verified
species['triglycerides'] = bst.Chemical('triglycerides', search_db=False, phase='s', formula="C76H52O46", MW=897.4168, Hf=-1546300)                          		#Triglycerides formula, update the triglycerides Hf

# Biomass inorganics
# This could be defined more explictly later
species['ash'] = bst.Chemical('CaO', phase='s', formula="CaO", Hf=5382027.245)                                         												#Ash, verify HF

# Pyrolysis solid products

species['char'] = bst.Chemical('carbon', phase='s', formula="C")        																							#Carbon, verified, properties from refs 1,2
species['hmwl'] = bst.Chemical('hmwl', search_db=False, phase='l', formula="C24H28O4", MW=380.48392, Hf=-958260)                                   					#Verified, High molecular weight lignin

# Fluidized bed media
species['sand'] = bst.Chemical('SiO2', formula="SiO2", Hf=-910860)                                            														#Sand, inert, properties from ref2

#metaplastic species
species['gco2'] = bst.Chemical('gco2',search_db=False, phase='s', formula="CO2", MW=44.01, Hf=-393510)                                            					#Trapped carbon dioxide, properties from ref2
species['gco'] = bst.Chemical('gco',search_db=False, phase='s', formula="CO",MW=28.01, Hf=-110520)                                               					#Trapped carbon monoxide, properties from ref2
species['gmeoh'] = bst.Chemical('gmeoh',search_db=False, phase='s', formula="CH3OH", MW=32.042, Hf=-238640)                                         				#Trapped methanol, verified, properties from ref2
species['gcoh2_l'] =bst.Chemical('gcoh2_l', search_db=False, phase='s', formula="CH2O",MW=30.026, Hf=-115890)                                        				#Trapped formaldehyde, properties from ref2, loose
species['gc2h6'] = bst.Chemical('gc2h6', search_db=False, phase='s', formula="C2H6", MW=30.07, Hf=-84667)                                           				#Trapped ethane, properties from ref2
species['gch4'] = bst.Chemical('gch4', search_db=False, phase='s', formula="CH4",MW=16.041, Hf=-74840)                                             					#Trapped methane, properties from ref2
species['gc2h4'] = bst.Chemical('gc2h4', search_db=False, phase='s', formula="C2H4", MW=28.052, Hf=52283)                                            				#Trapped ethylene, properties from ref2
species['gc6h5oh'] = bst.Chemical('gc6h5oh', search_db=False, phase='s', formula="C6H6O", MW=94.113, Hf=-96400)                                       				#Trapped phenol, properties from ref3
species['gcoh2_s'] = bst.Chemical('gcoh2_s', search_db=False, phase='s', formula="CH2O",MW=30.026, Hf=-115890)	                                      				#Trapped formaldehyde, properties from ref2, stiff
species['gh2'] = bst.Chemical('gh2', search_db=False, phase='s', formula="H2",MW=2.016, Hf=0)		                                              					#Trapped hydrogen, properties from ref2

#non-biomass species
species['argon'] = bst.Chemical('argon', formula="Ar")  
species['nitrogen'] = bst.Chemical('nitrogen', formula="N2")  
species['oxygen'] = bst.Chemical('oxygen', formula="O2")  

#Pyrolysis products
## Gas phase species
species['ethylene'] = bst.Chemical('ethylene', formula="C2H4", search_ID='ethylene', phase='g')                  
species['methane'] = bst.Chemical('methane', formula="CH4", search_ID='methane',  phase='g')                                   										#ref2
species['hydrogen'] = bst.Chemical('hydrogen', formula="H2", search_ID='hydrogen',  phase='g')                                     									#ref2
species['cacrbonmonoxide'] = bst.Chemical('carbon monoxide', formula="CO", search_ID='carbon monoxide', phase='g')                               					#ref2
species['carbondioxide'] = bst.Chemical('carbon dioxide', formula="CO2", search_ID='carbon dioxide', phase='g')                            							#ref2
species['formaldehyde'] = bst.Chemical('formaldehyde', formula="CH2O", search_ID='formaldehyde', phase='g')                            								#ref2

## Liquid phase species
species['acrolein'] = bst.Chemical('acrolein', formula="C3H4O", search_ID='acrolein', phase='l')      										#https://atct.anl.gov/Thermochemical%20Data/version%201.122e/species/?species_number=1209

species['sinapylaldehyde'] = bst.Chemical('sinapyaldehyde', search_db=False, phase='l', formula="C11H12O4", MW=208.21388, Hf=-483800)
species['methyllinoleate'] = bst.Chemical('methyllinoleate', phase='l', search_db=False, formula="C19H34O2", MW=294.4721, Hf=-604880)
species['glyoxal'] = bst.Chemical('glyoxal', search_ID='glyoxal', phase='l')							    														#NIST webbook
species['acetaldehyde'] = bst.Chemical('acetaldehyde', search_ID='acetaldehyde', phase='l')                       													#Acetaldehyde, properties from ref2 
#species['aceticacid'] = bst.Chemical('aceticacid', search_db=False, default=True, formula="CH3COOH", MW=60.05, Hf=-484930, Hvap=0, Psat=0, Tb=0)                   #Acetic acid, properties from ref2 
species['glycolaldehyde'] = bst.Chemical('glycolaldehyde', search_ID='glycolaldehyde', phase='l')                                    								#Source Pyrolysis of the Simplest Carbohydrate, Glycolaldehyde (CHO-CH2OH), and Glyoxal in a Heated Microreactor. Find additional sources
species['ethanol'] = bst.Chemical('ethanol', search_ID='ethanol', phase='l')                             															#Ethanol, properties from ref2 
species['propionaldehyde'] = bst.Chemical('propionaldehyde', search_ID='propionaldehyde', phase='l')                    							#Source: NIST webbook https://webbook.nist.gov/cgi/cbook.cgi?ID=C123386&Mask=2 
species['hydroxypropanal'] = bst.Chemical('hydroxypropanal', search_ID='hydroxypropanal',phase='l')  												#https://www.chemeo.com/cid/20-360-9/2%2C2-Dimethyl-3-hydroxypropionaldehyde
species['furfural'] = bst.Chemical('furfural', search_ID='furfural', phase='l')
species['xylosan'] = bst.Chemical('xylosan', phase='l', search_db=False, formula="C5H8O4", MW=132.115, Hf=-647000)                             				#University of Pittsburg Pitt Quantum Repository, DOI 10.17614/Q4SQ8QT7Q
species['levoglucosan'] = bst.Chemical('levoglucosan', search_ID='levoglucosan', phase='l')                      											#Gas phase Hf value, ref4.
species['phenol'] = bst.Chemical('phenol', search_ID='phenol', phase='l')                                 													#phenol, properties from ref5
species['hydroxymethylfurfural'] = bst.Chemical('hydroxymethylfurfural', search_db=False, formula="C6H6O3", MW=126.1100, Hf=-437400, phase='l')           	#Hydromethyl furfural DOI:10.1021/ie901012g
species['anisole'] = bst.Chemical('anisole', search_ID='anisole', phase='l')                                  												#DOI: 10.1021/jp507267f
species['pcoumarylalcohol'] = bst.Chemical('pcoumarylalcohol', search_db=False, formula="C9H10O2", MW=150.1772, Hf=-193500, phase='l')                 		#Unverified
species['formicacid'] = bst.Chemical('formicacid', search_db=False, formula="CH2O2", MW=46.0254, Hf=-431650, phase='l')                           			#ref: Heat of formation of liquid. Should it be of gas?
species['methanol'] = bst.Chemical('methanol', search_ID='methanol', phase='l')                               												#ref2
species['water'] = bst.Chemical('water', search_ID='water', phase='l')                                   													#ref2

species['propaniocacid'] = bst.Chemical('propaniocacid', search_db=False, formula="C3H6O2", Hf=-199481.24962000002, phase='l')
species['linalylpropioinate'] = bst.Chemical('linalylpropioinate', search_db=False, formula="C19H34O2", Hf=-194296.922694, phase='l')
species['vanillin'] = bst.Chemical('vanillin', search_ID='vanillin', phase='l')
species['cresol'] = bst.Chemical('cresol', search_db=False, formula="C7H8O", Hf=-34212.624394, phase='l')

# Create an empty Chemicals object and add each species to it
chemicals = tmo.Chemicals([])

for key, value in species.items():
    #print(key)
    value.default() # call the default method to populate all other properties not defined explicitly. Todo: some properties default to water's properties. Check individual chemical properties for accuracy once model is developed. 
    chemicals.append(value) # push a chemical object to back of the chemicals collection object

chemicals.compile()								# compile the loaded chemicals. Must do.
chemicals.set_synonym('water', 'moisture')		# set synonyms or common name for select species
chemicals.set_synonym('CaO', 'ash')

settings.set_thermo(chemicals, cache=False)				# Todo: read more set_thermo method.


# set raw material price (USD/kg)
price = {'biomass': 0.08, # @10% moisture content wet basis
         'water': 0.000353,
         'steam': 0.017,
         'waste': -0.33,
         'dry feedstock': 7.212,
         'gasoline': 0.756} # 2 USD/gal

#setting electructy price
settings.electricity_price = 0.1           # set electricity price to your desired rate

# set the name of your flowsheet
bst.main_flowsheet.set_flowsheet('fast_pyrolysis_biorefinery')

# define a chemical group with its composition in mass fraction here
chemicals.define_group('biomass',IDs = ['cellulose','xyhw','ligc','ligh','ligO','triglycerides','tannin','moisture','ash'], composition = [0.2347,0.1857,0.0915,0.0761,0.1290,0.0386,0.0176,0.10,0.0915], wt=True)


S101 = Stream('S101', biomass=333334.2, units='kg/hr', phase='s', price=price['biomass'])



# Vibrating Screen
S105 = Stream('S105', phase='s')
S106 = Stream('S106', phase='s')

UO104 = VibratingScreen('UO104', ins=S101, outs=(S105,S106), split=0.9)                     # Vibrating Screen


# Simulation of UO103 (Vibrating Screen)
UO104.simulate()
UO104.show(flow='kg/hr', composition = False, N=9)
r4 = UO104.results()
print(r4)
                                                                                           # Rotary Dryer
S107 = Stream('S107', phase='g') #Air
S108 = Stream('S108', phase='g') #Natural gas
S109 = Stream('S109', phase='s') #dry feedstock/solids
S110 = Stream('S110', phase='g') #hot air
S111 = Stream('S111', phase='g') #emissions
UO105 = units.DrumDryer('UO105', ins=(S105,S107,S108), outs=(S109,S110,S111),moisture_content = 0.01, split = (dict(moisture=0.9)) )
# Simulation of UO105 (Rotary dryer)
UO105.simulate()
UO105.show(flow='kg/hr', composition = False, N=10)
r5 = UO105.results()
print(r5)

S113 = Stream('S113', phase='s') #recycle stream
S114 = Stream('S114', phase='s', price=price['dry feedstock'])
UO106 = units.Mixer('UO106', ins=(S109, S113), outs=(S114))

# Hammermill
S112 = Stream('S112', phase='s') #crushed dry feedstock
S115 = Stream('S115', phase='s') 

UO107 = units.HammerMill('UO107', ins=(S114), outs=(S115))

#Vibrating Screen
S116 = Stream('S116', phase='s')
UO108 = VibratingScreen('UO108', ins=S115, outs=(S113,S116), split=0.1) 

Preprocessing = bst.main_flowsheet.create_system('Preprocessing')
Preprocessing.simulate()
Preprocessing.diagram('cluster', format='png')

UO106.show(flow='kg/hr', composition = False, N=9)
r6 = UO106.results()
print(r6)

UO107.show(flow='kg/hr', composition = False, N=9)
r7 = UO107.results()
print(r7)

UO108.show(flow='kg/hr', composition = False, N=9)
r8 = UO108.results()
print(r8)

class Preprocessing(bst.TEA):

    def __init__(self, system, IRR, duration, depreciation, income_tax,
                 operating_days, lang_factor, construction_schedule, WC_over_FCI,
                 labor_cost, fringe_benefits, property_tax,
                 property_insurance, supplies, maintenance, administration):
        # Huang et. al. does not take into account financing or startup
        # so these parameters are 0 by default
        super().__init__(system, IRR, duration, depreciation, income_tax,
                         operating_days, lang_factor, construction_schedule,
                         startup_months=0, startup_FOCfrac=0, startup_VOCfrac=0,
                         startup_salesfrac=0, finance_interest=0, finance_years=0,
                         finance_fraction=0, WC_over_FCI=WC_over_FCI)
        self.labor_cost = labor_cost
        self.fringe_benefits = fringe_benefits
        self.property_tax = property_tax
        self.property_insurance = property_insurance
        self.supplies= supplies
        self.maintenance = maintenance
        self.administration = administration


    def _DPI(self, installed_equipment_cost):
        return installed_equipment_cost


    def _TDC(self, DPI):
        return DPI

 
    def _FCI(self, TDC):
        return TDC

    
    def _FOC(self, FCI):
            return (FCI*(self.property_tax + self.property_insurance
                     + self.maintenance + self.administration)
                + self.labor_cost*(1+self.fringe_benefits+self.supplies))



sys = bst.main_flowsheet.create_system('sys')
tea =  Preprocessing(system=sys,
                  IRR=0.15,
                  duration=(2018, 2038),
                  depreciation='MACRS7',
                  income_tax=0.21, # Previously 35% in published study
                  operating_days=200,
                  lang_factor=3,
                  construction_schedule=(0.4, 0.6),
                  WC_over_FCI=0.05,
                  labor_cost=2.5e6,
                  fringe_benefits=0.4,
                  property_tax=0.11,
                  property_insurance=0.11,
                  supplies=0.20,
                  maintenance=0.01,
                  administration=0.005)


# # sys = main_flowsheet.stream
# sys.load_inlet_ports([S101], names={'feedstock': 0})
# sys.load_outlet_ports([S114], names={'product': 0})
# sys.show(data=False)
# sys.products[0]


tea.show()
# tea.get_cashflow_table()

feed = bst.main_flowsheet('S101')
products = [bst.main_flowsheet('S114')]
costs = tea.production_costs(products)# USD/yr
np.round(costs / 1e6) # million USD / yr


price = tea.solve_price(feed) # USD/kg
round(price, 6)

tea.IRR = tea.solve_IRR()
tea.show()
tea.get_cashflow_table()
#tea.save_report('report.xlsx')

#LCA

system=sys
feedstock=UO104.ins[0]
F_mass=feedstock.F_mass
GWP = 'GWP 100yr'
bst.settings.define_impact_indicator(key=GWP, units='kg*CO2e')
S101.set_CF(GWP, 0, basis='kg', units='kg*CO2e')
bst.settings.set_electricity_CF(GWP, 0, basis='kWhr', units='kg*CO2e')


GWP_preprocessing_displacement = sys.get_net_impact(key=GWP) / sys.get_mass_flow(S114)
#problem with this line
GWP_per_GGE = sys.get_property_allocated_impact(
    key=GWP, name='energy', basis='GGE', # Energy basis defaults to 'kJ'
    ignored={S106, S111}
)

GWP_preprocessing_energy = (
    GWP_per_GGE * S114.get_property('LHV', 'GGE/hr') / S114.F_mass
)

# GWP_per_USD = sys.get_property_allocated_impact(
#     key=GWP, name='revenue', basis='USD', # Revenue basis defaults to USD
# ) # kg-CO2e / USD

for s in sys.products: s.price = 50000000

GWP_per_USD = sys.get_property_allocated_impact(
    key=GWP, name='cost', basis='USD', # Revenue basis defaults to USD
) # kg-CO2e / USD

GWP_preprocessing_revenue = (
    GWP_per_USD * S114.price
) # kg-CO2e / kg-ethanol
print(
    "Cradle to biorefinery gate GWP of Preprocessing:\n "
    f"Displacement allocation: {GWP_preprocessing_displacement:.12f} [kg CO2e / kg]\n "
    f"Energy allocation: {GWP_preprocessing_energy:.12f} [kg CO2e / kg]\n "
    f"Economic allocation: {GWP_preprocessing_revenue:.2f} [kg CO2e / kg]"
)

Regards,
Uzair

[yoelcortes]

@Uzair1shah,

Yes, use the <Stream>.set_CF method to set the characterization factor of the natural gas inlet to the drum drier (<DrumDrier>.ins[2].set_CF). Don't forget to set characterization factors for your steam (using bst.settings.set_utility_agent_CF )if you are not using the BoilerTurbogenerator.

Thanks,