Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Radiation model update in StFlow.cpp #965

Closed
wants to merge 7 commits into from

Conversation

lavrenyukiv
Copy link

@lavrenyukiv lavrenyukiv commented Jan 25, 2021

Changes proposed in this pull request

  • Absorption coefficients for CO2, H2O, CO, CH4, OH molecules are added,
  • Spectrums from https://hitran.org/ are used for evaluated coefficients,
  • Polynomials are changed by linear interpolation in logarithmic coordinates between tabulated values,

@lavrenyukiv lavrenyukiv changed the title Update StFlow.cpp Radiation model update in StFlow.cpp Jan 25, 2021
@codecov
Copy link

codecov bot commented Jan 25, 2021

Codecov Report

Merging #965 (5db79c0) into main (e49d62b) will decrease coverage by 0.02%.
The diff coverage is 55.00%.

Impacted file tree graph

@@            Coverage Diff             @@
##             main     #965      +/-   ##
==========================================
- Coverage   71.26%   71.23%   -0.03%     
==========================================
  Files         377      377              
  Lines       46299    46371      +72     
==========================================
+ Hits        32993    33031      +38     
- Misses      13306    13340      +34     
Impacted Files Coverage Δ
src/oneD/StFlow.cpp 86.24% <55.00%> (-4.42%) ⬇️

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update e49d62b...89cc361. Read the comment docs.

Copy link
Member

@ischoegl ischoegl left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thanks for working on this. While others may have more comments, I have one observation to share.

@@ -94,9 +94,12 @@ StFlow::StFlow(ThermoPhase* ph, size_t nsp, size_t points) :
setupGrid(m_points, gr.data());

// Find indices for radiating species
m_kRadiating.resize(2, npos);
m_kRadiating.resize(5, npos);
Copy link
Member

@ischoegl ischoegl Jan 26, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

There has been some work to make the identification of species more flexible (which would mean that lowercase species names would be recognized). Rather than hard coding names in the subsequent lines, specifying the elemental composition will in most cases work (there usually is only one isomer). The function name is findIsomers

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I changed this part.

It would be great if the user could select the components involved in the radiation himself.

It would probably be better to load the absorption coefficient data from a separate file, so that the user can add his own dependencies.

Unfortunately, my experience with С++ is not enough.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I agree that it would be great to select species / import parameters instead of hard coding them. There are likely multiple ways of getting this done. I can probably give you pointers, but it’s up to @bryanwweber and @speth ...

Copy link
Member

@bryanwweber bryanwweber left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Thank you for contributing these changes! It's great to see work ongoing with the radiation models. Other than my style comment below, and the comment from @ischoegl, I have three other suggestions to help improve this:

  1. Please keep line lengths to less than 88 characters
  2. Please include the citation for where this data is derived from. Is there a specific license associated with the HITRAN data? I couldn't determine from a brief look at their website.
  3. Please add a test for the added code and/or an example that shows the effect of radiation for these species

Thanks again for your contribution! If there's anything I can do to help you with these comments, please let me know!


// Temperatures for Planck optical path length evaluation, K
const int OPL_table_size = 28;
const doublereal TemperatureOPL[OPL_table_size] = {300.0,400.0,500.0,600.0,700.0,800.0,900.0,1000.0,1100.0,1200.0,1300.0,1400.0,1500.0,1600.0,1700.0,1800.0,1900.0,2000.0,2100.0,2200.0,2300.0,2400.0,2500.0,2600.0,2700.0,2800.0,2900.0,3000.0};
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please don't use doublereal in new code, you can just use double. I also think this can be replaced with a std::vector. Finally, please indent your code with spaces, instead of tab characters.

Copy link
Author

@lavrenyukiv lavrenyukiv Jan 26, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

  1. Done
  2. I downloaded HAPI library from https://hitran.org/hapi/ (citation: R.V. Kochanov, I.E. Gordon, L.S. Rothman, P. Wcislo, C. Hill, J.S. Wilzewski, HITRAN Application Programming Interface (HAPI): A comprehensive approach to working with spectroscopic data, J. Quant. Spectrosc. Radiat. Transfer 177, 15-30 (2016), https://doi.org/10.1016/j.jqsrt.2016.03.005).
    Then I downloaded spectra for more important for combustion tasks molecules. It seems that the number of spectra available for download is limited (or I haven't figured it out yet). Then I computed spectrum averaged coefficients for optically thin medium (with Planck function) by myself. I visually compared my results with the results presented earlier(http://www.sandia.gov/TNF/radiation.html), but I am not completely sure.
  3. I prepared example *.py file. In which place it could be placed?

- replace tabs with spaces
- replace doublereal types with double
Copy link
Member

@ischoegl ischoegl left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

One more comment beyond my initial observation.

for (size_t n = 0; n <= 5; n++) {
k_P_H2O += c_H2O[n] * pow(1000 / T(x, j), (double) n);
for (int k = 0; k < OPL_table_size; k++) {
if (T(x, j) < TemperatureOPL[k]) {
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This implementation could be improved in terms of computational efficiency as it appears that the same temperature lookup loop is called for each of the species. Beyond, It should be possible to consolidate the species-specific code sections, which - as far as I can tell - are mostly identical except for the species index.

Copy link
Member

@ischoegl ischoegl Jan 28, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

PS: I am also not convinced that a linear interpolation within a lookup table is an improvement over a polynomial fit - unless the approximation is poor for the newly added curves. (Btw, how do new H20/CO2 curves compare to the previous implementation?)

Copy link
Author

@lavrenyukiv lavrenyukiv Jan 28, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I will make a comparison between old and new curves for H2O and CO2.
Linear interpolation is not an improvement, but rather a more general assumption. For other species we have not polynomial type assumption and polynomial coefficients. But we have spectra for it. Integration of spectra is a very expensive operation in terms of calculations. So tabulated data is accurate and fast.
The use of a polynomial can lead to negative values of the absorption coefficients.
I found that CO, CH4 and OH molecules spectra integrated before and presented by me are wrong, so I will change this data.

Copy link
Author

@lavrenyukiv lavrenyukiv Jan 28, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Plots in logarithmic coordinates are at the left part, the same plots in linear coordinates are at the right part.
Results for H2O look similar. Most likely the polynomial approximation was done on a linear scale. So there is a large difference in the logarithmic scale at high temperatures. So linear interpolation in logarithmic scale looks preferable.
The results for CO2 are very different. I will try to find the data from the original source.
models_check

If we compute S * a * (T^4-300^4), where S is Stefan - Boltzmann constant, we can make an estimate of the absorbed energy(W/m3/atm). Here the differences between approximations will be huge.
models_check_2

Copy link
Member

@ischoegl ischoegl Jan 28, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@lavrenyukiv ... thank you - this looks interesting. I was not involved in the original model, so it's hard to tell. I agree that the differences for energy are substantial. The original polynomial coefficients appear to be taken directly from the TNF link you had posted (where it states that "These curve fits were generated for temperatures between 300K and 2500K and may be very inaccurate outside this range."; the fit itself is for 1000/T so deviations may not look as dramatic?). I am not familiar with RADCAL and assumptions used therein, so some legitimate differences may exist.

PS: There appears to be an updated version of the radcal code on GitHub.

PPS: Incidentally, there was a post on the User Group about unphysically low temperatures but it's hard to tell whether these are related as no code was posted.

Copy link
Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Yes, it is a negative absorption coefficient at low temperatures.
models_check_3

Copy link
Author

@lavrenyukiv lavrenyukiv Jan 30, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I am trying to compare spectra.
Results for CO2 from https://github.com/firemodels/radcal/blob/master/Documentation/Figures/CO2_300K.png

image

and spectra from hapi

full range
co2_full

the RADCAL range
co2_crop

detailed view
co2_closest

Spectra look similar in the specified wavenumber range. The hapi has a greater resolution.

- Added coefficients for CO2, H2O, CO, CH4, OH, N2O, NH3, H2O2, H2CO, NO2, NO, HCN, C2H2, C2H4, C2H6, C2N2, C4H2, HCOOH, HNO3, O3
- Temprature range expanded to 200 - 3500 K
- Definitions for radiation model added
@ischoegl
Copy link
Member

ischoegl commented Feb 9, 2021

@lavrenyukiv ... I just noticed that there may be some related work on radiation properties posted in enhancements (Cantera/enhancements/issues/72). I believe it mainly addresses properties and not the implementation itself.

@lavrenyukiv
Copy link
Author

Cross validation results for HF, HCl, HBr, CO.
Points are from the article DOI: 10.1115/1.1416689
HF_validation

@speth
Copy link
Member

speth commented Feb 19, 2021

I think some improvements and generalizations to the way radiation is handled in Cantera's 1D code would be very welcome.

With the interest in extending the set of species handled by the radiation model, and the fact that data for so many more species are available from HITRAN, I think it would make sense to start reading this data in from a YAML input file, so users could provide their own data for additional species. The data could either be include as a new section within each species definition, or in a separate file. I'd also like to make sure that we keep the door open to allowing other radiation models, like the ones being suggested in Cantera/enhancements#72.

I think the choice of log-linear interpolation seems like a reasonable alternative to using polynomial fits. It's certainly familiar enough from its use in the PLOG rate constant calculations.

@ischoegl
Copy link
Member

@lavrenyukiv … while some time has passed since this was initially discussed, I believe there have been improvements in the YAML interface that should make the import of external data a lot simpler, especially if you interact from Python. Within C++, things are handled using an object named AnyMap that is fairly easy to use, and should work nicely for a more generic approach to radiative heat losses.

@bryanwweber bryanwweber removed their request for review February 10, 2022 20:48
@ischoegl
Copy link
Member

ischoegl commented Jun 8, 2022

@lavrenyukiv ... are you still pursuing this PR or should it be closed?

@ischoegl
Copy link
Member

ischoegl commented Nov 5, 2024

Superseded by #1799

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

Successfully merging this pull request may close these issues.

4 participants