resolved previous error related to rate object and implemented full o…

…utput of reaction rates (ORNL-Fusion#185)

xolotl/core/include/xolotl/core/network/impl/Reaction.tpp

@@ -1843,10 +1843,14 @@ ProductionReaction<TNetwork, TDerived>::computeTableOne(
 	ConcentrationsView concentrations, std::vector<std::vector<IndexType>> clusterBins, std::vector<std::vector<double>>& rates, IndexType gridIndex)
 {
 	// Check what bins our clusters are in
-	auto binProduct = (_products[0] == invalidIndex) ?  whichBin(clusterBins,_products[0]) : whichBin(clusterBins,_products[1]);
+	auto binProduct = (_products[0] == invalidIndex) ?  whichBin(clusterBins,_products[1]) : whichBin(clusterBins,_products[0]);
 	auto binReactant0 = whichBin(clusterBins,_reactants[0]);
 	auto binReactant1 = whichBin(clusterBins,_reactants[1]);
 
+	auto productId = (_products[0] == invalidIndex) ? (_products[1]) : (_products[0]);
+
+
+
 	auto rateTableOne = this->_rate(gridIndex) * concentrations[_reactants[1]] *
 			this->_coefs(0, 0, 0, 0) * concentrations[_reactants[0]];
 
@@ -1856,6 +1860,12 @@ ProductionReaction<TNetwork, TDerived>::computeTableOne(
 		rates[binProduct][binReactant0] += rateTableOne;
 		rates[binProduct][binReactant1] += rateTableOne;
 	}
+	/*
+	std::cout << "product ID " << productId <<" product bin " << binProduct<<std::endl
+		<< "reactant0 ID " << _reactants[0]<<" reactant0 bin " << binReactant0<<std::endl
+		<< "reactant1 ID " << _reactants[1]<<" reactant1 bin " << binReactant1<<std::endl
+		<< "rate " << rateTableOne<<std::endl
+		<< std::endl;*/
 
 	return 0.0;
 }
@@ -1866,7 +1876,7 @@ double
 ProductionReaction<TNetwork, TDerived>::computeTableTwo(
 	ConcentrationsView concentrations, std::vector<std::vector<IndexType>> clusterBins, std::vector<std::vector<double>>& rates, IndexType gridIndex)
 {
-	auto binProduct = (_products[0] == invalidIndex) ?  whichBin(clusterBins,_products[0]) : whichBin(clusterBins,_products[1]);
+	auto binProduct = (_products[0] == invalidIndex) ?  whichBin(clusterBins,_products[1]) : whichBin(clusterBins,_products[0]);
 	auto binReactant0 = whichBin(clusterBins,_reactants[0]);
 	auto binReactant1 = whichBin(clusterBins,_reactants[1]);
 

xolotl/core/include/xolotl/core/network/impl/ReactionNetwork.tpp

@@ -636,11 +636,12 @@ ReactionNetwork<TImpl>::getTableOne(
 {
 	// Set up vector for rates of TableOne reactions
 	std::vector<std::vector<double>> rates (8, std::vector<double> (8,0));
+	auto ratesRef = std::ref(rates);
 
 	// Loop on all the rates
 	_reactions.forEach(
 		"ReactionNetwork::getTableOne",DEVICE_LAMBDA(auto&& reaction) {
-			reaction.contributeTableOne(concentrations, clusterBins, rates, gridIndex);
+			reaction.contributeTableOne(concentrations, clusterBins, ratesRef, gridIndex);
 		});
 	Kokkos::fence();
 	return rates;
@@ -653,11 +654,12 @@ ReactionNetwork<TImpl>::getTableTwo(
 {
 	// Set up vector for rates of TableTwo reactions
 	std::vector<std::vector<double>> rates (8, std::vector<double> (8,0));
+	auto ratesRef = std::ref(rates);
 
 	// Loop on all the rates
 	_reactions.forEach(
 		"ReactionNetwork::getTableTwo",DEVICE_LAMBDA(auto&& reaction) {
-			reaction.contributeTableTwo(concentrations, clusterBins, rates, gridIndex);
+			reaction.contributeTableTwo(concentrations, clusterBins, ratesRef, gridIndex);
 		});
 	Kokkos::fence();
 	return rates;
@@ -670,11 +672,12 @@ ReactionNetwork<TImpl>::getTableThree(
 {
 	// Set up vector for rates of TableThree reactions
 	std::vector<std::vector<double>> rates (8, std::vector<double> (8,0));
+	auto ratesRef = std::ref(rates);
 
 	// Loop on all the rates
 	_reactions.forEach(
 		"ReactionNetwork::getTableThree",DEVICE_LAMBDA(auto&& reaction) {
-			reaction.contributeTableThree(concentrations, clusterBins, rates, gridIndex);
+			reaction.contributeTableThree(concentrations, clusterBins, ratesRef, gridIndex);
 		});
 	Kokkos::fence();
 	return rates;
@@ -687,13 +690,20 @@ ReactionNetwork<TImpl>::getTableFour(
 {
 	// Set up vector for rates of TableFour reactions
 	std::vector<std::vector<double>> rates (8, std::vector<double> (8,0));
+	auto ratesRef = std::ref(rates);
 
 	// Loop on all the rates
 	_reactions.forEach(
 		"ReactionNetwork::getTableFour",DEVICE_LAMBDA(auto&& reaction) {
-			reaction.contributeTableOne(concentrations, clusterBins, rates, gridIndex);
+			reaction.contributeTableFour(concentrations, clusterBins, ratesRef, gridIndex);
 		});
 	Kokkos::fence();
+	for (auto i = 0; i < rates.size(); i++){
+			for (auto j = 0; j < rates[i].size(); j++){
+				std::cout << rates[i][j]<<" ";
+			}
+			std::cout << std::endl;
+	}
 	return rates;
 }
 

xolotl/core/src/network/ZrReactionNetwork.cpp

@@ -381,9 +381,9 @@ ZrReactionNetwork::addRxnDataValues(Kokkos::View<const double*> conc,
 
 
 	auto tableOne = getTableOne(conc, clusterBins, 0); //8x8 vector of reaction rates for TableOne reaction 
-	auto tableTwo = getTableTwo(conc, clusterBins, 0);
-	auto tableThree = getTableOne(conc, clusterBins, 0);
-	auto tableFour = getTableOne(conc, clusterBins, 0);
+	auto tableTwo = getTableTwo(conc, clusterBins, 0); //8x8 vector of reaction rates for TableTwo reaction 
+	auto tableThree = getTableThree(conc, clusterBins, 0); //8x8 vector of reaction rates for TableThree reaction 
+	auto tableFour = getTableFour(conc, clusterBins, 0); //8x8 vector of reaction rates for TableFour reaction 
 
 	for (auto i=0;i<totalVals.size();i++){
 		totalVals[i][0] = i;
@@ -396,7 +396,20 @@ ZrReactionNetwork::addRxnDataValues(Kokkos::View<const double*> conc,
 
 		}
 	} 
-
+	// Open the output file
+	/*
+	const int outputPrecision = 5;
+		std::fstream outputFile;
+		outputFile.open(
+			getRxnOutputFileName(), std::fstream::out | std::fstream::app);
+		outputFile << std::setprecision(outputPrecision);
+		
+		for (auto i = 0; i < clusterBins.size(); i++){
+			for (auto j = 0; j < clusterBins[i].size(); j++){
+				outputFile << clusterBins[i][j]<<" ";
+			}
+			outputFile << std::endl;
+		}*/
 }
 
 void
@@ -424,7 +437,7 @@ ZrReactionNetwork::writeRxnDataLine(
 	*/
 		// Set the output precision
 		const int outputPrecision = 5;
-
+		
 		// Open the output file
 		std::fstream outputFile;
 		outputFile.open(
@@ -441,8 +454,8 @@ ZrReactionNetwork::writeRxnDataLine(
 					   << iSize <<" "<< localData[8*i+3] << " "
 					   << localData[8*i+4] << " "<< localData[8*i+5] << " "
 						 << localData[8*i+6]<<" "<<localData[8*i+7]  << " "<< std::endl;
-		}*/
-
+		}//
+		*/
 		for (auto i = 0; i < localData.size(); i++){
 			for (auto j = 0; j < localData[i].size(); j++){
 				outputFile << localData[i][j] <<" ";