Mentions légales du service

Skip to content
Snippets Groups Projects
Commit 4b333838 authored by Quentin Khan's avatar Quentin Khan
Browse files

Update Rotation Kernel interface to use the new data organisation

 - Divide the FRotationCell inner data layout into two sub-types:
   multipole_t and local_expansion_t. Two class attribute are accessible
   through the `getMultipoledata` and `getLocalExpansionData` methods.

 - Change FMM operators signature to take advantage of the new
   layout. The operators only take as parameters the data they may need
   instead of the whole cell. For instance (simplified):

   void M2M(CellClass* parent, CellClass** children);

   becomes

   void M2M(multipole_t*  parent_m,    symbolic_data_t* parent_s,
            multipole_t** children_ms, symbolic_data_t* children_ss);
parent 4af27174
No related branches found
No related tags found
No related merge requests found
......@@ -4,13 +4,13 @@
// This software is a computer program whose purpose is to compute the FMM.
//
// This software is governed by the CeCILL-C and LGPL licenses and
// abiding by the rules of distribution of free software.
//
// abiding by the rules of distribution of free software.
//
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
// GNU General Public and CeCILL-C Licenses for more details.
// "http://www.cecill.info".
// "http://www.cecill.info".
// "http://www.gnu.org/licenses".
// ===================================================================================
#ifndef FROTATIONCELL_HPP
......@@ -37,75 +37,68 @@
template <class FReal, int P>
class FRotationCell : public FBasicCell, public FAbstractSendable {
protected:
//< Size of multipole vector
static const int MultipoleSize = ((P+2)*(P+1))/2; // Artimethique suite (n+1)*n/2
//< Size of local vector
static const int LocalSize = ((P+2)*(P+1))/2; // Artimethique suite (n+1)*n/2
//< Multipole vector (static memory)
FComplex<FReal> multipole_exp[MultipoleSize]; //< For multipole extenssion
//< Local vector (static memory)
FComplex<FReal> local_exp[LocalSize]; //< For local extenssion
template<std::size_t S, class Tag>
struct expansion_impl {
enum {Size = S};
FComplex<FReal> exp[Size];
FComplex<FReal>* get() noexcept {
return exp;
}
const FComplex<FReal>* get() const noexcept {
return exp;
}
int getSize() const noexcept {
return Size;
}
void reset() noexcept {
for(int idx = 0; idx < Size; ++idx) {
exp[idx].setRealImag(FReal(0.0), FReal(0.0));
}
}
template<class BufferWriterClass>
void serialize(BufferWriterClass& buffer) const {
buffer.write(exp, Size);
}
template<class BufferReaderClass>
void deserialize(BufferReaderClass& buffer) {
buffer.fillArray(exp, Size);
}
};
public:
/** Default constructor
* Put 0 in vectors
*/
FRotationCell(){
}
/** Copy constructor
* Copy the value in the vectors
*/
FRotationCell(const FRotationCell& other){
(*this) = other;
}
using multipole_t = expansion_impl<((P+2)*(P+1))/2, class multipole_tag>;
using local_expansion_t = expansion_impl<((P+2)*(P+1))/2, class local_expansion_tag>;
/** Default destructor */
virtual ~FRotationCell(){
}
protected:
/** Copy operator
* copies only the value in the vectors
*/
FRotationCell& operator=(const FRotationCell& other) {
FMemUtils::copyall(multipole_exp, other.multipole_exp, MultipoleSize);
FMemUtils::copyall(local_exp, other.local_exp, LocalSize);
return *this;
}
multipole_t m_data;
local_expansion_t l_data;
/** Get Multipole array */
const FComplex<FReal>* getMultipole() const {
return multipole_exp;
}
/** Get Local array */
const FComplex<FReal>* getLocal() const {
return local_exp;
}
public:
/** Get Multipole array */
FComplex<FReal>* getMultipole() {
return multipole_exp;
const multipole_t& getMultipoleData() const noexcept {
return m_data;
}
/** Get Local array */
FComplex<FReal>* getLocal() {
return local_exp;
multipole_t& getMultipoleData() {
return m_data;
}
int getArraySize() const
{
return MultipoleSize;
const local_expansion_t& getLocalExpansionData() const noexcept {
return l_data;
}
local_expansion_t& getLocalExpansionData() {
return l_data;
}
/** Make it like the begining */
void resetToInitialState(){
for(int idx = 0 ; idx < MultipoleSize ; ++idx){
multipole_exp[idx].setRealImag(FReal(0.0), FReal(0.0));
}
for(int idx = 0 ; idx < LocalSize ; ++idx){
local_exp[idx].setRealImag(FReal(0.0), FReal(0.0));
}
m_data.reset();
l_data.reset();
}
///////////////////////////////////////////////////////
......@@ -113,28 +106,28 @@ public:
///////////////////////////////////////////////////////
template <class BufferWriterClass>
void serializeUp(BufferWriterClass& buffer) const{
buffer.write(multipole_exp, MultipoleSize);
m_data.serialize(buffer);
}
template <class BufferReaderClass>
void deserializeUp(BufferReaderClass& buffer){
buffer.fillArray(multipole_exp, MultipoleSize);
m_data.deserialize(buffer);
}
template <class BufferWriterClass>
void serializeDown(BufferWriterClass& buffer) const{
buffer.write(local_exp, LocalSize);
l_data.serialize(buffer);
}
template <class BufferReaderClass>
void deserializeDown(BufferReaderClass& buffer){
buffer.fillArray(local_exp, LocalSize);
l_data.deserialize(buffer);
}
FSize getSavedSizeUp() const {
return ((FSize) sizeof(FComplex<FReal>)) * (MultipoleSize);
return ((FSize) sizeof(FComplex<FReal>)) * (multipole_t::Size);
}
FSize getSavedSizeDown() const {
return ((FSize) sizeof(FComplex<FReal>)) * (LocalSize);
return ((FSize) sizeof(FComplex<FReal>)) * (local_expansion_t::Size);
}
///////////////////////////////////////////////////////
......@@ -143,18 +136,18 @@ public:
template <class BufferWriterClass>
void save(BufferWriterClass& buffer) const{
FBasicCell::save(buffer);
buffer.write(multipole_exp, MultipoleSize);
buffer.write(local_exp, LocalSize);
m_data.serialize();
l_data.serialize();
}
template <class BufferReaderClass>
void restore(BufferReaderClass& buffer){
FBasicCell::restore(buffer);
buffer.fillArray(multipole_exp, MultipoleSize);
buffer.fillArray(local_exp, LocalSize);
m_data.deserialize();
l_data.deserialize();
}
FSize getSavedSize() const {
return FSize(((int) sizeof(FComplex<FReal>)) * (MultipoleSize + LocalSize)
return FSize(((int) sizeof(FComplex<FReal>)) * (multipole_t::Size + local_expansion_t::Size)
+ FBasicCell::getSavedSize());
}
};
......
......@@ -909,13 +909,17 @@ public:
* \omega (q,a) = q \frac{a^{l}}{(l+|m|)!} P_{lm}(cos( \alpha ) )e^{-im \beta}
* \f]
*/
void P2M(CellClass* const inPole, const ContainerClass* const inParticles ) override {
template<class SymbolicData>
void P2M(typename CellClass::multipole_t* const LeafMultipole,
const SymbolicData* const LeafSymb,
const ContainerClass* const inParticles )
{
const FReal i_pow_m[4] = {0, FMath::FPiDiv2<FReal>(), FMath::FPi<FReal>(), -FMath::FPiDiv2<FReal>()};
// w is the multipole moment
FComplex<FReal>* FRestrict const w = inPole->getMultipole();
FComplex<FReal>* FRestrict const w = LeafMultipole->get();
// Copying the position is faster than using cell position
const FPoint<FReal> cellPosition = getLeafCenter(inPole->getCoordinate());
const FPoint<FReal> cellPosition = getLeafCenter(LeafSymb->getCoordinate());
// We need a legendre array
FReal legendre[SizeArray];
......@@ -972,17 +976,23 @@ public:
* then transfer using the formula
* and finaly rotate back.
*/
void M2M(CellClass* const FRestrict inPole, const CellClass*const FRestrict *const FRestrict inChildren, const int inLevel) override {
template<class SymbolicData>
void M2M(typename CellClass::multipole_t* const FRestrict ParentMultipole,
const SymbolicData* const ParentSymb,
const typename CellClass::multipole_t*const FRestrict *const FRestrict ChildMultipoles,
const SymbolicData* const /*ChildSymbs*/[])
{
int level = static_cast<int>(ParentSymb->getLevel());
// Get the translation coef for this level (same for all child)
const FReal*const coef = M2MTranslationCoef[inLevel];
const FReal*const coef = M2MTranslationCoef[level];
// A buffer to copy the source w allocated once
FComplex<FReal> source_w[SizeArray];
// For all children
for(int idxChild = 0 ; idxChild < 8 ; ++idxChild){
// if child exists
if(inChildren[idxChild]){
if(ChildMultipoles[idxChild]){
// Copy the source
FMemUtils::copyall(source_w, inChildren[idxChild]->getMultipole(), SizeArray);
FMemUtils::copyall(source_w, ChildMultipoles[idxChild]->get(), SizeArray);
// rotate it forward
RotationZVectorsMul(source_w,rotationExpMinusImPhi[idxChild]);
......@@ -1012,7 +1022,7 @@ public:
RotationZVectorsMul(target_w,rotationExpImPhi[idxChild]);
// Sum the result
FMemUtils::addall( inPole->getMultipole(), target_w, SizeArray);
FMemUtils::addall( ParentMultipole->get(), target_w, SizeArray);
}
}
}
......@@ -1028,16 +1038,20 @@ public:
* then transfer using the formula
* and finaly rotate back.
*/
void M2L(CellClass* const FRestrict inLocal, const CellClass* inInteractions[], const int /*inSize*/, const int inLevel) {
void M2L(typename CellClass::local_expansion_t* const FRestrict NodeLocalExpansion,
const typename CellClass::multipole_t* const FRestrict NeighborMultipoles[],
const int /*inSize*/,
const int inLevel)
{
// To copy the multipole data allocated once
FComplex<FReal> source_w[SizeArray];
// For all children
for(int idxNeigh = 0 ; idxNeigh < 343 ; ++idxNeigh){
// if interaction exits
if(inInteractions[idxNeigh]){
if(NeighborMultipoles[idxNeigh]){
const FReal*const coef = M2LTranslationCoef[inLevel][idxNeigh];
// Copy multipole data into buffer
FMemUtils::copyall(source_w, inInteractions[idxNeigh]->getMultipole(), SizeArray);
FMemUtils::copyall(source_w, NeighborMultipoles[idxNeigh]->get(), SizeArray);
// Rotate
RotationZVectorsMul(source_w,rotationM2LExpMinusImPhi[idxNeigh]);
......@@ -1070,22 +1084,29 @@ public:
RotationZVectorsMul(target_u,rotationM2LExpMinusImPhi[idxNeigh]);
// Sum
FMemUtils::addall(inLocal->getLocal(), target_u, SizeArray);
FMemUtils::addall(NodeLocalExpansion->get(), target_u, SizeArray);
}
}
}
void M2L(CellClass* const FRestrict inLocal, const CellClass* inInteractions[],
const int neighborPositions[], const int inSize, const int inLevel) override {
template<class SymbolicData>
void M2L(typename CellClass::local_expansion_t * const FRestrict TargetExpansion,
const SymbolicData* const TargetSymb,
const typename CellClass::multipole_t * const FRestrict SourceMultipoles[],
const SymbolicData* const FRestrict /*SourceSymbs*/[],
const int neighborPositions[],
const int inSize)
{
int level = static_cast<int>(TargetSymb->getLevel());
// To copy the multipole data allocated once
FComplex<FReal> source_w[SizeArray];
// For all children
for(int idxExistingNeigh = 0 ; idxExistingNeigh < inSize ; ++idxExistingNeigh){
const int idxNeigh = neighborPositions[idxExistingNeigh];
// if interaction exits
const FReal*const coef = M2LTranslationCoef[inLevel][idxNeigh];
const FReal*const coef = M2LTranslationCoef[level][idxNeigh];
// Copy multipole data into buffer
FMemUtils::copyall(source_w, inInteractions[idxExistingNeigh]->getMultipole(), SizeArray);
FMemUtils::copyall(source_w, SourceMultipoles[idxExistingNeigh]->get(), SizeArray);
// Rotate
RotationZVectorsMul(source_w,rotationM2LExpMinusImPhi[idxNeigh]);
......@@ -1118,7 +1139,7 @@ public:
RotationZVectorsMul(target_u,rotationM2LExpMinusImPhi[idxNeigh]);
// Sum
FMemUtils::addall(inLocal->getLocal(), target_u, SizeArray);
FMemUtils::addall(TargetExpansion->get(), target_u, SizeArray);
}
}
......@@ -1133,17 +1154,23 @@ public:
* then transfer using the formula
* and finaly rotate back.
*/
void L2L(const CellClass* const FRestrict inLocal, CellClass* FRestrict *const FRestrict inChildren, const int inLevel) override {
template<class SymbolicData>
void L2L(const typename CellClass::local_expansion_t * const FRestrict ParentLocalExpansion,
const SymbolicData* const ParentSymb,
typename CellClass::local_expansion_t * FRestrict * const FRestrict ChildLocalExpansion,
const SymbolicData* const /*ChildSymbs*/[])
{
int level = static_cast<int>(ParentSymb->getLevel());
// Get the translation coef for this level (same for all chidl)
const FReal*const coef = L2LTranslationCoef[inLevel];
const FReal*const coef = L2LTranslationCoef[level];
// To copy the source local allocated once
FComplex<FReal> source_u[SizeArray];
// For all children
for(int idxChild = 0 ; idxChild < 8 ; ++idxChild){
// if child exists
if(inChildren[idxChild]){
if(ChildLocalExpansion[idxChild]){
// Copy the local data into the buffer
FMemUtils::copyall(source_u, inLocal->getLocal(), SizeArray);
FMemUtils::copyall(source_u, ParentLocalExpansion->get(), SizeArray);
// Rotate
RotationZVectorsMul(source_u,rotationExpImPhi[idxChild]);
......@@ -1172,7 +1199,7 @@ public:
RotationZVectorsMul(target_u,rotationExpMinusImPhi[idxChild]);
// Sum in child
FMemUtils::addall(inChildren[idxChild]->getLocal(), target_u, SizeArray);
FMemUtils::addall(ChildLocalExpansion[idxChild]->get(), target_u, SizeArray);
}
}
}
......@@ -1195,13 +1222,17 @@ public:
* F_{ \phi } = -\frac{1}{r sin \phi} \sum_{j=0}^P \sum_{k=1}^j{(-2k) Im(u_{j,k} I_{j,k}(r, \theta, \phi)) }
* \f]
*/
void L2P(const CellClass* const inLocal, ContainerClass* const inParticles) override {
template<class SymbolicData>
void L2P(const typename CellClass::local_expansion_t* const LeafLocalExpansion,
const SymbolicData* const LeafSymb,
ContainerClass* const inParticles)
{
const FReal i_pow_m[4] = {0, FMath::FPiDiv2<FReal>(), FMath::FPi<FReal>(), -FMath::FPiDiv2<FReal>()};
// Take the local value from the cell
const FComplex<FReal>* FRestrict const u = inLocal->getLocal();
const FComplex<FReal>* FRestrict const u = LeafLocalExpansion->get();
// Copying the position is faster than using cell position
const FPoint<FReal> cellPosition = getLeafCenter(inLocal->getCoordinate());
const FPoint<FReal> cellPosition = getLeafCenter(LeafSymb->getCoordinate());
// For all particles in the leaf box
const FReal*const physicalValues = inParticles->getPhysicalValues();
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment