/***************************************
	Auteur : Pierre Aubert
	Mail : aubertp7@gmail.com
	Licence : CeCILL-C
****************************************/

#include "phoenix_intrinsics.h"

#include <algorithm>
#include <execution>

#include "intrinsics_propagation.h"

#include "intrinsics_propagation_block.h"

///Propagate the U and V species in the matVecVecU and matVecV for blocks
/**	@param[out] outMatVecU : updated matrix U version (with vectorial neighbours)
 * 	@param[out] outMatVecV : updated matrix V version (with vectorial neighbours)
 * 	@param matVecVecU : input of matrix U (with vectorial neighbours)
 * 	@param matVecV : input of matrix V (with vectorial neighbours)
 * 	@param matBroadcastDeltaSquare : matrix of the delta square values (with broadcast neighbours)
 * 	@param nbStencilRow : number of rows of the matrix matBroadcastDeltaSquare
 * 	@param nbStencilCol : number of columns of the matrix matBroadcastDeltaSquare
 * 	@param diffusionRateU : diffusion rate of the U specie
 * 	@param diffudionRateV : diffusion rate of the V specie
 * 	@param feedRate : rate of the process which feeds U and drains U, V and P
 * 	@param killRate : rate of the process which converts V into P
 * 	@param dt : time interval between two steps
*/
void grayscott_propagation_compute_block(PBlock<float> & outMatVecU, PBlock<float> & outMatVecV, const PBlock<float> & matVecVecU, const PBlock<float> & matVecVecV,
					const float * matBroadcastDeltaSquare, long nbStencilRow, long nbStencilCol,
					float diffusionRateU, float diffusionRateV, float feedRate, float killRate, float dt)
{
	size_t nbRow(outMatVecU.getFullNbRow()), nbCol(outMatVecU.getNbCol());
	grayscott_propagation(outMatVecU.getData(), outMatVecV.getData(), matVecVecU.getData(), matVecVecV.getData(), nbRow, nbCol,
				matBroadcastDeltaSquare, nbStencilRow, nbStencilCol,
				diffusionRateU, diffusionRateV, feedRate, killRate, dt);
}


void grayscott_propagation_compute_block_dummy(PBlock<float> & outMatVecU, PBlock<float> & outMatVecV,
// 					const PBlock<float> & matVecVecU, const PBlock<float> & matVecVecV,
					const float * matBroadcastDeltaSquare, long nbStencilRow, long nbStencilCol,
					float diffusionRateU, float diffusionRateV, float feedRate, float killRate, float dt)
{
	
	
}

///Propagate the U and V species in the matVecVecU and matVecV
/**	@param[out] outMatVecU : updated matrix U version (with vectorial neighbours)
 * 	@param[out] outMatVecV : updated matrix V version (with vectorial neighbours)
 * 	@param matVecVecU : input of matrix U (with vectorial neighbours)
 * 	@param matVecV : input of matrix V (with vectorial neighbours)
 * 	@param[out] vecBlockOutU : vector of output blocks for U
 * 	@param[out] vecBlockOutV : vector of output blocks for V
 * 	@param[out] vecBlockOutU : vector of input blocks for U
 * 	@param[out] vecBlockOutV : vector of input blocks for V
 * 	@param blockNbRow : maximum number of rows per block
 * 	@param blockNbCol : maximum number of columns per block
 * 	@param matBroadcastDeltaSquare : matrix of the delta square values (with broadcast neighbours)
 * 	@param nbStencilRow : number of rows of the matrix matBroadcastDeltaSquare
 * 	@param nbStencilCol : number of columns of the matrix matBroadcastDeltaSquare
 * 	@param diffusionRateU : diffusion rate of the U specie
 * 	@param diffudionRateV : diffusion rate of the V specie
 * 	@param feedRate : rate of the process which feeds U and drains U, V and P
 * 	@param killRate : rate of the process which converts V into P
 * 	@param dt : time interval between two steps
*/
void grayscott_propagation_block(PTensor<float> & outMatVecU, PTensor<float> & outMatVecV, const PTensor<float> & matVecVecU, const PTensor<float> & matVecVecV,
				std::vector<PBlock<float> > & vecBlockOutU, std::vector<PBlock<float> > & vecBlockOutV,
				std::vector<PBlock<float> > & vecBlockInU, std::vector<PBlock<float> > & vecBlockInV,
				size_t blockNbRow, size_t blockNbCol,
				const float * matBroadcastDeltaSquare, long nbStencilRow, long nbStencilCol,
				float diffusionRateU, float diffusionRateV, float feedRate, float killRate, float dt)
{
	outMatVecU.splitBlock(vecBlockOutU, blockNbRow, blockNbCol, 1lu);
	outMatVecV.splitBlock(vecBlockOutV, blockNbRow, blockNbCol, 1lu);
	matVecVecU.splitBlock(vecBlockInU, blockNbRow, blockNbCol, 1lu);
	matVecVecV.splitBlock(vecBlockInV, blockNbRow, blockNbCol, 1lu);
	
	std::vector<size_t> vecIndex;
	for(size_t i(0lu); i < vecBlockOutU.size(); ++i){vecIndex.push_back(i);}

	std::for_each(std::execution::seq, vecIndex.begin(), vecIndex.end(),
			[&](size_t i){
				grayscott_propagation_compute_block((PBlock<float>&)vecBlockOutU[i], (PBlock<float>&)vecBlockOutV[i],
								vecBlockInU[i], vecBlockInV[i],
								matBroadcastDeltaSquare, nbStencilRow, nbStencilCol,
								diffusionRateU, diffusionRateV, feedRate, killRate, dt);
			});

	
	outMatVecU.mergeBlock(vecBlockOutU, 1lu);
	outMatVecV.mergeBlock(vecBlockOutV, 1lu);
}