Julia vs Python vs Analytical solutions

[Khaleeh]

I ran each version using the following variables. (I'm not 100% confident that consts is correct)

Variables = [A, D, R, AD, RD, ADR]
concs = [1.0, 0.01, 0.0001, 0.0, 0.0, 0.0]
keqs = [1.0e6, 100, 1, 100] #A+D, D+R,A+DR, AD+R
consts = [[1.0e6, 100, 1, 100] [100, 1, 100, 1.0e6]]
stoich = [[-1, -1, 0, 1, 0, 0],
[ 0, -1, -1, 0, 1, 0],
[ 0, 0, -1, -1, 0, 1],
[-1, 0, 0, 0, -1, 1]]

here are the results.

julia

nekmc

[0.99000, 2.3876e-16, 0.0001, 0.009999, 0.0, 0.0]

kmc

[0.99000, 2.3876e-16, 0.0001, 0.009999, 0.0, 3.3881e-21]

python

nekmc

[9.90e-01 2.00e-06 1.00e-04 9.998e-03 0.0e+00 0.0e+00]

exact

[9.900e-01 1.99e-06 1.00e-04 9.998e-03 -7.209e-14 7.209e-14]

analytical

1.0010701238772311e-14 non-coop [ADR]

These results are interesting but I think they might be edge cases. @PaulWAyers @msricher what do you think?

[Khaleeh]

note: I reduced the number of decimal places to make it easier to read.

[PaulWAyers]

I would guess that it isn't converged, since it seems to me that $R$ should be depleted, but it doesn't seem to be depleted in this example. I expect the final approximations to have (almost) no $D$ and no $R$ left as they are the limiting reagents for reaction 1 and 3.

You may need to run longer then. I don't understand why the Python results are different, but it seems like maybe the Python simulations didn't run for as long. Does the final time print out too?

I don't know how to interpret the analytic result; I expected a list of final concentrations.