diff --git a/cmlm/utils/input_file.py b/cmlm/utils/input_file.py index 695055f..c87a2bf 100644 --- a/cmlm/utils/input_file.py +++ b/cmlm/utils/input_file.py @@ -8,6 +8,14 @@ import tomlkit +def scalar_to_list(val): + """Leave lists as is, convert scalars to 1 element lists.""" + if not isinstance(val, (list, tuple)): + return [val] + else: + return val + + def recursively_update_dict(base, new): """ Update items and subitems in one nested dict-like object based on another. @@ -143,6 +151,9 @@ def __init__( self.base = base self.output = kwargs.get("output", None) if self.output is not None: + # if output ends with .toml, its a file, otherwise its a directory + if not self.output.endswith(".toml"): + self.output = os.path.join(self.output, "config.toml") self.output_dir = os.path.split(self.output)[0] if self.base is None: if len(self.output_dir) > 0 and not os.path.exists(self.output_dir): @@ -202,7 +213,7 @@ def parse_file(cls, file_name=None, additional_args=None, **kwargs): return cls(data, name=name, base=None, **kwargs) @classmethod - def parse_args(cls, description=None, infile=None): + def parse_args(cls, description=None, infile=None, require_output=False): """ Parse command line arguments specifying file and arguments to create a TPP. @@ -212,6 +223,8 @@ def parse_args(cls, description=None, infile=None): Short description of program for which config is being loaded infile: str, optional Default TOML input file to use + require_output: bool, optional + if true, the `-o` command line argument is required. Default False. Returns ------- @@ -235,7 +248,8 @@ def parse_args(cls, description=None, infile=None): "-o", "--output", default=None, - help="File in which to write used inputs/outputs", + required=require_output, + help="File (.toml) or directory in which to write used inputs/outputs", ) parser.add_argument( "-t", @@ -346,7 +360,7 @@ def __setitem__(self, item_name, value): self.data[prefix] = tomlkit.document() self[prefix][suffix] = value - def get(self, item_name, default=None, doc=None): + def get(self, item_name, default=None, doc=None, choices=None): """ Retrieve a leaf or subtable form the TomlParmParse table. @@ -363,6 +377,9 @@ def get(self, item_name, default=None, doc=None): value to use if item_name is not found in table doc: optional string to add as a comment in the TOML file + choices: optional + list or tuple of allowable options for input parameter. If specified, + an error will be raised if the specified value is not in the list. Returns ------- @@ -381,9 +398,19 @@ def get(self, item_name, default=None, doc=None): ) self[item_name] = retval + if choices is not None: + if retval not in choices: + raise ValueError( + f"In TomlParmParse object {self.name}:\n" + f" Invalid value specified for item <{item_name}> (doc: {doc})\n" + f" Choices are: {choices}" + ) + if doc is not None and self.output_type == "doc": if default is not None: doc += f" | optional, default: {default}" + if choices is not None: + doc += f" | choices are: {choices}" self[item_name].comment(doc) return retval diff --git a/docs/conf.py b/docs/conf.py index 48a4811..1b22405 100644 --- a/docs/conf.py +++ b/docs/conf.py @@ -6,6 +6,7 @@ import sys sys.path.insert(0, os.path.abspath("../run_scripts/ctable")) +sys.path.insert(0, os.path.abspath("../run_scripts/cantera")) # -- Project information ----------------------------------------------------- diff --git a/docs/run-scripts.rst b/docs/run-scripts.rst index a0f0da5..e741c3e 100644 --- a/docs/run-scripts.rst +++ b/docs/run-scripts.rst @@ -3,12 +3,60 @@ Run Scripts =========== +These scripts all use the TomlParmParse utility to manage configuration/inputs. Sample +inputs with descriptions are included in ``.toml`` files with the same base file name +as the script. All scripts can be run using:: + + python script_name.py script_name.toml + +Adding the ``-h`` flag will give the additional command line inputs that are associated +with the TomlParmParse utility:: + + positional arguments: + infile Input file to parse inputs from + + options: + -h, --help show this help message and exit + -o OUTPUT, --output OUTPUT + File (.toml) or directory in which to write used inputs/outputs + -t OUTPUT_TYPE, --output_type OUTPUT_TYPE + Type of output file: *doc*: include comments documenting used inputs *clean*: no comments, only used inputs *original*: input file in original formatting + -a ARGS, --args ARGS Override arguments, as a toml string + -w, --no_overwrite Disallow overwriting entries once they have been used or set + -e, --error_unused Raise error if there are unused inputs + -l, --live_update Update output file continuously as code runs + +Notably, specifying ``-o`` will result in a new toml file with the configuration that +was used, including documentation for each input by default, and ``-a`` allows inputs +that are usually specified in the input file to be specified on the command line. + +---- + +cantera +------- + +These python scripts use Cantera to generate 1D flame solutions that will then be +compiled into tabulated chemistry tables. + +---- + +compute_premixed_flamelets.py +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + +.. automodule:: + compute_premixed_flamelets + +---- + +ctable +------ + These python scripts use the CMLM package to generate tabulated chemistry data. ---- create_dummy_table_nd.py ------------------------- +~~~~~~~~~~~~~~~~~~~~~~~~ .. automodule:: create_dummy_table_nd @@ -16,7 +64,7 @@ create_dummy_table_nd.py ---- create_spray_table_nd.py ------------------------- +~~~~~~~~~~~~~~~~~~~~~~~~ .. automodule:: create_spray_table_nd diff --git a/run_scripts/cantera/compute_premixed_flamelets.py b/run_scripts/cantera/compute_premixed_flamelets.py new file mode 100644 index 0000000..42989c0 --- /dev/null +++ b/run_scripts/cantera/compute_premixed_flamelets.py @@ -0,0 +1,239 @@ +""" +Compute 1D unstrained premixed flamelets using Cantera. + +Flames can be computed across a range of conditions including +pressures, temperatures, and compositions. Three methods are +possible for specifying composition sweeps: equivalence ratio, +mixture fraction, and directly specifying a composition. + +Uses TOML format input files (see example input file for more details) +and requires specification of an output directory on the command line:: + + python compute_premixed_flamelets.py compute_premixed_flamelets.toml -o output +""" + +if __name__ == "__main__": + + import itertools + import os + + import cantera as ct + import numpy as np + import pandas as pd + + from cmlm.utils import TomlParmParse + from cmlm.utils.input_file import scalar_to_list + + # ------------------ Parse relevant inputs ----------------------------# + pp = TomlParmParse.parse_args( + description="A tool to compute planar premixed flames using Cantera", + require_output=True, + ) + outdir = pp.output_dir + + # Conditions + ppc = pp["conditions"].doc("Physical conditions for flames") + pressures = scalar_to_list(ppc.get("pressures", doc="Single pressure or list, atm")) + + composition_type = ppc.get( + "composition_type", + choices=["zmixs", "phis", "single"], + doc="Method of specifying compositions to run", + ) + if composition_type == "zmixs": + fuel_comp = ppc.get( + "fuel_comp", doc="Fuel stream cantera composition, mass basis" + ) + fuel_temp = ppc.get("fuel_temp", doc="K") + oxid_comp = ppc.get( + "oxid_comp", doc="Oxidizer stream cantera composition, mass basis" + ) + oxid_temp = ppc.get("oxid_temp", doc="K") + Zvalues = scalar_to_list(ppc.get("Zvalues", doc="Single Z value or list")) + cond_labels = ["p{:.4f}", "zmix{:.4f}"] # noqa : FS003 + cond_iterator_global = list(itertools.product(pressures, Zvalues)) + elif composition_type == "phis": + fuel_comp = ppc.get( + "fuel_comp", doc="Fuel stream cantera composition, mass basis" + ) + oxid_comp = ppc.get( + "oxid_comp", doc="Oxidizer stream cantera composition, mass basis" + ) + temperatures = scalar_to_list( + ppc.get("temperatures", doc="Single temperature or list, K") + ) + phis = scalar_to_list(ppc.get("phis", doc="Single equivalence ratio or list")) + cond_labels = ["p{:.4f}", "T{:.1f}", "phi{:.4f}"] # noqa : FS003 + cond_iterator_global = list(itertools.product(pressures, temperatures, phis)) + else: + comp = ppc.get("comp", doc="cantera composition, mass basis") + temperatures = scalar_to_list( + ppc.get("temperatures", doc="Single temperature or list, K") + ) + cond_labels = ["p{:.4f}", "T{:.1f}"] # noqa : FS003 + cond_iterator_global = list(itertools.product(pressures, temperatures)) + flame_width = ppc.get( + "dom_width", default=0.1, doc="Domain width for flame simulations (m)" + ) + + # Models + ppm = pp["models"].doc("Physical models used by solver") + mechanism = ppm.get( + "mechanism", doc="Chemical mechanism to use (Cantera yaml format)" + ) + transport = ppm.get( + "transport", + default="mixture-averaged", + choices=[ + "mixture-averaged", + "unity-Lewis-number", + "multicomponent", + "high-pressure", + "high-pressure-Chung", + ], + doc="Transport property model from Cantera to use, " + "see Cantera documentation for more information", + ) + eos = ppm.get( + "eos", + default="", + doc="Equation of State, must be a phase option in your mechanism yaml file. " + "If unspecified, the default phase from the mechanism file is used", + ) + + # Flame numerics + ppf = pp["numerics"].doc( + "Numerics for premixed flame solve, " + "see Cantera documentation for more information" + ) + loglevel = ppf.get("loglevel", default=0) + ratio = ppf.get("ratio", default=2) + slope = ppf.get("slope", default=0.2) + curve = ppf.get("curve", default=0.2) + prune = ppf.get("prune", default=0.1) + max_points = ppf.get("max_points", default=10000) + + # Options + ppo = pp["options"].doc("Additional options for script") + use_mpi = ppo.get("use_mpi", default=False, doc="Using MPI requires mpi4py") + if use_mpi: + from mpi4py import MPI + + comm = MPI.COMM_WORLD + rank = comm.Get_rank() + nprocs = comm.Get_size() + else: + rank = 0 + nprocs = 1 + metadata_file = ppo.get( + "metadata_file", + default="", + doc="If specified, save chemtable metadata to this file", + ) + cp_fuel_species = ppo.get( + "cp_fuel_species", + default="", + doc="If specified, Cantera composition for fuel to compute cp_fuel in " + "flame output, or 'fuel_comp' to use that input if available", + ) + if cp_fuel_species == "fuel_comp" and composition_type != "single": + cp_fuel_species = fuel_comp + + # ------------------ Solve Flames ---------------------------- # + cond_iterator = cond_iterator_global[rank::nprocs] + gas = ct.Solution(mechanism, eos) + + output = pd.DataFrame( + data=np.nan, + index=pd.MultiIndex.from_tuples( + cond_iterator, names=[lab[: lab.index("{")] for lab in cond_labels] + ), + columns=["s_L", "T_ad", "l_f", "Nx", "dx_min"], + ) + for cond in cond_iterator: + label = "_".join(condlabel.format(c) for condlabel, c in zip(cond_labels, cond)) + print(f"Rank {rank} - Computing flame: {label}", flush=True) + + # Set up conditions + if composition_type == "zmixs": + gas.TPY = oxid_temp, cond[0] * ct.one_atm, oxid_comp + oxstream = ct.Quantity(gas, constant="HP") + oxstream.mass = 1.0 - cond[1] + gas.TPY = fuel_temp, cond[0] * ct.one_atm, fuel_comp + fustream = ct.Quantity(gas, constant="HP") + fustream.mass = cond[1] + mix = fustream + oxstream + gas.TPY = mix.T, mix.P, mix.Y + elif composition_type == "phis": + gas.TP = cond[1], cond[0] * ct.one_atm + gas.set_equivalence_ratio(cond[2], fuel_comp, oxid_comp, basis="mass") + else: + gas.TPY = cond[1], cond[0] * ct.one_atm, comp + + # Solve Flame + flame = ct.FreeFlame(gas, width=flame_width) + flame.set_refine_criteria(ratio=ratio, slope=slope, curve=curve, prune=prune) + flame.set_max_grid_points(1, max_points) + flame.transport_model = transport + flame.solve(loglevel=loglevel) + + # Extract outputs from flame + max_dtempdx = np.max(np.diff(flame.T) / np.diff(flame.grid)) + flame_thickness = (np.max(flame.T) - np.min(flame.T)) / max_dtempdx + flame_temp = flame.T[-1] + flame_speed = flame.velocity[0] + flame_grid = len(flame.T) + dx_min = np.min(np.diff(flame.grid)) + output.loc[cond] = flame_speed, flame_temp, flame_thickness, flame_grid, dx_min + + # Save flame solution - default Cantera MKS units + data = pd.DataFrame() + data["X"] = flame.grid + data["T"] = flame.T + data["VEL"] = flame.velocity + data["RHO"] = flame.density_mass + data["DIFF"] = flame.thermal_conductivity / flame.cp_mass + data["VISC"] = flame.viscosity + data["LAMBDA"] = flame.thermal_conductivity + data["CP"] = flame.cp_mass + data["MW"] = flame.mean_molecular_weight + if cp_fuel_species != "": + sa = flame.to_array() + sa.TPY = sa.T, sa.P, cp_fuel_species + data["CP_FUEL"] = sa.cp_mass + spec_names = [f"Y-{spec}" for spec in gas.species_names] + spec_y_data = pd.DataFrame(flame.Y.T, columns=spec_names) + rr_names = [f"SRC_{spec}" for spec in gas.species_names] + spec_rr_data = pd.DataFrame( + flame.net_production_rates.T * list(gas.molecular_weights), columns=rr_names + ) + data = pd.concat([data, spec_y_data, spec_rr_data], axis=1) + data.to_csv(os.path.join(outdir, f"flame_{label}.csv")) + + # We're finished with this flame + print( + f"Rank {rank} - Finished flame: {label}. " + f"sL={flame_speed:7.4f} Tad={flame_temp:7.1f} l_f={flame_thickness:10.3e} " + f"N={flame_grid:5n} dx_min={dx_min:10.3e}", + flush=True, + ) + + # Collect data and save + print(f"Rank {rank} - Completed all required tasks") + if use_mpi: + all_output = comm.gather(output) + else: + all_output = [output] + if rank == 0: + all_output = pd.concat(all_output).sort_index() + all_output.to_csv(os.path.join(outdir, "flamestats.csv")) + print("All Ranks Completed.") + print(all_output) + if metadata_file != "": + with open(metadata_file, "w") as fi: + fi.write("manifold.has_species_mw = true\n") + for i, spec in enumerate(gas.species_names): + fi.write(f"manifold.{spec}_mw = {gas.molecular_weights[i]}\n") + if len(pressures) != 1: + raise RuntimeError("Can only save metadata for a single pressure") + fi.write(f"manifold.nominal_pressure_cgs = {pressures[0] * 10.0}\n") diff --git a/run_scripts/cantera/compute_premixed_flamelets.toml b/run_scripts/cantera/compute_premixed_flamelets.toml new file mode 100644 index 0000000..d299104 --- /dev/null +++ b/run_scripts/cantera/compute_premixed_flamelets.toml @@ -0,0 +1,75 @@ +# ==================================================================================== # +# Physical conditions for flames +# ==================================================================================== # +[conditions] + +# Single pressure or list, atm +pressures = 1.0 + +# Method of specifying compositions to run | choices are: ['zmixs', 'phis', 'single'] +composition_type = "phis" + +# Fuel stream cantera composition, mass basis +fuel_comp = "H2:1" + +# Oxidizer stream cantera composition, mass basis +oxid_comp = "O2:1, N2:4" + +# Single temperature or list, K +temperatures = [300.0, 400.0] + +# Single equivalence ratio or list +phis = [1.0,1.1,1.2,1.3] + +# Domain width for flame simulations (m) | optional, default: 0.1 +dom_width = 0.1 + + + +# ==================================================================================== # +# Physical models used by solver +# ==================================================================================== # +[models] + +# Chemical mechanism to use (Cantera yaml format) +mechanism = 'h2o2.yaml' + +# Transport property model from Cantera to use, see Cantera documentation for more +# information | optional, default: mixture-averaged | choices are: ['mixture- +# averaged', 'unity-Lewis-number', 'multicomponent', 'high-pressure', 'high-pressure- +# Chung'] +transport = 'unity-Lewis-number' + +# Equation of State, must be a phase option in your mechanism yaml file. If +# unspecified, the default phase from the mechanism file is used | optional, default: +eos = "" + + + +# ==================================================================================== # +# Numerics for premixed flame solve, see Cantera documentation for more information +# ==================================================================================== # +[numerics] + +loglevel = 0 +ratio = 2 +slope = 0.05 +curve = 0.05 +prune = 0.02 +max_points = 10000 + + +# ==================================================================================== # +# Additional options for script +# ==================================================================================== # +[options] + +# Using MPI requires mpi4py | optional, default: False +use_mpi = true + +# If specified, save chemtable metadata to this file | optional, default: +metadata_file = "" + +# If specified, Cantera composition for fuel to compute cp_fuel in flame output, or +# 'fuel_comp' to use that input if available | optional, default: +cp_fuel_species = "" diff --git a/run_scripts/ctable/create_dummy_table_nd.py b/run_scripts/ctable/create_dummy_table_nd.py index 6d4cc35..7254378 100644 --- a/run_scripts/ctable/create_dummy_table_nd.py +++ b/run_scripts/ctable/create_dummy_table_nd.py @@ -1,22 +1,10 @@ """ Create sample chemtables with dummy data for testing purposes. -Usage ------ +Uses TOML format input files (see example input file for more details):: -Invoke on the command line:: + python create_dummy_table_nd.py create_dummy_table_nd.toml - python create_dummy_table_nd.py - -Input File ----------- - - The input file is a TOML format file with single section ``[table]`` that - specifies the following: - - - ``ndim`` (int): number of dimensions in table - - ``ngrid`` (int): number of grid points in each dimension - - ``outfile_prefix`` (string): path to save files """ if __name__ == "__main__": diff --git a/run_scripts/ctable/create_spray_table_nd.py b/run_scripts/ctable/create_spray_table_nd.py index 586b50a..89bd645 100644 --- a/run_scripts/ctable/create_spray_table_nd.py +++ b/run_scripts/ctable/create_spray_table_nd.py @@ -15,23 +15,9 @@ mechanism. The `liquid_fuels_nonreacting` mechanism that comes with PelePhysics is likely a good choice. -Usage ------ +Uses TOML format input files (see example input file for more details):: -Invoke on the command line:: - - python create_spray_table_nd.py - -To see additional runtime options you can run:: - - python create_spray_table_nd.py -h - -Input File ----------- - -The input file is a TOML format file with sections ``[phys]`` and ``[table]``, -corresponding to physical/BC inputs and parameters for the table. See the -sample input file for more details. + python create_spray_table_nd.py create_spray_table_nd.toml """ if __name__ == "__main__": @@ -58,7 +44,8 @@ pressure = ppp.get("pressure", doc="ambient pressure, Pa") liq_temp_fuel = ppp.get("liq_temp_fuel", doc="liquid temps for each fuel, K") X_fuel = ppp.get("X_fuel", doc="Cantera composition string") - species_list = ["O2"] + [spec.split(":")[0] for spec in X_fuel] + fuel_species_list = [spec.split(":")[0] for spec in X_fuel] + species_list = ["O2"] + fuel_species_list delta_h_vap = ppp.get("delta_h_vap", doc="Latent heats for each fuel, J/kg") ox = ct.Solution(mechanism) @@ -157,7 +144,7 @@ ", ".join( [ f"{spec}:{max(mixture.Y[mixture.species_index(spec)],eps)}" - for spec in species_list + for spec in fuel_species_list ] ), )