Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Simulating a radioactive pile

  1. Dec 6, 2008 #1
    I'm looking to simulate a radioactive pile in the context of a simple computer game. The problem is, however, that for various reasons I can't use a straightforward tick-by-tick simulation.

    Instead I'm looking for an equation that can be solved based on time and would give me the amount of stuff remaining in the pile. And this is quite easy:

    [tex]R\left(t\right)=R_0\left(\frac{1}{2}\right)^{\frac{t}{t_h}}[/tex]

    Where [tex]t_h[/tex] is the half-life of the element and [tex]R_0[/tex] is the initial amount at [tex]t=0[/tex]

    But there's a complication. Extra material can be added to the pile at any [tex]t[/tex] (For example, 5 units added at [tex]t=15[/tex]) Which is still quite simple. Each added amount could be modeled as a separate pile and the [tex]R\left(t\right)[/tex] of all piles could be added together:

    [tex]R_{total}\left(t\right)=\sum_{i=0}^{n}R_{i}\left(t\right)[/tex]

    (Without taking piles into account if [tex]t < t_{i_{0}}[/tex])

    But there's one more complication and this is the one I'm unable to crack.

    Extra material is to be added at a steady rate over several time units. (For example, 5 units added between [tex]t_0=17[/tex] and [tex]t_1=19[/tex]) So there's a function which grows a pile between the pile's [tex]t_0[/tex] and [tex]t_1[/tex] and looks something like:

    [tex]R(t)=\left(\frac{m}{a} t - t_0\right)[/tex] if [tex]t\in\left[t_0,t_1\right][/tex] and m = the total amount added to create the pile and [tex]a=t_1 - t_0[/tex] or the time it takes to create the pile.

    The problem is that I can't figure out how to start applying the decay function while the pile is still being created. At [tex]t_1[/tex] the amount of stuff in the pile will be equal to m, which shouldn't happen as some of the stuff should have decayed already.

    So my question would be: how do I model the lifetime of a decaying pile in such a way that takes into account the fact that it can take several time units to initially create the pile with material being added at a steady rate?

    What I'd like to end up with is a graph showing the past, present, and future combined size of all piles and the ability to add events at any point in time that contribute new material to this combined pile.
     
  2. jcsd
  3. Dec 6, 2008 #2

    CompuChip

    User Avatar
    Science Advisor
    Homework Helper

    So let T be the half-value time, and C(t) the function that indicates how much you add at time t. Then you could do something like:

    [tex]R(t + \Delta t) = R(t) \left( \frac12 \right)^{\Delta t / T} + N(t, t + \Delta t)[/tex]
    or
    [tex]R(t + \Delta t) = (R(t) + N(t, t + \Delta t)) \left( \frac12 \right)^{\Delta t / T}[/tex]
    (depending on whether you first add and then decay, or vice versa), and with N denoting the amount added during the step:
    [tex]N(t, t + \Delta t) = \int_{t}^{t + \Delta t} C(t') dt'[/tex]
    Then you set R(0) and solve iteratively (calculating each next step from the previous), or plug in some explicit forms of C(t) and derive a differential equation :smile:

    Does that help any? If not, I apologize, I'm not quite sure what you want/
     
  4. Dec 6, 2008 #3
    If I understand you correctly (and it's been over two years since I last did maths on this level, so I'm sorry if I get something wrong) your solution requires a step-by-step (an iterative) calculation. And that's a somewhat valid solution, but I would much prefer to have a solution that doesn't require a step-by-step calculation.

    I did a quick step-by-step simulation using a small Python script and plotted the results with OpenOffice:

    [​IMG]

    What I'd like to do is end up with a formula that takes [tex]t_0[/tex] (in this case 0) and [tex]t_1[/tex] (in this case 500) and the half-life of the particles and the rate at which they are added and finally, when plotted, spit out the blue line above.
     
  5. Dec 7, 2008 #4

    CompuChip

    User Avatar
    Science Advisor
    Homework Helper

    Yes, you understood me correctly.

    As I said, it is much easier to write down a recursion. Just to show you what I have in mind, let me give you a relatively easy example. I hope you have your calculus ready :smile:
    Suppose that every time-unit (second) we add a constant amount of n particles, and the half-life is T. Then the recursion equation reads
    [tex]R(t + \Delta t) = \left( \frac12 \right)^{\Delta t / T} R(t) + n \Delta.[/tex]
    We can Taylor-expand the power around [itex]\Delta t = 0[/itex], giving
    [tex]R(t + \Delta t) \simeq \left(1 - \frac{\log 2}{T} \Delta \right) R(t) + n \Delta + \mathcal O(\Delta^2)[/tex].
    Therefore,
    [tex]\frac{R(t + \Delta t) - R(t)}{\Delta t} \simeq - \frac{\log 2}{T} R(t) + n + \mathcal O(\Delta)[/tex].
    If we now take the limit [itex]\Delta t \to 0[/itex] we get the differential equation
    [tex]R'(t) = -\frac{\log 2}{T} R(t) + n[/tex].
    If you solve this with an initial condition R(0) = R0, for example, you get
    [tex]R(t) = \frac{n T}{\log 2} \left( 1 - 2^{-t/T} \right) + R_0 2^{-t/T}[/tex].
    Of course, if you want to have an arbitrary time-dependent function n(t) to add the particles, it becomes little more complicated, especially if n(t) will be discontinuous (e.g. a constant > 0 up to some time and 0 from then on). Or you could take multiple of this solutions and patch them up (for example, take n = 100 up to t = 50 with R0 = 0; then from t = 50 take another such solution with n = 0 where you pick R0 to match the end value of the first part, etc.)

    You decide :smile:
    Didn't mean to scare you off, but this is the simplest I could think of... so I suggest the recursion equation.
     
  6. Dec 7, 2008 #5
    Wow. Okay, that is quite scary. I've done the things you describe in the past, so I can see where you're going with this, but it's been so long, I don't remember enough to actually apply any of that anymore.

    I'm gonna try finding a way to simulate this instead of trying to solve it as an equation.
     
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook

Have something to add?