# Bateman Equation in C++

1. Jun 13, 2007

### lapo3399

Hello,

I am attempting to create a program that will evaluate the Bateman equation for radiactive decay series, given several decay constants stored in a text file. I am not sure where this is going wrong, but no matter what chain member I enter in as the desired nuclide to find the mass remaining after 1 million years, it outputs NAN - even when I choose the nuclide that I am certain to have a perceivable value. The equation I am emulating later in the program is shown below, please help to explain why this is not working! All help is appreciated!

http://img515.imageshack.us/img515/5097/screenshotml6.jpg [Broken]

Code (Text):

#include <iostream>
#include <math.h>
#include <fstream>
#include <string>

using namespace std;

int main(int argc, char** argv)
{
string line; //for file input
int iter; //stores which nuclide will be calculated for

double dcs[23]; //array for decay constants
int anum = 1; //counting var

ifstream deconst;
deconst.open ("decayconstants.txt", ios::in); //file with decay constants
if (deconst.is_open())
{
while (! deconst.eof())
{
getline(deconst, line);
dcs[anum]=strtod(line.c_str(), NULL); //store as double in decay constant array
anum++;
}

}

std::cout << "What nuclide in the chain should be evaluated?";
std::cin >> iter;

//I separated the equation into multiple parts

double consts=1; //part one: the multiplied decay constants at the beginning of the equation
double init = 100; //the intial mass, N(O)
double time = 1000000; //the variable t, one million years
double sum=0; // the second half of the equation, the summation
double totalmass; //the product of consts,init, and sum

for (int i = 1; i < iter; i++)
{
consts=consts*dcs[i]; //multiply all the necessary decay constants
}

for (int q = 1; q < iter; q++) //the summation
{
double expo = exp(-dcs[q]*time); //determine the numerator in the summation
double denom = 1;

for (int w = 1; w < iter; w++) //the product
{
if (w==iter) //k cant be j
{
break;
}
else
{
denom = denom * (dcs[w]-dcs[q]); //multiply together to get the product
}

}

sum = sum + (expo/denom); //add up all of the summations

}

totalmass = sum*consts*init;

std::cout << totalmass;

return 0;
}

Thanks again!

Last edited by a moderator: May 2, 2017
2. Jun 13, 2007

### chroot

Staff Emeritus
Well, I can't run your program without the file it needs to read, but I'd suggest putting some cout statements throughout the code to print out the intermediate values.

- Warren

3. Jun 13, 2007

### lapo3399

I've attempted that, and I finally realized what the problems where this morning. Firstly, the first for loop should read "i<(iter-1)", and I shouldn't have included the break statement in the other for loop.

4. Nov 1, 2007

### mehafroozs

Hi,
I see that you said your program worked. I am working with a few students to develope a detector but i also needed them to study the decan chains of Radon. I think your program will be useful, do you think i could get a copy of it?

5. Nov 1, 2007

### KingMike

have you checked all of the values prior to use?
If so, are you trying to calculate a null with an int/double?

6. Dec 3, 2010

### ustoehlker

Hi,
would it be possible to get a copy of the source code?