SSUMES - Quick Start step-5

[Top] < Prev Step | Next Step >

SSUMES - Quick Start step-5

Example-5: Chebyshev-polynomial fit

The last example is the Chebyshev-polynomial fitting of the rate coefficients calculated by RRKM theory. The example rate coefficients used here were calculated by carate program with the input rc2h2ph_m_ca_cheb.inp. This is essentially the same input as rc2h2ph_m_ca.inp used in the Example problem-3 but with different temperature and pressure specifications suitable for Chebyshev fitting as shown below.

# phenyl-c2h2 quadruple-well model - chemical activation (Gaus-Chebyshev grid)

tempGauChebGrd 400 2000 10
pressUnit atm
pressGauChebGrd 0.1 1000 10

Input

A part of the content of the input file, rc2h2ph_m_ca_fit_w1c3.inp, is shown below.

# Chebyshev-polynomial fit for rc2h2ph_m_ca_cheb well-1 channel-3

tempCheb 400 2000
pressCheb 0.1 1000

numPars 0 0
! numPars 6 3

TPk_Table{
 ! T, P, k
401.98   0.105833553   8.21632E-16
418.234  0.105833553   1.28085E-15
453.082  0.105833553   3.1334E-15
   :
1951.94  944.8802632   6.68191E-12
}

The keys 'tempCheb' and 'pressCheb' are used to specify the ranges of temperature and pressure, respectively. Unlike the other solver programs in SSUMES, the unit of the pressure is 'atm' which is used by the Chemkin program. The 'numPars' is used to specify the numbers of parameters in the temperature and pressure axes. We will first run the fitcheb program with 'numPars 0 0' which tells the program to scan the number of parameters. The 'TPk_Table' block contains the list of rate coefficients, containing one set of T [K], P [atm], and k in a line.

Parameter scan

Run fitcheb as:
fitcheb rc2h2ph_m_ca_fit_w1c3
With 'numPars 0 0', the program fitcheb scans the number of parameters and a diagnostic is written in rc2h2ph_m_ca_fit_w1c3_fitcheb_diag.csv. The contents of the diagnostic output looks like:

The first table lists the standard deviation (STD) for each pair of nParT (number of parameters in temperature axis) and nParP (number of parameters in pressure axis). For example, the cell (nParT, nParP) = (6, 3) shows the standard deviation (0.041048) of the common logarithm of the input rate coefficient, log₁₀k, from the Chebyshev fitting with six parameters in T-axis and three parameters in P-axis. Note that this standard deviation corresponding to the factor of 10^0.041048 = 1.099 deviation, that is, about 10% error as the standard deviation.
The second and third tables contains the mean absolute deviation (MAD) and maximum absolute deviation (MAXAD), respectively, of log₁₀k. This output is helpful to determine the most cost-effective selection of (nParT, nParP). Note that the increasing the number of parameters results in the high-cost (or slow computation) in the evaluation of the rate coefficients in programs such as Chemkin.

Chebyshev coefficients

Edit the rc2h2ph_m_ca_fit_w1c3.inp as:
```
! numPars 0 0
numPars 6 3
```
if you select (nParT, nParP) = (6, 3).
Then, rerun fitcheb as:
fitcheb rc2h2ph_m_ca_fit_w1c3
This time, the program fitcheb generates a file, rc2h2ph_m_ca_fit_w1c3_fitcheb_out.csv, containing the 6×3 coefficients for a Chebyshev-polynomial as below:

Input for kcheb

The program kcheb can be used to evaluate the rate coefficients, k(T, P) from a set of Chebyshev coefficients. The sample input, rc2h2ph_m_ca_rate_w1c3.inp, is shown below:

# Chebyshev-polynomial rate for rc2h2ph_m_ca_cheb well-1 channel-3

tempCheb 400 2000
pressCheb 0.1 1000

tempRecipRange 100000 55 245 5
tempList 1999.9999 400.00001
pressUnit atm
pressLog10Range 0 2 1

coefTable{
 6 3
 -13.983,-1.09038,-0.103847
 2.9531,0.96261,0.0334369
 0.0865576,0.300817,0.063038
 -0.134253,-0.0188854,0.0366201
 -0.0763344,-0.097919,0.00221055
 -0.0350195,-0.073087,-0.0154463
}

The arguments to 'tempCheb' and 'pressCheb' should be the same as in those used in fitcheb input (and note that the unit of pressure is 'atm'). Temperature and pressure range inputs are similar to those in solvers (carate etc.). The 'coefTable' block contains the number of parameters and coefficients.

Rate coefficients calculation from Chebyshev coefficients

Run kcheb as:
kcheb rc2h2ph_m_ca_rate_w1c3
The output is written into rc2h2ph_m_ca_rate_w1c3_kcheb_out.csv which looks like:

[Top] < Prev Step | Next Step >