ver-1gc

Series Chemical Reactions

Problem set up

This verification problem is taken from Ambrosek and Longhurst (2008) and builds on the capabilities verified in ver-1g for simple chemical reactions. This case is simulated in (test/tests/ver-1gc/ver-1gc.i).

This problem models a set of chemical reactions in series with three species: $var element = document.getElementById("moose-equation-3d0a7ac3-b88e-47c7-b99c-5310a6101b85");katex.render("A", element, {displayMode:false,throwOnError:false});$ , $var element = document.getElementById("moose-equation-642439b8-90d0-4138-9430-9b3608e200b6");katex.render("B", element, {displayMode:false,throwOnError:false});$ , and $var element = document.getElementById("moose-equation-9e75befa-85c4-45f7-9499-ed44a6adab9d");katex.render("C", element, {displayMode:false,throwOnError:false});$ . The system was configured so that the enclosure initially contained only species $var element = document.getElementById("moose-equation-3b90ee45-d45e-48e1-afe6-0814a954f5ce");katex.render("A", element, {displayMode:false,throwOnError:false});$ . At time $var element = document.getElementById("moose-equation-a27e3bfa-87b7-4069-85f0-8bb73d1fbcf0");katex.render("t \\geq 0", element, {displayMode:false,throwOnError:false});$ s, the reactions were allowed to proceed. The reactions that were modeled are

$(1)var element = document.getElementById("moose-equation-7b70f6ab-90ab-40a3-9eb1-9a4f637a1f3e");katex.render(" A \\xrightarrow{\\mathit{k_1}} B \\xrightarrow{\\mathit{k_2}} C,", element, {displayMode:true,throwOnError:false});$ with $var element = document.getElementById("moose-equation-66d0f3ed-e5e2-4169-b956-46474a99ac24");katex.render("k_1", element, {displayMode:false,throwOnError:false});$ and $var element = document.getElementById("moose-equation-f0091786-29de-4661-ade1-8fc20ebde7b1");katex.render("k_2", element, {displayMode:false,throwOnError:false});$ the reaction rates.

The concentration of each species is therefore described as $(2)var element = document.getElementById("moose-equation-5cb57bef-3824-47c3-8e4d-b3b747672739");katex.render(" \\frac{dc_A}{dt} = - k_1 c_A,", element, {displayMode:true,throwOnError:false});$ $(3)var element = document.getElementById("moose-equation-421e0ea0-97da-4589-a107-25ef543e32e6");katex.render(" \\frac{dc_B}{dt} = k_1 c_A - k_2 c_B,", element, {displayMode:true,throwOnError:false});$ $(4)var element = document.getElementById("moose-equation-545d24da-4a24-4085-bba5-6b252918ed66");katex.render(" \\frac{dc_C}{dt} = k_2 c_B,", element, {displayMode:true,throwOnError:false});$ with $var element = document.getElementById("moose-equation-30b58343-c017-4a78-b0a2-e097a5a9fcc0");katex.render("c_i", element, {displayMode:false,throwOnError:false});$ the concentration of species $var element = document.getElementById("moose-equation-d785a6dd-98c6-4712-b137-2de79d940650");katex.render("i", element, {displayMode:false,throwOnError:false});$ .

Analytical solution

Fogler (1999) provides the analytical equations for the time evolution of the concentrations of $var element = document.getElementById("moose-equation-8a73cbbd-cb76-4717-969c-aa733768f94b");katex.render("A", element, {displayMode:false,throwOnError:false});$ and $var element = document.getElementById("moose-equation-7bc46cd8-2231-4cb5-8321-c69a409937f1");katex.render("B", element, {displayMode:false,throwOnError:false});$ as $(5)var element = document.getElementById("moose-equation-6dbfc0dc-4834-4503-bc97-68c5b515a429");katex.render(" c_A(t) = c_{A0} \\exp\\left( -k_1 t\\right),", element, {displayMode:true,throwOnError:false});$ and $(6)var element = document.getElementById("moose-equation-3a0e2142-3c32-47be-8800-2d3169eb3b55");katex.render(" c_B(t) = k_1 c_{A0} \\frac{\\exp\\left( -k_1 t\\right) - \\exp\\left( -k_2 t\\right)}{k_2-k_1},", element, {displayMode:true,throwOnError:false});$ where $var element = document.getElementById("moose-equation-ee9bfbd8-6078-4787-877e-7f6b89118445");katex.render("t", element, {displayMode:false,throwOnError:false});$ is the time in s, $var element = document.getElementById("moose-equation-1f5fc412-f3f2-4ab5-a342-92c5d539ba95");katex.render("c_{A0} = 2.415 \\times 10^{14}", element, {displayMode:false,throwOnError:false});$ atoms/m $var element = document.getElementById("moose-equation-0aabfcae-3432-4f67-8492-5396aa8f364a");katex.render("^3", element, {displayMode:false,throwOnError:false});$ is the initial concentration of species $var element = document.getElementById("moose-equation-a164fcbf-caec-44b4-80cb-5496be828396");katex.render("A", element, {displayMode:false,throwOnError:false});$ , $var element = document.getElementById("moose-equation-da9f21d9-b555-4335-9776-3e878894821b");katex.render("k_1 = 0.0125", element, {displayMode:false,throwOnError:false});$ s $var element = document.getElementById("moose-equation-584151ff-75cc-4b2f-a89d-c11a133a491a");katex.render("^{-1}", element, {displayMode:false,throwOnError:false});$ , and $var element = document.getElementById("moose-equation-93a46c61-a345-48c3-b6ee-0a343ee1e2fe");katex.render("k_2 = 0.0025", element, {displayMode:false,throwOnError:false});$ s $var element = document.getElementById("moose-equation-2f818a6c-b72d-4353-a3e3-ccb0699fe035");katex.render("^{-1}", element, {displayMode:false,throwOnError:false});$ .

The concentration of $var element = document.getElementById("moose-equation-aedecc0a-fe04-4ca2-834b-4c365adb3ea3");katex.render("C", element, {displayMode:false,throwOnError:false});$ was found by applying a mass balance over the system in Ambrosek and Longhurst (2008). From the stoichiometry of this reaction it was found that $(7)var element = document.getElementById("moose-equation-123ffdf1-bc79-4b6e-857e-f6a6ed8230d6");katex.render(" c_C(t) = c_{A0} - c_A(t) - c_B(t).", element, {displayMode:true,throwOnError:false});$ The time evolution of the species concentrations from the analytical solution is provided in Figure 1.

Results and comparison against analytical solution

The comparison of TMAP8 results against the analytical solution is shown in Figure 1. The match between TMAP8's predictions and the analytical solution is satisfactory, with root mean square percentage errors (RMSPE) of RMSPE = 4.78 % for species $var element = document.getElementById("moose-equation-ac073ac9-4e5c-4662-b25d-b07afe212d54");katex.render("A", element, {displayMode:false,throwOnError:false});$ , RMSPE = 3.20 % for species $var element = document.getElementById("moose-equation-188dda1b-3258-412f-ab04-77945c6561a1");katex.render("B", element, {displayMode:false,throwOnError:false});$ , and RMSPE = 0.38 % for species $var element = document.getElementById("moose-equation-51dfe708-4fc5-4aef-9191-1ca2fb96634a");katex.render("C", element, {displayMode:false,throwOnError:false});$ . The larger values of RMSPE for species $var element = document.getElementById("moose-equation-f64b2e07-e5bb-47a3-ad5b-84d0ac350b50");katex.render("A", element, {displayMode:false,throwOnError:false});$ and $var element = document.getElementById("moose-equation-80f1f75d-a1e2-4d13-9f36-29e04bf0ff7e");katex.render("B", element, {displayMode:false,throwOnError:false});$ are due to the small average values of these concentrations over time, and do not reflect on the accuracy of the TMAP8 solve itself.

Figure 1: Comparison of partial pressures of species in a series reaction predicted by TMAP8 and provided by the analytical solution. The RMSPE is the root mean square percent error between the analytical solution and TMAP8 predictions.

Input files

The input file for this case can be found at (test/tests/ver-1gc/ver-1gc.i), which is also used as test in TMAP8 at (test/tests/ver-1gc/tests).

References

James Ambrosek and GR Longhurst. Verification and Validation of TMAP7. Technical Report INEEL/EXT-04-01657, Idaho National Engineering and Environmental Laboratory, December 2008.

@techreport{ambrosek2008verification,
    author = "Ambrosek, James and Longhurst, GR",
    title = "{Verification and Validation of TMAP7}",
    year = "2008",
    number = "INEEL/EXT-04-01657",
    month = "December",
    institution = "Idaho National Engineering and Environmental Laboratory"
}

H. S. Fogler. Elements of Chemical Reaction Engineering. Prentice Hall PTR, 3rd edition, 1999.

@book{Fogler1999,
    author = "Fogler, H. S.",
    city = "Upper Saddle River, NJ",
    edition = "3rd",
    publisher = "Prentice Hall PTR",
    title = "Elements of Chemical Reaction Engineering",
    year = "1999"
}

(test/tests/ver-1gc/ver-1gc.i)

concentration_A_0 = ${units 2.415e14 at/m^3 -> at/mum^3} # atoms/microns^3 initial concentration of species A
k_1 = 0.0125 # 1/s reaction rate for first reaction
k_2 = 0.0025 # 1/s reaction rate for second reaction
end_time = 1500 # s

[Mesh]
  type = GeneratedMesh
  dim = 1
[]

[Variables]
  [c_A]
    initial_condition = ${concentration_A_0}
    scaling = 1e4
  []
  [c_B]
    scaling = 1e4
  []
  [c_C]
    scaling = 1e4
  []
[]

[Kernels]
  # Equation for species A
  [timeDerivative_c_A]
    type = ADTimeDerivative
    variable = c_A
  []
  [MatReaction_A]
    type = ADMatReactionFlexible
    variable = c_A
    vs = 'c_A'
    coeff = -1
    reaction_rate_name = K1
  []
  # Equation for species B
  [timeDerivative_c_B]
    type = TimeDerivative
    variable = c_B
  []
  [MatReaction_B_1]
    type = ADMatReactionFlexible
    variable = c_B
    vs = 'c_A'
    coeff = 1
    reaction_rate_name = K1
  []
  [MatReaction_B_2]
    type = ADMatReactionFlexible
    variable = c_B
    vs = 'c_B'
    coeff = -1
    reaction_rate_name = K2
  []
  # Equation for species C
  [timeDerivative_c_C]
    type = TimeDerivative
    variable = c_C
  []
  [MatReaction_C]
    type = ADMatReactionFlexible
    variable = c_C
    vs = 'c_B'
    coeff = 1
    reaction_rate_name = K2
  []
[]

[Materials]
  [K1]
    type = ADParsedMaterial
    property_name = 'K1'
    expression = ${k_1}
  []
  [K2]
    type = ADParsedMaterial
    property_name = 'K2'
    expression = ${k_2}
  []
[]

[BCs]
  [c_A_neumann] # No flux on the sides
    type = NeumannBC
    variable = c_A
    boundary = 'left right'
    value = 0
  []
  [c_B_neumann] # No flux on the sides
    type = NeumannBC
    variable = c_B
    boundary = 'left right'
    value = 0
  []
  [c_C_neumann] # No flux on the sides
    type = NeumannBC
    variable = c_C
    boundary = 'left right'
    value = 0
  []
[]

[Postprocessors]
  [concentration_A]
    type = ElementAverageValue
    variable = c_A
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [concentration_B]
    type = ElementAverageValue
    variable = c_B
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [concentration_C]
    type = ElementAverageValue
    variable = c_C
    execute_on = 'INITIAL TIMESTEP_END'
  []
[]

[Executioner]
  type = Transient
  scheme = bdf2
  nl_rel_tol = 1e-11
  nl_abs_tol = 1e-50
  l_tol = 1e-10
  solve_type = 'NEWTON'
  petsc_options = '-snes_ksp_ew'
  petsc_options_iname = '-pc_type'
  petsc_options_value = 'lu'
  end_time = ${end_time}
  dtmax = 50
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 1.0
    optimal_iterations = 15
    growth_factor = 1.25
    cutback_factor = 0.8
  []
[]

[Outputs]
  exodus = true
  csv = true
[]

(test/tests/ver-1gc/ver-1gc.i)

concentration_A_0 = ${units 2.415e14 at/m^3 -> at/mum^3} # atoms/microns^3 initial concentration of species A
k_1 = 0.0125 # 1/s reaction rate for first reaction
k_2 = 0.0025 # 1/s reaction rate for second reaction
end_time = 1500 # s

[Mesh]
  type = GeneratedMesh
  dim = 1
[]

[Variables]
  [c_A]
    initial_condition = ${concentration_A_0}
    scaling = 1e4
  []
  [c_B]
    scaling = 1e4
  []
  [c_C]
    scaling = 1e4
  []
[]

[Kernels]
  # Equation for species A
  [timeDerivative_c_A]
    type = ADTimeDerivative
    variable = c_A
  []
  [MatReaction_A]
    type = ADMatReactionFlexible
    variable = c_A
    vs = 'c_A'
    coeff = -1
    reaction_rate_name = K1
  []
  # Equation for species B
  [timeDerivative_c_B]
    type = TimeDerivative
    variable = c_B
  []
  [MatReaction_B_1]
    type = ADMatReactionFlexible
    variable = c_B
    vs = 'c_A'
    coeff = 1
    reaction_rate_name = K1
  []
  [MatReaction_B_2]
    type = ADMatReactionFlexible
    variable = c_B
    vs = 'c_B'
    coeff = -1
    reaction_rate_name = K2
  []
  # Equation for species C
  [timeDerivative_c_C]
    type = TimeDerivative
    variable = c_C
  []
  [MatReaction_C]
    type = ADMatReactionFlexible
    variable = c_C
    vs = 'c_B'
    coeff = 1
    reaction_rate_name = K2
  []
[]

[Materials]
  [K1]
    type = ADParsedMaterial
    property_name = 'K1'
    expression = ${k_1}
  []
  [K2]
    type = ADParsedMaterial
    property_name = 'K2'
    expression = ${k_2}
  []
[]

[BCs]
  [c_A_neumann] # No flux on the sides
    type = NeumannBC
    variable = c_A
    boundary = 'left right'
    value = 0
  []
  [c_B_neumann] # No flux on the sides
    type = NeumannBC
    variable = c_B
    boundary = 'left right'
    value = 0
  []
  [c_C_neumann] # No flux on the sides
    type = NeumannBC
    variable = c_C
    boundary = 'left right'
    value = 0
  []
[]

[Postprocessors]
  [concentration_A]
    type = ElementAverageValue
    variable = c_A
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [concentration_B]
    type = ElementAverageValue
    variable = c_B
    execute_on = 'INITIAL TIMESTEP_END'
  []
  [concentration_C]
    type = ElementAverageValue
    variable = c_C
    execute_on = 'INITIAL TIMESTEP_END'
  []
[]

[Executioner]
  type = Transient
  scheme = bdf2
  nl_rel_tol = 1e-11
  nl_abs_tol = 1e-50
  l_tol = 1e-10
  solve_type = 'NEWTON'
  petsc_options = '-snes_ksp_ew'
  petsc_options_iname = '-pc_type'
  petsc_options_value = 'lu'
  end_time = ${end_time}
  dtmax = 50
  [TimeStepper]
    type = IterationAdaptiveDT
    dt = 1.0
    optimal_iterations = 15
    growth_factor = 1.25
    cutback_factor = 0.8
  []
[]

[Outputs]
  exodus = true
  csv = true
[]

(test/tests/ver-1gc/tests)

[Tests]
  design = 'ADMatReactionFlexible.md'
  issues = '#12 #104'
  verification = 'ver-1gc.md'
  [ver-1gc]
    type = Exodiff
    input = 'ver-1gc.i'
    exodiff = ver-1gc_out.e
    requirement = 'The system shall be able to model a series of chemical reactions involving three species and calculate the concentrations of each species as a function of time.'
  []
  [ver-1gc_csv]
    type = CSVDiff
    input = 'ver-1gc.i'
    should_execute = False # this test relies on the output files from ver-1gc, so it shouldn't be run twice
    csvdiff = ver-1gc_out.csv
    requirement = 'The system shall be able to model a series of chemical reactions involving three species and calculate the concentrations of each species as a function of time and to generate CSV data for use in comparisons with the analytic solution over time.'
    prereq = ver-1gc
  []
  [ver-1gc_comparison]
    type = RunCommand
    command = 'python3 comparison_ver-1gc.py'
    requirement = 'The system shall be able to generate comparison plots between the analytical solution and simulated solution of a series of chemical reactions involving three species and calculate the concentrations of each species as a function of time and to generate CSV data for use in comparisons with the analytic solution over time.'
    required_python_packages = 'matplotlib numpy pandas scipy os'
  []
[]

Problem set up
Analytical solution
Results and comparison against analytical solution
Input files
References