Commit dc8ecaff authored by COULAUD Olivier's avatar COULAUD Olivier

Merge branch 'master' of git+ssh://scm.gforge.inria.fr//gitroot/scalfmm/scalfmm

# By Berenger Bramas (5) and others
# Via piacibel (2) and Berenger Bramas (1)
* 'master' of git+ssh://scm.gforge.inria.fr//gitroot/scalfmm/scalfmm:
  Changes relative to int -> FSize propagated inside API
  Add a test function in order to test huge case without small verifications
  Add an allreduce to fmpi
  Add the inheritance of FAbstractSendable on each cell and update the taylor cell to support mpi stuffs
  Update the quick start
  update mpi part to work with fsize and add some assertions
  Provided new distribution of particles (hyperbolic paraboloid Z=C*(AX*AX-BY*BY)).
parents 8ff30c57 40313cfb
......@@ -16,6 +16,8 @@
#ifndef CKERNELAPI_H
#define CKERNELAPI_H
typedef long long int FSize;
/**
* @file
* This file defines the API for the C USER. The objective is to
......@@ -362,7 +364,7 @@ void scalfmm_algorithm_config(scalfmm_handle Handle,scalfmm_algorithm config);
* @param particleIndexes indexes of particles currently computed
* @param userData datas specific to the user's kernel
*/
typedef void (*Callback_P2M)(void* leafCell, FSize nbParticles, const int* particleIndexes, void* userData);
typedef void (*Callback_P2M)(void* leafCell, FSize nbParticles, const FSize* particleIndexes, void* userData);
/**
* @brief Function to be filled by user's M2M
......@@ -411,7 +413,7 @@ typedef void (*Callback_L2L)(int level, void* parentCell, int childPosition, voi
* @param particleIndexes indexes of particles currently computed
* @param userData datas specific to the user's kernel
*/
typedef void (*Callback_L2P)(void* leafCell, FSize nbParticles,const int* particleIndexes, void* userData);
typedef void (*Callback_L2P)(void* leafCell, FSize nbParticles,const FSize* particleIndexes, void* userData);
/**
* @brief Function to be filled by user's P2P
......@@ -421,7 +423,7 @@ typedef void (*Callback_L2P)(void* leafCell, FSize nbParticles,const int* partic
* @param sourceParticleIndexes indexes of cource particles currently computed
* @param userData datas specific to the user's kernel
*/
typedef void (*Callback_P2P)(FSize nbParticles, const int* particleIndexes, int nbSourceParticles, const int* sourceParticleIndexes, void* userData);
typedef void (*Callback_P2P)(FSize nbParticles, const FSize* particleIndexes, FSize nbSourceParticles, const FSize* sourceParticleIndexes, void* userData);
/**
* @brief Function to be filled by user's P2P
......@@ -432,8 +434,8 @@ typedef void (*Callback_P2P)(FSize nbParticles, const int* particleIndexes, int
* @param sourceNbPart array containing the number of part in each neighbors
* @param userData datas specific to the user's kernel
*/
typedef void (*Callback_P2PFull)(FSize nbParticles, const int* particleIndexes,
const int * sourceParticleIndexes[27],int sourceNbPart[27], void* userData);
typedef void (*Callback_P2PFull)(FSize nbParticles, const FSize* particleIndexes,
const FSize * sourceParticleIndexes[27],FSize sourceNbPart[27], void* userData);
/**
......@@ -442,7 +444,7 @@ typedef void (*Callback_P2PFull)(FSize nbParticles, const int* particleIndexes,
* @param particleIndexes indexes of particles currently computed
* @param userData datas specific to the user's kernel
*/
typedef void (*Callback_P2PInner)(FSize nbParticles, const int* particleIndexes, void* userData);
typedef void (*Callback_P2PInner)(FSize nbParticles, const FSize* particleIndexes, void* userData);
......
......@@ -113,7 +113,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
sources->getPhysicalValues()[idxPart] = physicalValues[indexes[idxPart]];
}
......@@ -127,7 +127,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
int iterPart = 0;
bool notFoundYet = true;
......@@ -149,7 +149,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
physicalValues[indexes[idxPart]] = sources->getPhysicalValues()[idxPart];
}
......@@ -161,7 +161,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
int iterPart = 0;
bool notFoundYet = true;
......@@ -182,7 +182,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
forcesToFill[indexes[idxPart]*3+0] = sources->getForcesX()[idxPart];
forcesToFill[indexes[idxPart]*3+1] = sources->getForcesY()[idxPart];
......@@ -195,7 +195,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
int iterPart = 0;
bool notFoundYet = true;
......@@ -218,7 +218,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
fX[indexes[idxPart]] = sources->getForcesX()[idxPart];
fY[indexes[idxPart]] = sources->getForcesY()[idxPart];
......@@ -231,7 +231,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
int iterPart = 0;
bool notFoundYet = true;
......@@ -255,7 +255,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
sources->getForcesX()[idxPart] = forcesToRead[indexes[idxPart]*3+0];
sources->getForcesY()[idxPart] = forcesToRead[indexes[idxPart]*3+1];
......@@ -267,7 +267,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
int iterPart = 0;
bool notFoundYet = true;
......@@ -289,7 +289,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
sources->getForcesX()[idxPart] = fX[indexes[idxPart]];
sources->getForcesY()[idxPart] = fY[indexes[idxPart]];
......@@ -301,7 +301,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
int iterPart = 0;
bool notFoundYet = true;
......@@ -325,7 +325,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
sources->getPotentials()[idxPart] = potentialsToRead[indexes[idxPart]];
}
......@@ -339,7 +339,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
int iterPart = 0;
bool notFoundYet = true;
......@@ -361,7 +361,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
potentialsToFill[indexes[idxPart]] = sources->getPotentials()[idxPart];
}
......@@ -373,7 +373,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
int iterPart = 0;
bool notFoundYet = true;
......@@ -394,7 +394,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
positionsToFill[indexes[idxPart]*3+0] = sources->getPositions()[0][idxPart];
positionsToFill[indexes[idxPart]*3+1] = sources->getPositions()[1][idxPart];
......@@ -407,7 +407,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
int iterPart = 0;
bool notFoundYet = true;
......@@ -430,7 +430,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
X[indexes[idxPart]] = sources->getPositions()[0][idxPart];
Y[indexes[idxPart]] = sources->getPositions()[1][idxPart];
......@@ -443,7 +443,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
int iterPart = 0;
bool notFoundYet = true;
......@@ -469,7 +469,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
sources->getWPositions()[0][idxPart] += updatedXYZ[indexes[idxPart]*3+0];
sources->getWPositions()[1][idxPart] += updatedXYZ[indexes[idxPart]*3+1];
......@@ -482,7 +482,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
sources->getWPositions()[0][idxPart] += X[indexes[idxPart]];
sources->getWPositions()[1][idxPart] += Y[indexes[idxPart]];
......@@ -496,7 +496,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
sources->getWPositions()[0][idxPart] = updatedXYZ[indexes[idxPart]*3+0];
sources->getWPositions()[1][idxPart] = updatedXYZ[indexes[idxPart]*3+1];
......@@ -509,7 +509,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
sources->getWPositions()[0][idxPart] = X[indexes[idxPart]];
sources->getWPositions()[1][idxPart] = Y[indexes[idxPart]];
......
......@@ -222,8 +222,8 @@ struct ScalFmmCoreHandle {
struct ScalFmmCoreConfig {
// Read/Write parameter
int treeHeight; // Number of level in the octree
FReal boxWidth; // Simulation box size (root level)
FPoint<FReal> boxCenter; // Center position of the box simulation(FReal[3])
double boxWidth; // Simulation box size (root level)
FPoint<double> boxCenter; // Center position of the box simulation(FReal[3])
};
ScalFmmCoreConfig config;
......@@ -391,5 +391,27 @@ extern "C" void scalfmm_user_kernel_config(scalfmm_handle Handle, Scalfmm_Kernel
((ScalFmmCoreHandle * ) Handle)->engine->user_kernel_config(userKernel,userDatas);
}
/**
* These functions are just translating functions.
*/
//< This function fill the childFullPosition[3] with [-1;1] to know the position of a child relatively to
//< its position from its parent
extern "C" void scalfmm_utils_parentChildPosition(int childPosition, int* childFullPosition){
childFullPosition[2] = (childPosition%2 ? 1 : -1);
childFullPosition[1] = ((childPosition/2)%2 ? 1 : -1);
childFullPosition[0] = ((childPosition/4)%2 ? 1 : -1);
}
//< This function fill the childFullPosition[3] with [-3;3] to know the position of a interaction
//< cell relatively to its position from the target
extern "C" void scalfmm_utils_interactionPosition(int interactionPosition, int* srcPosition){
srcPosition[2] = interactionPosition%7 - 3;
srcPosition[1] = (interactionPosition/7)%7 - 3;
srcPosition[0] = (interactionPosition/49)%7 - 3;
}
#endif
......@@ -9,6 +9,7 @@ extern "C" {
extern "C" scalfmm_handle scalfmm_init(/*int TreeHeight,double BoxWidth,double* BoxCenter, */scalfmm_kernel_type KernelType){
ScalFmmCoreHandle * handle = new ScalFmmCoreHandle();
typedef double FReal;
switch(KernelType){
case 0:
handle->engine = new FUserKernelEngine<FReal>(/*TreeHeight, BoxWidth, BoxCenter, */KernelType);
......@@ -24,16 +25,16 @@ extern "C" scalfmm_handle scalfmm_init(/*int TreeHeight,double BoxWidth,double*
handle->engine = new FInterEngine<FReal,ChebCell,ChebKernel>(/*TreeHeight,BoxWidth,BoxCenter, */KernelType);
break;
// case 2:
// //TODO typedefs
// typedef FP2PParticleContainerIndexed<FReal> ContainerClass;
// typedef FUnifCell<7> UnifCell;
// case 2:
// //TODO typedefs
// typedef FP2PParticleContainerIndexed<FReal> ContainerClass;
// typedef FUnifCell<7> UnifCell;
// typedef FInterpMatrixKernelR<FReal> MatrixKernelClass;
// typedef FUnifKernel<UnifCell,ContainerClass,MatrixKernelClass,7> UnifKernel;
// typedef FInterpMatrixKernelR<FReal> MatrixKernelClass;
// typedef FUnifKernel<UnifCell,ContainerClass,MatrixKernelClass,7> UnifKernel;
// handle->engine = new FInterEngine<UnifCell,UnifKernel>(/*TreeHeight,BoxWidth,BoxCenter, */KernelType);
// break;
// handle->engine = new FInterEngine<UnifCell,UnifKernel>(/*TreeHeight,BoxWidth,BoxCenter, */KernelType);
// break;
default:
std::cout<< "Kernel type unsupported" << std::endl;
......@@ -126,15 +127,15 @@ extern "C" void ChebKernelStruct_free(void *inKernel){
}
extern "C" void ChebKernel_P2M(void * leafCell, FSize nbParticles, const int *particleIndexes, void * inKernel){
extern "C" void ChebKernel_P2M(void * leafCell, FSize nbParticles, const FSize *particleIndexes, void * inKernel){
//make temporary array of parts
FP2PParticleContainerIndexed<double>* tempContainer = new FP2PParticleContainerIndexed<double>();
tempContainer->reserve(nbParticles);
FPoint<double> pos;
for(int i=0 ; i<nbParticles ; ++i){
pos = FPoint<double>(reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3 ],
reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3+1],
reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3+2]);
reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3+1],
reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3+2]);
double Phi = reinterpret_cast<UserData *>(inKernel)->myPhyValues[particleIndexes[i]];
tempContainer->push(pos,particleIndexes[i],Phi);
......@@ -168,8 +169,8 @@ extern "C" void ChebKernel_M2M(int level, void* parentCell, int childPosition,
//Get the kernel
ChebKernelStruct * inKernelStruct = reinterpret_cast<UserData*>(inKernel)->kernelStruct;
inKernelStruct->kernel[id_thread]->getPtrToInterpolator()->applyM2M(childPosition,
childChebCell->getMultipole(0),
parentChebCell->getMultipole(0));
childChebCell->getMultipole(0),
parentChebCell->getMultipole(0));
}
extern "C" void ChebKernel_M2L(int level, void* targetCell, void* sourceCell[343], void* inKernel){
......@@ -198,7 +199,7 @@ extern "C" void ChebKernel_M2L(int level, void* targetCell, void* sourceCell[34
}
extern "C" void ChebKernel_L2L(int level, void* parentCell, int childPosition, void* childCell, void* inKernel){
//Get our structures
//Get our structures
ChebCellStruct * parentCellStruct = reinterpret_cast<ChebCellStruct *>(parentCell);
ChebCellStruct * childCellStruct = reinterpret_cast<ChebCellStruct *>(childCell);
//get real cheb cell
......@@ -211,22 +212,22 @@ extern "C" void ChebKernel_L2L(int level, void* parentCell, int childPosition, v
//Get the kernel
ChebKernelStruct * inKernelStruct = reinterpret_cast<UserData*>(inKernel)->kernelStruct;
inKernelStruct->kernel[id_thread]->getPtrToInterpolator()->applyL2L(childPosition,
parentChebCell->getLocal(0),
childChebCell->getLocal(0));
parentChebCell->getLocal(0),
childChebCell->getLocal(0));
}
extern "C" void ChebKernel_L2P(void* leafCell, FSize nbParticles, const int* particleIndexes, void* inKernel){
extern "C" void ChebKernel_L2P(void* leafCell, FSize nbParticles, const FSize* particleIndexes, void* inKernel){
//Create temporary FSimpleLeaf
FP2PParticleContainerIndexed<double>* tempContainer = new FP2PParticleContainerIndexed<double>();
tempContainer->reserve(nbParticles);
FPoint<double> pos;
for(int i=0 ; i<nbParticles ; ++i){
pos = FPoint<double>(reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3 ],
reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3+1],
reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3+2]);
double Phi = reinterpret_cast<UserData *>(inKernel)->myPhyValues[particleIndexes[i]];
tempContainer->push(pos,particleIndexes[i],Phi);
}
pos = FPoint<double>(reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3 ],
reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3+1],
reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3+2]);
double Phi = reinterpret_cast<UserData *>(inKernel)->myPhyValues[particleIndexes[i]];
tempContainer->push(pos,particleIndexes[i],Phi);
}
//Get our structures
ChebCellStruct * leafCellStruct = reinterpret_cast<ChebCellStruct *>(leafCell);
//get real cheb cell
......@@ -242,28 +243,28 @@ extern "C" void ChebKernel_L2P(void* leafCell, FSize nbParticles, const int* par
//Then retrieve the results
double * forcesToFill = reinterpret_cast<UserData *>(inKernel)->forcesComputed[id_thread];
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
const FVector<FSize>& indexes = tempContainer->getIndexes();
for(FSize idxPart = 0 ; idxPart<nbParticles ; ++idxPart){
forcesToFill[indexes[idxPart]*3+0] += tempContainer->getForcesX()[idxPart];
forcesToFill[indexes[idxPart]*3+1] += tempContainer->getForcesY()[idxPart];
forcesToFill[indexes[idxPart]*3+2] += tempContainer->getForcesZ()[idxPart];
}
forcesToFill[indexes[idxPart]*3+0] += tempContainer->getForcesX()[idxPart];
forcesToFill[indexes[idxPart]*3+1] += tempContainer->getForcesY()[idxPart];
forcesToFill[indexes[idxPart]*3+2] += tempContainer->getForcesZ()[idxPart];
}
delete tempContainer;
tempContainer=nullptr;
}
void ChebKernel_P2P(FSize nbParticles, const int* particleIndexes,
const int * sourceParticleIndexes[27],int sourceNbPart[27],void* inKernel){
void ChebKernel_P2P(FSize nbParticles, const FSize* particleIndexes,
const FSize * sourceParticleIndexes[27],FSize sourceNbPart[27],void* inKernel){
//Create temporary FSimpleLeaf for target
FP2PParticleContainerIndexed<double>* tempContTarget = new FP2PParticleContainerIndexed<double>();
tempContTarget->reserve(nbParticles);
for(int i=0 ; i<nbParticles ; ++i){
FPoint<double> pos = FPoint<double>(reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3 ],
reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3+1],
reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3+2]);
reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3+1],
reinterpret_cast<UserData *>(inKernel)->insertedPositions[particleIndexes[i]*3+2]);
double Phi = reinterpret_cast<UserData *>(inKernel)->myPhyValues[particleIndexes[i]];
tempContTarget->push(pos,particleIndexes[i],Phi);
}
......@@ -277,12 +278,12 @@ extern "C" void ChebKernel_L2P(void* leafCell, FSize nbParticles, const int* par
//Allocate memory
tempContSources[idSource]->reserve(sourceNbPart[idSource]);
//Store a ptr to the indices of that source leaf
const int * indSource = sourceParticleIndexes[idSource];
const FSize * indSource = sourceParticleIndexes[idSource];
//Then, for each part in this source
for(int i=0 ; i<sourceNbPart[idSource] ; ++i){
FPoint<double> pos = FPoint<double>(reinterpret_cast<UserData *>(inKernel)->insertedPositions[indSource[i]*3 ],
reinterpret_cast<UserData *>(inKernel)->insertedPositions[indSource[i]*3+1],
reinterpret_cast<UserData *>(inKernel)->insertedPositions[indSource[i]*3+2]);
reinterpret_cast<UserData *>(inKernel)->insertedPositions[indSource[i]*3+1],
reinterpret_cast<UserData *>(inKernel)->insertedPositions[indSource[i]*3+2]);
double Phi = reinterpret_cast<UserData *>(inKernel)->myPhyValues[indSource[i]];
tempContSources[idSource]->push(pos,indSource[i],Phi);
}
......@@ -306,7 +307,7 @@ extern "C" void ChebKernel_L2P(void* leafCell, FSize nbParticles, const int* par
//get back forces
double * forcesToFill = reinterpret_cast<UserData *>(inKernel)->forcesComputed[id_thread];
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
const FVector<FSize>& indexes = tempContTarget->getIndexes();
for(FSize idxPart = 0 ; idxPart<nbParticles ; ++idxPart){
forcesToFill[indexes[idxPart]*3+0] += tempContTarget->getForcesX()[idxPart];
forcesToFill[indexes[idxPart]*3+1] += tempContTarget->getForcesY()[idxPart];
......@@ -316,12 +317,12 @@ extern "C" void ChebKernel_L2P(void* leafCell, FSize nbParticles, const int* par
//Note that sources are also modified.
//get back sources forces
for(int idSource = 0 ; idSource < 27 ; ++idSource){
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
const int nbPartInSource = sourceNbPart[idSource];
const FVector<FSize>& indexesSource = tempContSources[idSource]->getIndexes();
const FSize nbPartInSource = sourceNbPart[idSource];
for(int idxSourcePart = 0; idxSourcePart < nbPartInSource ; ++idxSourcePart){
forcesToFill[indexes[idxSourcePart]*3+0] += tempContSources[idSource]->getForcesX()[idxSourcePart];
forcesToFill[indexes[idxSourcePart]*3+1] += tempContSources[idSource]->getForcesY()[idxSourcePart];
forcesToFill[indexes[idxSourcePart]*3+2] += tempContSources[idSource]->getForcesZ()[idxSourcePart];
forcesToFill[indexesSource[idxSourcePart]*3+0] += tempContSources[idSource]->getForcesX()[idxSourcePart];
forcesToFill[indexesSource[idxSourcePart]*3+1] += tempContSources[idSource]->getForcesY()[idxSourcePart];
forcesToFill[indexesSource[idxSourcePart]*3+2] += tempContSources[idSource]->getForcesZ()[idxSourcePart];
}
}
......
......@@ -87,6 +87,7 @@ Scalfmm_Cell_Descriptor CoreCell::user_cell_descriptor;
*/
template< class CellClass, class ContainerClass>
class CoreKernel : public FAbstractKernels<CellClass,ContainerClass> {
Scalfmm_Kernel_Descriptor kernel;
void* userData;
......@@ -165,8 +166,8 @@ public:
if(kernel.p2p_full){
//Create the arrays of size and indexes
int nbPartPerNeighbors[27];
const int * indicesPerNeighbors[27];
FSize nbPartPerNeighbors[27];
const FSize * indicesPerNeighbors[27];
for(int idx=0 ; idx<27 ; ++idx){
if(neighbors[idx]){
nbPartPerNeighbors[idx] = neighbors[idx]->getNbParticles();
......@@ -277,7 +278,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
positionsToFill[indexes[idxPart]*3+0] = sources->getPositions()[0][idxPart];
positionsToFill[indexes[idxPart]*3+1] = sources->getPositions()[1][idxPart];
......@@ -290,7 +291,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
int iterPart = 0;
bool notFoundYet = true;
......@@ -313,7 +314,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
X[indexes[idxPart]] = sources->getPositions()[0][idxPart];
Y[indexes[idxPart]] = sources->getPositions()[1][idxPart];
......@@ -326,7 +327,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
int iterPart = 0;
bool notFoundYet = true;
......@@ -353,7 +354,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
sources->getWPositions()[0][idxPart] += updatedXYZ[indexes[idxPart]*3+0];
sources->getWPositions()[1][idxPart] += updatedXYZ[indexes[idxPart]*3+1];
......@@ -366,7 +367,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
sources->getWPositions()[0][idxPart] += X[indexes[idxPart]];
sources->getWPositions()[1][idxPart] += Y[indexes[idxPart]];
......@@ -380,7 +381,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
sources->getWPositions()[0][idxPart] = updatedXYZ[indexes[idxPart]*3+0];
sources->getWPositions()[1][idxPart] = updatedXYZ[indexes[idxPart]*3+1];
......@@ -393,7 +394,7 @@ public:
octree->forEachLeaf([&](LeafClass* leaf){
ContainerClass * sources = leaf->getSrc();
const FVector<FSize>& indexes = leaf->getTargets()->getIndexes();
int nbPartThere = sources->getNbParticles();
FSize nbPartThere = sources->getNbParticles();
for(FSize idxPart = 0 ; idxPart<nbPartThere ; ++idxPart){
sources->getWPositions()[0][idxPart] = X[indexes[idxPart]];