deng-mit/arrhenius.jl

Differentiable Reacting Flow Modeling Software

arrhenius.jl's Issues

Use NN to solve PSR

Use a NN to produce S-Curve

Determine the choice of sparse matrix based on the matrix size

Need a careful benchmarking to optimize the choice of whether using sparse matrix according to the size of the matrix and its sparsity

Set up documentation

https://sparkbox.com/foundry/github_wiki_tutorial_for_technical_wiki_documentation

Fail to load CH4_Kazakov_without_NOx.yaml in PSR.jl

hi,

I try just to compile the PSR.jl resulting in the error below. Can you recompile the DRM22 mech (CH4_Kazakov_without_NOx.yaml) / or fix it ?

Thanks, Dmitry

julia> gas = CreateSolution("CH4_Kazakov_without_NOx.yaml")
ERROR: MethodError: no method matching Array{Float64,2}(::Array{Float64,3})
Closest candidates are:
Array{Float64,2}(::AbstractArray{S,N}) where {T, N, S} at array.jl:562
Array{Float64,2}(::UndefInitializer, ::Int64, ::Int64) where T at boot.jl:408
Array{Float64,2}(::UndefInitializer, ::Int64...) where {T, N} at boot.jl:412
...
Stacktrace:
[1] convert(::Type{Array{Float64,2}}, ::Array{Float64,3}) at .\array.jl:554
[2] Arrhenius.Transport(::Int64, ::Array{Float64,2}, ::Array{Float64,2}, ::Array{Float64,3}) at C:\Users\lysenkod.julia\packages\Arrhenius\BOKjU\src\DataStructure.jl:23
[3] CreateSolution(::String) at C:\Users\lysenkod.julia\packages\Arrhenius\BOKjU\src\Solution.jl:142
[4] top-level scope at REPL[31]:1

Implement a PSR in the example folder

For example
https://github.com/DENG-MIT/NN-PSR/blob/main/node_pyrolysis.jl

Interpret kinetic data from yaml directly and bypass .npz

Arrhenius.jl/src/Solution.jl

Line 56 in 51a2f3d

product_stoich_coeffs = sparse(npz["product_stoich_coeffs"])

Basically, following the code in ReacTorch. So far, it works well, and this task is not urgent at all.

Investigate if it is beneficial to use sparse matrix for stoichiometric coeffcients

Arrhenius.jl/src/Solution.jl

Line 56 in 51a2f3d

product_stoich_coeffs = sparse(npz["product_stoich_coeffs"])

Perhaps testing one small mechanism of GRI3.0, and a large one for n-heptane.

Implement an equilibrium solver

This is very essential for a lot of the cases

Arrhenius.jl for Rapid Compression Machine

Code is almost ready, research can come online soon

Prepare a list of resources for learning Julia

As Julia is not widely adopted in the Combustion community, it will be great to have instructions for beginners from different backgrounds. For example,

• without any programming experience
• knowing Fortran/Chemkin well
• knowing Cantera's python/Matlab interface well

Wrong definition of X2C

The current equation for X2C(gas, X, ρ_mass)

C = X * ρ_mass ./ gas.MW

seems incorrect to me. It is the same as for Y2C, which can't be the case.
Instead it should read

C = X * ρ_mass / sum(X.*gas.MW)

Arrhenius.jl for shock tubes

Already done

Divide the interpretation of kinetic model into three parts

Arrhenius.jl/src/Solution.jl

Line 9 in 51a2f3d

function CreateSolution(mech)

Currently, the interpretation of kinetic, thermodynamic, and transport are coded together. It will be great to spill them into three parts, to make the code more modulized.

Determine if we need those @inbounds

If those @inbounds are not necessary, we shall remove them to improve the readiness of the source code.

Arrhenius.jl/src/Kinetics.jl

Lines 2 to 53 in 03e31ea

	function wdot_func(reaction, T, C, S0, h_mole; get_qdot=false)
	@inbounds _kf = @. @view(reaction.Arrhenius_coeffs[:, 1]) * exp(
	@view(reaction.Arrhenius_coeffs[:, 2]) * log(T) -
	@view(reaction.Arrhenius_coeffs[:, 3]) * (4184.0 / R / T),
	)
	for i in reaction.index_three_body
	@inbounds _kf[i] *= dot(@view(reaction.efficiencies_coeffs[:, i]), C)
	end
	for (j, i) in enumerate(reaction.index_falloff)
	@inbounds A0, b0, Ea0 = reaction.Arrhenius_0[j, :]
	@inbounds k0 = A0 * exp(b0 * log(T) - Ea0 * 4184.0 / R / T)
	@inbounds Pr =
	k0 * dot(@view(reaction.efficiencies_coeffs[:, i]), C) / _kf[i]
	lPr = log10(Pr)
	_kf[i] *= (Pr / (1 + Pr))
	# reference:
	# http://web.mit.edu/2.62/cantera/doc/html/classCantera_1_1Troe4.html#a38aa787421d426dfd0a587fd6fc8108e
	if reaction.index_falloff_Troe[j] > 0
	k = reaction.index_falloff_Troe[j]
	@inbounds F_cent =
	(1 - reaction.Troe_[k, 1]) * exp(-T / reaction.Troe_[k, 4]) +
	reaction.Troe_[k, 1] * exp(-T / reaction.Troe_[k, 2]) +
	exp(-reaction.Troe_[k, 3] / T)
	lF_cent = log10(F_cent)
	_C = -0.4 - 0.67 * lF_cent
	N = 0.75 - 1.27 * lF_cent
	@inbounds f1 = (lPr + _C) / (N - 0.14 * (lPr + _C))
	@inbounds _kf[i] *= exp(log(10.0) * lF_cent / (1 + f1^2))
	end
	end
	@inbounds ΔS_R = reaction.vk' * S0 / R
	@inbounds ΔH_RT = reaction.vk' * h_mole / (R * T)
	@inbounds Keq =
	@. exp(ΔS_R - ΔH_RT + log(one_atm / R / T) * reaction.vk_sum)
	@inbounds _kr = @. _kf / Keq * reaction.is_reversible
	for i = 1:reaction.n_reactions
	@inbounds for j in reaction.i_reactant[i]
	@inbounds _kf[i] *= C[j]^reaction.reactant_orders[j, i]
	end
	if reaction.is_reversible[i]
	@inbounds for j in reaction.i_product[i]
	@inbounds _kr[i] *= C[j]^reaction.product_stoich_coeffs[j, i]
	end
	end
	end

Ideally, we could come up with a guideline on when we need explicitly to use the keyword, and when the compiler will automatically do the code optimization.

Implement a multiplier to rate constants

Pass a vector that contains factors to create the solution, and multiply it in wdot_func

Handle falloff reaction without troe coefficients

Arrhenius.jl/src/Kinetics.jl

Line 23 in cb0b2ee

if reaction.Troe_A[j] > 1.e-12

Need better approach to handle such situations. Currently, I have to manually go to the yaml file and put zero Troe coefficients.

Discretization for one-dimensional flame

Goal:

implement the discretization for one-dimensional flame, then one can compute the residuals and the sensitivity for a given solution. Further, one can also solve the equations.

Equations:

see Chemkin manual Sec. 12.6

Notes:

similar to the one did for sensBVP for ignition, but more complex with the central difference.

TagBot trigger issue

This issue is used to trigger TagBot; feel free to unsubscribe.

If you haven't already, you should update your TagBot.yml to include issue comment triggers.
Please see this post on Discourse for instructions and more details.

If you'd like for me to do this for you, comment TagBot fix on this issue.
I'll open a PR within a few hours, please be patient!

Determine the tolerance for combustion simulations

References

• ReactionMechanismGenerator/ReactionMechanismSimulator.jl#87

add elemental matrix into struct gas

something like

E_C = [10, 0, 0, 1, 1, 3, 2, 3, 5, 6, 7, 0, 0, 0, 0];
E_H = [16, 1, 2, 3, 4, 5, 4, 6, 6, 6, 8, 0, 0, 1, 2];
E_N = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0];
E_O = [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 2, 2, 2];

This is especially useful for computing elemental conversation subspace.

Add Unit Testing for development

As the package becomes more and more complex, it is important to adopt unit testing. We can following the example of RMS.jl at https://github.com/ReactionMechanismGenerator/ReactionMechanismSimulator.jl/blob/ac2a3bd88504b86d5bc60dedc09533839a9c9cbb/src/Calculators/TestThermo.jl#L1-L60.

Some thoughts:

• We could put the test files inside the src folder, with a consistent prefix, so that one can easily spot the actual source file.
• We can test the functions against Cantera via PyCall. This should work well for most of the functionality.
• @test Cps ≈ Cpexplist rtol=1e-3 where 1e-3 is a good default value for tolerance.