Commit 824bbc34 authored by DUFOYER Benjamin's avatar DUFOYER Benjamin

Factorisation of validation method

parent 94de3af9
#ifndef _VALIDATION_METHOD_
#define _VALIDATION_METHOD_
namespace validation_methods
{
/**
* validate_result function to check the current group tree
* @author benjamin.dufoyer@inria.fr
* @param loader loader of the file used in
* the programm
* @param TreeHeight height of the group tree
* @param SubTreeHeight size of the subTree
* @param myParticles container of particle
* @param mpiComm MPI comm
* @param groupedTree groupTree
*/
template
< class ContainerClass,
class LeafClass,
class CellClass,
class OctreeClass,
class KernelClass,
class FmmClass,
class LoaderClass,
class ParticleContainer,
class MatrixKernelClass,
class GroupTreeClass,
class GroupCellSymbClass,
class GroupCellUpClass,
class GroupCellDownClass>
void validate_result(
LoaderClass loader,
int TreeHeight,
int SubTreeHeight,
ParticleContainer myParticles,
FMpi mpiComm,
GroupTreeClass groupedTree
)
{
const FReal epsi = 1E-10;
OctreeClass treeCheck(TreeHeight, SubTreeHeight,loader.getBoxWidth(),loader.getCenterOfBox());
for(FSize idxPart = 0 ; idxPart < myParticles.getSize() ; ++idxPart){
// put in tree
treeCheck.insert(myParticles[idxPart].position,
myParticles[idxPart].physicalValue);
}
MatrixKernelClass MatrixKernel;
KernelClass kernels(TreeHeight, loader.getBoxWidth(), loader.getCenterOfBox(), &MatrixKernel);
FmmClass algorithm(mpiComm.global(),&treeCheck, &kernels);
algorithm.execute();
std::cout << "Algo is over" << std::endl;
groupedTree.forEachCellWithLevel(
[&](GroupCellSymbClass* gsymb ,
GroupCellUpClass* gmul,
GroupCellDownClass* gloc,
const int level)
{
const CellClass* cell = treeCheck.getCell(gsymb->getMortonIndex(), level);
if(cell == nullptr){
std::cout << "[Empty] Error cell should exist " << gsymb->getMortonIndex() << "\n";
}
else {
FMath::FAccurater<FReal> diffUp;
diffUp.add(cell->getMultipoleData().get(0), gmul->get(0), gmul->getVectorSize());
if(diffUp.getRelativeInfNorm() > epsi || diffUp.getRelativeL2Norm() > epsi){
std::cout << "[Up] Up is different at index " << gsymb->getMortonIndex() << " level " << level << " is " << diffUp << "\n";
}
FMath::FAccurater<FReal> diffDown;
diffDown.add(cell->getLocalExpansionData().get(0), gloc->get(0), gloc->getVectorSize());
if(diffDown.getRelativeInfNorm() > epsi || diffDown.getRelativeL2Norm() > epsi){
std::cout << "[Up] Down is different at index " << gsymb->getMortonIndex() << " level " << level << " is " << diffDown << "\n";
}
}
});
groupedTree.forEachCellLeaf(
[&](GroupCellSymbClass* gsymb ,
GroupCellUpClass* /* gmul */,
GroupCellDownClass* /* gloc */,
FP2PGroupParticleContainer<FReal> * leafTarget)
{
const ContainerClass* targets = treeCheck.getLeafSrc(gsymb->getMortonIndex());
if(targets == nullptr){
std::cout << "[Empty] Error leaf should exist " << gsymb->getMortonIndex() << "\n";
}
else{
const FReal*const gposX = leafTarget->getPositions()[0];
const FReal*const gposY = leafTarget->getPositions()[1];
const FReal*const gposZ = leafTarget->getPositions()[2];
const FSize gnbPartsInLeafTarget = leafTarget->getNbParticles();
const FReal*const gforceX = leafTarget->getForcesX();
const FReal*const gforceY = leafTarget->getForcesY();
const FReal*const gforceZ = leafTarget->getForcesZ();
const FReal*const gpotential = leafTarget->getPotentials();
const FReal*const posX = targets->getPositions()[0];
const FReal*const posY = targets->getPositions()[1];
const FReal*const posZ = targets->getPositions()[2];
const FSize nbPartsInLeafTarget = targets->getNbParticles();
const FReal*const forceX = targets->getForcesX();
const FReal*const forceY = targets->getForcesY();
const FReal*const forceZ = targets->getForcesZ();
const FReal*const potential = targets->getPotentials();
if(gnbPartsInLeafTarget != nbPartsInLeafTarget){
std::cout << "[Empty] Not the same number of particles at " << gsymb->getMortonIndex()
<< " gnbPartsInLeafTarget " << gnbPartsInLeafTarget << " nbPartsInLeafTarget " << nbPartsInLeafTarget << "\n";
}
else{
FMath::FAccurater<FReal> potentialDiff;
FMath::FAccurater<FReal> fx, fy, fz;
for(FSize idxPart = 0 ; idxPart < nbPartsInLeafTarget ; ++idxPart){
if(gposX[idxPart] != posX[idxPart] || gposY[idxPart] != posY[idxPart]
|| gposZ[idxPart] != posZ[idxPart]){
std::cout << "[Empty] Not the same particlea at " << gsymb->getMortonIndex() << " idx " << idxPart
<< gposX[idxPart] << " " << posX[idxPart] << " " << gposY[idxPart] << " " << posY[idxPart]
<< " " << gposZ[idxPart] << " " << posZ[idxPart] << "\n";
}
else{
potentialDiff.add(potential[idxPart], gpotential[idxPart]);
fx.add(forceX[idxPart], gforceX[idxPart]);
fy.add(forceY[idxPart], gforceY[idxPart]);
fz.add(forceZ[idxPart], gforceZ[idxPart]);
}
}
if(potentialDiff.getRelativeInfNorm() > epsi || potentialDiff.getRelativeL2Norm() > epsi){
std::cout << "[Up] potentialDiff is different at index " << gsymb->getMortonIndex() << " is " << potentialDiff << "\n";
}
if(fx.getRelativeInfNorm() > epsi || fx.getRelativeL2Norm() > epsi){
std::cout << "[Up] fx is different at index " << gsymb->getMortonIndex() << " is " << fx << "\n";
}
if(fy.getRelativeInfNorm() > epsi || fy.getRelativeL2Norm() > epsi){
std::cout << "[Up] fy is different at index " << gsymb->getMortonIndex() << " is " << fy << "\n";
}
if(fz.getRelativeInfNorm() > epsi || fz.getRelativeL2Norm() > epsi){
std::cout << "[Up] fz is different at index " << gsymb->getMortonIndex() << " is " << fz << "\n";
}
}
}
});
std::cout << "Comparing is over" << std::endl;
}
}
#endif //_VALIDATION_METHOD_
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment