Simcenter STAR-CCM+ Understanding User Coding for Modified Internal Rates of Individual Reactions - with Example

Reaction rates play a crucial role in understanding and predicting chemical reactions. They determine how fast a reaction occurs and are essential in various fields, such as chemistry, physics, and engineering. Simcenter STAR-CCM+ allows the import of all CHEMKIN defined reaction rates by default. In some cases this is still too limited. For these situations  the user-coded Modified Internal Rates of Individual Reactions was introduced in Simcenter STAR-CCM+ 2310 allowing researchers and engineers to customize and modify reaction rates based on their requirements. With user code you can program freely in c, c++ or Fortran - using external libraries if necessary. 

Please find a simple example simulation and code in c and c++ attached.

For the latest information on this method please conduct the latest version of the User Guide: Home > Physics Simulation > Reacting Flow > Reacting Flow General Reference > Reacting Model Reference > User-Defined Reaction Rate Source Terms

Accessible Data Inputs

User Coding for Reaction Rates provides a flexible approach to adjust reaction rates by accessing relevant data through a structured interface and even utilizing user-defined field functions. This method offers these main inputs:

1. zi: A pointer to a vector containing gas phase specific mole fractions, surface site fractions, bulk activities, and cell temperature.
2. struct UserAccessibleData: A structured data object that holds important information, such as the number of reactions, loaded field functions, and absolute pressure. The contents of the struct can change, always use the latest information given in the user guide.
3. The wall temperature is directly available as T, whereas the cell temperature is accessible as the last element in zi using zi[nSpe].
4. You can also load your own field functions and retrieve them from the struct. To include field functions in your user code select them below the user defined reaction source in the reacting system properties settings. 
   This is also where you enable the method in general.
5. reacRateForward and reacRateBackward are in and outputs at the same time. These pointers can be used to retrieve, adjust and set the reaction rates.

Handling Indices

One challenge with User Coding for Reaction Rates is the handling of indices. The indices of species, reactions, and field functions can change as new elements are added or removed. To address this issue, a workaround involves defining an enum to assign names to species, reactions, and field functions. This approach enables the use of names instead of indices, making the code more readable and maintainable.

For example: here, we adjust the (forward) reaction rate of "CH4+2PT(S)=>CH3(S)+H(S)" by using the partial pressure of CH4, with the forward reaction rate reacRateForward and the partial pressure of Methane pi[CH4].

reacRateForward[_CH42PTStoCH3SHS] *= 1 / (1 + pi[CH4]); 

Next example: here, we modify the reaction rate of all reactions based on a user field function f[_onOff]. For testing purposes, the user function dynamically deactivates certain cells by multiplying by either 0 or 1. The User Field Function $_onOff and the resulting Chemistry Heat Release Rate is shown in the first two animations of the example surface chemistry simulation.

for (unsigned int i = 0; i < nReactions; ++i) { 
	reacRateForward[i] *= f[_onOff]; 
	reacRateBackward[i] *= f[_onOff]; 
} 

The enum indices can be either read manually from the interface or created using user code. Simply define the VERBOSE output option in the code and run once. This creates the indices automatically in the output, copy paste them into the enum definition on top. 

typedef enum {  
	CH4, CO, CO2, H2, H2O, N2, O2, OH,				// fluid phase species, 
	CS, CHS, CH2SS, CH3S, COS, CO2S, HS, H2OS, OS, OHS, PtS,	// surface site species, 
	PtB								// bulk species IN CORRECT ORDER
} SpeciesOrder;  

Debugging: Exporting Field Functions from Reaction Rate User Coding

Another challenge with the reaction rate user coding is the limited options for debugging output. Be cautious with printf. The user code is executed for each CVODE iteration (plus numerical linearization executions) and for each cell. Therefore, each line of the user code is executed VERY often and the output using printf is massive.
A workaround is build into the provided example. You can now export field functions from the reaction rate User Coding.

This is done via global vectors storing the data during the reaction rate calculations. The stored function is recalled by a normal User-Coded Field Function and can be displayed as such. Have a look at the bottom two animations above. Here the previously used expression 1/(1+pi[CH4]) is plotted once using a User Field function '_1/(1+pi)' and the User-Coded Field Function 'User transferData' for comparison with the same result.
Be aware this method works best in serial simulations. For parallel simulations the used 'Cell Index' function is not uniquely defined.

To enable this feature please define the ENABLE_DATA_TRANSFER in the provided user code. In the line 'transferData[index] =' you can define what data you want to be transferred and shown in 'User transferData' for debugging.