search

ukaea / decay-rate-ode-solver

Solver for a hyper-stiff set of ODEs
GNU Lesser General Public License v3.0
1 stars 1 forks source link

Are the expected answers really correct? #4

Open jwscook opened 5 years ago

jwscook commented 5 years ago

@thomasms I found this odd, but wanted to double check with you to make sure.

In decay_inventory_full_1.000000e+15.txt the result for Fe56

301,260560,Fe56,1.000000000000000e+20,3.500000000000000e+21

Is this right? All the nuclides decay down to the most stable iron, which accumulates to exactly 3.5e21? This number seems to neat and tidy to me. Are we sure it's right?

thomasms commented 5 years ago

OK - so this was the one file I did not correct in the last pull request. I will add it in another pull request.

The result for Fe56 changes to 6.999999999999996e+20 (7e20) in the correct run (LSODES), an almost exact amount again. This makes some sense because the data set used to create this is EAF2010, and in this dataset the grandparents of Fe56 are only V56, V57, and Cu56.

The decay paths are then: V56 -> Cr56 -> Mn56 -> Fe56 V57 -> Cr56 -> Mn56 -> Fe56 Cu56 -> Ni56 -> Co56 -> Fe56

All of which have half lives much, much less than the time step (1e15 secs). Therefore, all will decay to Fe56 completely (apart from small numerical issues), and we should expect all contributions from V56 + V57 + Cr56 + Mn56 + Ni56 + Co56 (7 nuclides) to go to Fe56. So, if we include the original 1e20 for Fe56 I would expect (7+1)*1e20 = 8e20 atoms for Fe56 at the end, not 7e20. Have I miscounted? This difference might be due to the nuclear data, since some channels do not exist in some data sets, and it could be that the V57 -> Cr56 channel is omitted in the data, and just decays to the sink. This would bring it to the correct 7e20 atoms, whilst correctly decaying V57.

For your reference, the decay paths can be seen here: http://www.periodictable.com/Isotopes/026.56/index2.full.dm.wt.html

With other nuclear decay libraries, which contain a lot more nuclides, other pathways are enabled, which match with the channels on that link. However for the matrix here we used EAF 2010 data which only contains the nuclides and channels mentioned above.

jwscook commented 5 years ago

The result for Fe56 changes to 6.999999999999996e+20 (7e20) in the correct run (LSODES), an almost exact amount again. This makes some sense because the data set used to create this is EAF2010, and in this dataset the grandparents of Fe56 are only V56, V57, and Cu56.

Ok, I'll take your word for it.

V56 + V57 + Cr56 + Mn56 + Ni56 + Co56 (6 nuclides) to go to Fe56. So, if we include the original 1e20 for Fe56 I would expect (6+1)1e20 = 7*e20 atoms for Fe56 at the end, not 7e20.

There was an out by one error, corrected a above