Skip to content

Commit

Permalink
exclude species feature
Browse files Browse the repository at this point in the history
  • Loading branch information
Nicholaswogan committed Oct 17, 2024
1 parent fa8b32e commit 754c517
Show file tree
Hide file tree
Showing 4 changed files with 41 additions and 22 deletions.
12 changes: 8 additions & 4 deletions photochem/equilibrate.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,16 +3,18 @@
from .utils._format import FormatReactions_main, yaml, Loader, MyDumper, mechanism_dict_with_atoms
from photochem_clima_data import DATA_DIR

def generate_zahnle_earth_thermo(outfile='zahnle_earth_thermo.yaml', atoms_names=None):
def generate_zahnle_earth_thermo(outfile='zahnle_earth_thermo.yaml', atoms_names=None, exclude_species=[]):
"""Generates a thermodynamic file for equilibrium solving that includes
condensible species (e.g., H2O condensate).
Parameters
----------
outfile : str, optional
Name of the output file, by default 'zahnle_earth_thermo.yaml'
atoms_names: list, optional
atoms_names : list, optional
List of atoms to keep. By default all atoms in the mechanism are kept
exclude_species : list, optional
List of species to exclude.
"""

rx_folder = DATA_DIR+'/reaction_mechanisms/'
Expand All @@ -32,8 +34,10 @@ def generate_zahnle_earth_thermo(outfile='zahnle_earth_thermo.yaml', atoms_names
for i,sp in enumerate(dat1['species']):
dat['species'].append(sp)

if atoms_names is not None:
dat = mechanism_dict_with_atoms(dat, atoms_names)
if atoms_names is None:
atoms_names = [a['name'] for a in dat['atoms']]

dat = mechanism_dict_with_atoms(dat, atoms_names, exclude_species)

dat = FormatReactions_main(dat)

Expand Down
7 changes: 5 additions & 2 deletions photochem/extensions/gasgiants.py
Original file line number Diff line number Diff line change
Expand Up @@ -955,6 +955,7 @@ def generate_photochem_rx_and_thermo_files(
atoms_names=['H','He','N','O','C','S'],
rxns_filename='photochem_rxns.yaml',
thermo_filename='photochem_thermo.yaml',
exclude_species=[],
remove_particles=False,
remove_reaction_particles=True
):
Expand All @@ -968,6 +969,8 @@ def generate_photochem_rx_and_thermo_files(
Name of output reactions file, by default 'photochem_rxns.yaml'
thermo_filename : str, optional
Name of output thermodynamic file, by default 'photochem_thermo.yaml'
exclude_species : list, optional
List of species to exclude.
remove_particles : bool, optional
If True, then particles will be removed, by default False.
remove_reaction_particles : bool, optional
Expand All @@ -976,11 +979,11 @@ def generate_photochem_rx_and_thermo_files(
# Kinetics
with open(zahnle_earth,'r') as f:
rxns = yaml.load(f,Loader=yaml.Loader)
rxns = mechanism_dict_with_atoms(rxns, atoms_names, remove_particles, remove_reaction_particles)
rxns = mechanism_dict_with_atoms(rxns, atoms_names, exclude_species, remove_particles, remove_reaction_particles)
rxns = FormatReactions_main(rxns)
with open(rxns_filename,'w') as f:
yaml.dump(rxns,f,Dumper=MyDumper,sort_keys=False,width=70)

# Thermodynamics
equilibrate.generate_zahnle_earth_thermo(thermo_filename, atoms_names)
equilibrate.generate_zahnle_earth_thermo(thermo_filename, atoms_names, exclude_species)

42 changes: 27 additions & 15 deletions photochem/utils/_format.py
Original file line number Diff line number Diff line change
Expand Up @@ -177,7 +177,12 @@ def FormatSettings_main(data):

return data

def mechanism_dict_with_atoms(dat, atoms_names, remove_particles=False, remove_reaction_particles=False):
def species_in_reaction(rx):
rx1 = rx.replace('<=>','=>').replace('(','').replace(')','')
react, prod = [a.replace(' ','').split('+') for a in rx1.split('=>')]
return react + prod

def mechanism_dict_with_atoms(dat, atoms_names, exclude_species=[], remove_particles=False, remove_reaction_particles=False):

atoms = []
exclude_atoms = []
Expand All @@ -199,6 +204,10 @@ def mechanism_dict_with_atoms(dat, atoms_names, remove_particles=False, remove_r
if tmp in exclude_atoms:
exclude = True
break

if sp['name'] in exclude_species:
exclude = True

if not exclude:
species.append(sp)

Expand All @@ -215,26 +224,24 @@ def mechanism_dict_with_atoms(dat, atoms_names, remove_particles=False, remove_r

if remove_reaction_particles and sp['formation'] == "reaction":
exclude = True

if sp['name'] in exclude_species:
exclude = True
if sp['formation'] == "saturation":
if sp['gas-phase'] in exclude_species:
exclude = True
elif not remove_reaction_particles and sp['formation'] == "reaction":
rx_sp = species_in_reaction(sp['equation'])
if any(a in exclude_species for a in rx_sp):
exclude = True

if not exclude:
particles.append(sp)

if "reactions" in dat:
reactions = []
for i,rx in enumerate(dat['reactions']):
eq = rx['equation']
eq = eq.replace('(','').replace(')','')
if '<=>' in eq:
split_str = '<=>'
else:
split_str = '=>'

a,b = eq.split(split_str)
a = a.split('+')
b = b.split('+')
a = [a1.strip() for a1 in a]
b = [b1.strip() for b1 in b]
sp = a + b
sp = species_in_reaction(rx['equation'])

exclude = False
for s in sp:
Expand All @@ -244,6 +251,10 @@ def mechanism_dict_with_atoms(dat, atoms_names, remove_particles=False, remove_r
break
if exclude:
break

if any(a in exclude_species for a in sp):
exclude = True

if not exclude:
reactions.append(rx)

Expand All @@ -257,14 +268,15 @@ def mechanism_dict_with_atoms(dat, atoms_names, remove_particles=False, remove_r

return out

def resave_mechanism_with_atoms(infile, outfile, atoms_names, remove_particles=False, remove_reaction_particles=False):
def resave_mechanism_with_atoms(infile, outfile, atoms_names, exclude_species=[], remove_particles=False, remove_reaction_particles=False):

with open(infile,'r') as f:
dat = yaml.load(f,Loader=Loader)

out = mechanism_dict_with_atoms(
dat,
atoms_names,
exclude_species,
remove_particles,
remove_reaction_particles
)
Expand Down
2 changes: 1 addition & 1 deletion setup.py
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,7 @@

setup(
name="photochem",
packages=['photochem','photochem.utils'],
packages=['photochem','photochem.utils','photochem.extensions'],
python_requires='>=3.6',
version=version,
license="GNU General Public License v3.0",
Expand Down

0 comments on commit 754c517

Please sign in to comment.