diff --git a/Src/Kernels/Rotation/FRotationCell.hpp b/Src/Kernels/Rotation/FRotationCell.hpp
index 593f96535964fed23c3476b51d939a7eb0bc001b..7585e1c896a69e54112a5b7c3a6d493b23cd3b84 100644
--- a/Src/Kernels/Rotation/FRotationCell.hpp
+++ b/Src/Kernels/Rotation/FRotationCell.hpp
@@ -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());
}
};
diff --git a/Src/Kernels/Rotation/FRotationKernel.hpp b/Src/Kernels/Rotation/FRotationKernel.hpp
index 513a7ab6b78f6fdf968282b841700937f8c4d9cb..7b0163706a86d86abc34c6c2acbde79ba16f8bd2 100644
--- a/Src/Kernels/Rotation/FRotationKernel.hpp
+++ b/Src/Kernels/Rotation/FRotationKernel.hpp
@@ -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();