Skip to content
GitLab
Projects
Groups
Snippets
Help
Loading...
Help
What's new
7
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Open sidebar
solverstack
ScalFMM
Commits
c93d2bb7
Commit
c93d2bb7
authored
Jan 14, 2015
by
PIACIBELLO Cyrille
Browse files
Options
Browse Files
Download
Email Patches
Plain Diff
Operator Timers added to ChebyshevInterpolation example
parent
5d4e332e
Changes
1
Hide whitespace changes
Inline
Side-by-side
Showing
1 changed file
with
133 additions
and
124 deletions
+133
-124
Examples/ChebyshevInterpolationFMM.cpp
Examples/ChebyshevInterpolationFMM.cpp
+133
-124
No files found.
Examples/ChebyshevInterpolationFMM.cpp
View file @
c93d2bb7
...
...
@@ -68,7 +68,7 @@ int main(int argc, char* argv[])
FParameterDefinitions
::
NbThreads
);
const
std
::
string
defaultFile
(
/*SCALFMMDataPath+*/
"../Data/
test20k.
fma"
);
const
std
::
string
defaultFile
(
/*SCALFMMDataPath+*/
"../Data/
unitCubeXYZQ100.b
fma"
);
const
std
::
string
filename
=
FParameters
::
getStr
(
argc
,
argv
,
FParameterDefinitions
::
InputFile
.
options
,
defaultFile
.
c_str
());
const
unsigned
int
TreeHeight
=
FParameters
::
getValue
(
argc
,
argv
,
FParameterDefinitions
::
OctreeHeight
.
options
,
5
);
const
unsigned
int
SubTreeHeight
=
FParameters
::
getValue
(
argc
,
argv
,
FParameterDefinitions
::
OctreeSubHeight
.
options
,
2
);
...
...
@@ -76,133 +76,142 @@ int main(int argc, char* argv[])
#ifdef _OPENMP
omp_set_num_threads
(
NbThreads
);
std
::
cout
<<
"
\n
>> Using "
<<
omp_get_max_threads
()
<<
" threads.
\n
"
<<
std
::
endl
;
omp_set_num_threads
(
NbThreads
);
std
::
cout
<<
"
\n
>> Using "
<<
omp_get_max_threads
()
<<
" threads.
\n
"
<<
std
::
endl
;
#else
std
::
cout
<<
"
\n
>> Sequential version.
\n
"
<<
std
::
endl
;
std
::
cout
<<
"
\n
>> Sequential version.
\n
"
<<
std
::
endl
;
#endif
//
std
::
cout
<<
"Parameters "
<<
std
::
endl
<<
" Octree Depth "
<<
TreeHeight
<<
std
::
endl
<<
" SubOctree depth "
<<
SubTreeHeight
<<
std
::
endl
<<
" Input file name: "
<<
filename
<<
std
::
endl
<<
" Thread number: "
<<
NbThreads
<<
std
::
endl
<<
std
::
endl
;
//
// init timer
FTic
time
;
// open particle file
////////////////////////////////////////////////////////////////////
//
FFmaGenericLoader
loader
(
filename
);
//
////////////////////////////////////////////////////////////////////
// begin Chebyshev kernel
// accuracy
const
unsigned
int
ORDER
=
7
;
// typedefs
typedef
FP2PParticleContainerIndexed
<>
ContainerClass
;
typedef
FSimpleLeaf
<
ContainerClass
>
LeafClass
;
typedef
FChebCell
<
ORDER
>
CellClass
;
typedef
FOctree
<
CellClass
,
ContainerClass
,
LeafClass
>
OctreeClass
;
//
typedef
FInterpMatrixKernelR
MatrixKernelClass
;
//
std
::
cout
<<
"Parameters "
<<
std
::
endl
<<
" Octree Depth "
<<
TreeHeight
<<
std
::
endl
<<
" SubOctree depth "
<<
SubTreeHeight
<<
std
::
endl
<<
" Input file name: "
<<
filename
<<
std
::
endl
<<
" Thread number: "
<<
NbThreads
<<
std
::
endl
<<
std
::
endl
;
//
// init timer
FTic
time
;
// open particle file
////////////////////////////////////////////////////////////////////
//
FFmaGenericLoader
loader
(
filename
);
//
////////////////////////////////////////////////////////////////////
// begin Chebyshev kernel
// accuracy
const
unsigned
int
ORDER
=
7
;
// typedefs
typedef
FP2PParticleContainerIndexed
<>
ContainerClass
;
typedef
FSimpleLeaf
<
ContainerClass
>
LeafClass
;
typedef
FChebCell
<
ORDER
>
CellClass
;
typedef
FOctree
<
CellClass
,
ContainerClass
,
LeafClass
>
OctreeClass
;
//
typedef
FInterpMatrixKernelR
MatrixKernelClass
;
const
MatrixKernelClass
MatrixKernel
;
typedef
FChebSymKernel
<
CellClass
,
ContainerClass
,
MatrixKernelClass
,
ORDER
>
KernelClass
;
//
typedef
FChebSymKernel
<
CellClass
,
ContainerClass
,
MatrixKernelClass
,
ORDER
>
KernelClass
;
//
#ifdef _OPENMP
typedef
FFmmAlgorithmThread
<
OctreeClass
,
CellClass
,
ContainerClass
,
KernelClass
,
LeafClass
>
FmmClass
;
typedef
FFmmAlgorithmThread
<
OctreeClass
,
CellClass
,
ContainerClass
,
KernelClass
,
LeafClass
>
FmmClass
;
#else
typedef
FFmmAlgorithm
<
OctreeClass
,
CellClass
,
ContainerClass
,
KernelClass
,
LeafClass
>
FmmClass
;
typedef
FFmmAlgorithm
<
OctreeClass
,
CellClass
,
ContainerClass
,
KernelClass
,
LeafClass
>
FmmClass
;
#endif
// init oct-tree
OctreeClass
tree
(
TreeHeight
,
SubTreeHeight
,
loader
.
getBoxWidth
(),
loader
.
getCenterOfBox
());
{
// -----------------------------------------------------
std
::
cout
<<
"Creating & Inserting "
<<
loader
.
getNumberOfParticles
()
<<
" particles ..."
<<
std
::
endl
;
std
::
cout
<<
"
\t
Height : "
<<
TreeHeight
<<
"
\t
sub-height : "
<<
SubTreeHeight
<<
std
::
endl
;
time
.
tic
();
//
FPoint
position
;
FReal
physicalValue
=
0.0
;
//
for
(
int
idxPart
=
0
;
idxPart
<
loader
.
getNumberOfParticles
()
;
++
idxPart
){
//
// Read particle per particle from file
loader
.
fillParticle
(
&
position
,
&
physicalValue
);
//
// put particle in octree
tree
.
insert
(
position
,
idxPart
,
physicalValue
);
}
time
.
tac
();
std
::
cout
<<
"Done "
<<
"(@Creating and Inserting Particles = "
<<
time
.
elapsed
()
<<
" s) ."
<<
std
::
endl
;
}
// -----------------------------------------------------
{
// -----------------------------------------------------
std
::
cout
<<
"
\n
Chebyshev FMM (ORDER="
<<
ORDER
<<
") ... "
<<
std
::
endl
;
time
.
tic
();
//
KernelClass
kernels
(
TreeHeight
,
loader
.
getBoxWidth
(),
loader
.
getCenterOfBox
(),
&
MatrixKernel
);
//
FmmClass
algorithm
(
&
tree
,
&
kernels
);
//
algorithm
.
execute
();
// Here the call of the FMM algorithm
//
time
.
tac
();
std
::
cout
<<
"Done "
<<
"(@Algorithm = "
<<
time
.
elapsed
()
<<
" s) ."
<<
std
::
endl
;
}
// -----------------------------------------------------
//
// Some output
//
//
{
// -----------------------------------------------------
long
int
N1
=
0
,
N2
=
loader
.
getNumberOfParticles
()
/
2
,
N3
=
loader
.
getNumberOfParticles
()
-
1
;
;
FReal
energy
=
0.0
;
//
// Loop over all leaves
//
std
::
cout
<<
std
::
endl
<<
" &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& "
<<
std
::
endl
;
std
::
cout
<<
std
::
scientific
;
std
::
cout
.
precision
(
10
)
;
tree
.
forEachLeaf
([
&
](
LeafClass
*
leaf
){
const
FReal
*
const
posX
=
leaf
->
getTargets
()
->
getPositions
()[
0
];
const
FReal
*
const
posY
=
leaf
->
getTargets
()
->
getPositions
()[
1
];
const
FReal
*
const
posZ
=
leaf
->
getTargets
()
->
getPositions
()[
2
];
const
FReal
*
const
potentials
=
leaf
->
getTargets
()
->
getPotentials
();
const
FReal
*
const
forcesX
=
leaf
->
getTargets
()
->
getForcesX
();
const
FReal
*
const
forcesY
=
leaf
->
getTargets
()
->
getForcesY
();
const
FReal
*
const
forcesZ
=
leaf
->
getTargets
()
->
getForcesZ
();
const
int
nbParticlesInLeaf
=
leaf
->
getTargets
()
->
getNbParticles
();
const
FReal
*
const
physicalValues
=
leaf
->
getTargets
()
->
getPhysicalValues
();
const
FVector
<
int
>&
indexes
=
leaf
->
getTargets
()
->
getIndexes
();
for
(
int
idxPart
=
0
;
idxPart
<
nbParticlesInLeaf
;
++
idxPart
){
const
int
indexPartOrig
=
indexes
[
idxPart
];
if
((
indexPartOrig
==
N1
)
||
(
indexPartOrig
==
N2
)
||
(
indexPartOrig
==
N3
)
)
{
std
::
cout
<<
"Index "
<<
indexPartOrig
<<
" potential "
<<
potentials
[
idxPart
]
<<
" Pos "
<<
posX
[
idxPart
]
<<
" "
<<
posY
[
idxPart
]
<<
" "
<<
posZ
[
idxPart
]
<<
" Forces: "
<<
forcesX
[
idxPart
]
<<
" "
<<
forcesY
[
idxPart
]
<<
" "
<<
forcesZ
[
idxPart
]
<<
std
::
endl
;
}
energy
+=
potentials
[
idxPart
]
*
physicalValues
[
idxPart
]
;
}
});
std
::
cout
<<
std
::
endl
<<
"Energy: "
<<
energy
<<
std
::
endl
;
std
::
cout
<<
std
::
endl
<<
" &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& "
<<
std
::
endl
<<
std
::
endl
;
}
// -----------------------------------------------------
return
0
;
// init oct-tree
OctreeClass
tree
(
TreeHeight
,
SubTreeHeight
,
loader
.
getBoxWidth
(),
loader
.
getCenterOfBox
());
{
// -----------------------------------------------------
std
::
cout
<<
"Creating & Inserting "
<<
loader
.
getNumberOfParticles
()
<<
" particles ..."
<<
std
::
endl
;
std
::
cout
<<
"
\t
Height : "
<<
TreeHeight
<<
"
\t
sub-height : "
<<
SubTreeHeight
<<
std
::
endl
;
time
.
tic
();
//
FPoint
position
;
FReal
physicalValue
=
0.0
;
//
for
(
int
idxPart
=
0
;
idxPart
<
loader
.
getNumberOfParticles
()
;
++
idxPart
){
//
// Read particle per particle from file
loader
.
fillParticle
(
&
position
,
&
physicalValue
);
//
// put particle in octree
tree
.
insert
(
position
,
idxPart
,
physicalValue
);
}
time
.
tac
();
std
::
cout
<<
"Done "
<<
"(@Creating and Inserting Particles = "
<<
time
.
elapsed
()
<<
" s) ."
<<
std
::
endl
;
}
// -----------------------------------------------------
{
// -----------------------------------------------------
std
::
cout
<<
"
\n
Chebyshev FMM (ORDER="
<<
ORDER
<<
") ... "
<<
std
::
endl
;
time
.
tic
();
//
KernelClass
kernels
(
TreeHeight
,
loader
.
getBoxWidth
(),
loader
.
getCenterOfBox
(),
&
MatrixKernel
);
//
FmmClass
algo
(
&
tree
,
&
kernels
);
//
algo
.
execute
();
// Here the call of the FMM algorithm
//
time
.
tac
();
std
::
cout
<<
"Timers Far Field
\n
"
<<
"P2M "
<<
algo
.
getTime
(
FAlgorithmTimers
::
P2MTimer
)
<<
" seconds
\n
"
<<
"M2M "
<<
algo
.
getTime
(
FAlgorithmTimers
::
M2MTimer
)
<<
" seconds
\n
"
<<
"M2L "
<<
algo
.
getTime
(
FAlgorithmTimers
::
M2LTimer
)
<<
" seconds
\n
"
<<
"L2L "
<<
algo
.
getTime
(
FAlgorithmTimers
::
L2LTimer
)
<<
" seconds
\n
"
<<
"P2P and L2P "
<<
algo
.
getTime
(
FAlgorithmTimers
::
NearTimer
)
<<
" seconds
\n
"
<<
std
::
endl
;
std
::
cout
<<
"Done "
<<
"(@Algorithm = "
<<
time
.
elapsed
()
<<
" s) ."
<<
std
::
endl
;
}
// -----------------------------------------------------
//
// Some output
//
//
{
// -----------------------------------------------------
long
int
N1
=
0
,
N2
=
loader
.
getNumberOfParticles
()
/
2
,
N3
=
loader
.
getNumberOfParticles
()
-
1
;
;
FReal
energy
=
0.0
;
//
// Loop over all leaves
//
std
::
cout
<<
std
::
endl
<<
" &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& "
<<
std
::
endl
;
std
::
cout
<<
std
::
scientific
;
std
::
cout
.
precision
(
10
)
;
tree
.
forEachLeaf
([
&
](
LeafClass
*
leaf
){
const
FReal
*
const
posX
=
leaf
->
getTargets
()
->
getPositions
()[
0
];
const
FReal
*
const
posY
=
leaf
->
getTargets
()
->
getPositions
()[
1
];
const
FReal
*
const
posZ
=
leaf
->
getTargets
()
->
getPositions
()[
2
];
const
FReal
*
const
potentials
=
leaf
->
getTargets
()
->
getPotentials
();
const
FReal
*
const
forcesX
=
leaf
->
getTargets
()
->
getForcesX
();
const
FReal
*
const
forcesY
=
leaf
->
getTargets
()
->
getForcesY
();
const
FReal
*
const
forcesZ
=
leaf
->
getTargets
()
->
getForcesZ
();
const
int
nbParticlesInLeaf
=
leaf
->
getTargets
()
->
getNbParticles
();
const
FReal
*
const
physicalValues
=
leaf
->
getTargets
()
->
getPhysicalValues
();
const
FVector
<
int
>&
indexes
=
leaf
->
getTargets
()
->
getIndexes
();
for
(
int
idxPart
=
0
;
idxPart
<
nbParticlesInLeaf
;
++
idxPart
){
const
int
indexPartOrig
=
indexes
[
idxPart
];
if
((
indexPartOrig
==
N1
)
||
(
indexPartOrig
==
N2
)
||
(
indexPartOrig
==
N3
)
)
{
std
::
cout
<<
"Index "
<<
indexPartOrig
<<
" potential "
<<
potentials
[
idxPart
]
<<
" Pos "
<<
posX
[
idxPart
]
<<
" "
<<
posY
[
idxPart
]
<<
" "
<<
posZ
[
idxPart
]
<<
" Forces: "
<<
forcesX
[
idxPart
]
<<
" "
<<
forcesY
[
idxPart
]
<<
" "
<<
forcesZ
[
idxPart
]
<<
std
::
endl
;
}
energy
+=
potentials
[
idxPart
]
*
physicalValues
[
idxPart
]
;
}
});
std
::
cout
<<
std
::
endl
<<
"Energy: "
<<
energy
<<
std
::
endl
;
std
::
cout
<<
std
::
endl
<<
" &&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&& "
<<
std
::
endl
<<
std
::
endl
;
}
// -----------------------------------------------------
return
0
;
}
Write
Preview
Markdown
is supported
0%
Try again
or
attach a new file
.
Attach a file
Cancel
You are about to add
0
people
to the discussion. Proceed with caution.
Finish editing this message first!
Cancel
Please
register
or
sign in
to comment