Skip to content
GitLab
Projects
Groups
Snippets
/
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Menu
Open sidebar
SWIFT
SWIFTsim
Commits
f8be5419
Commit
f8be5419
authored
Jun 23, 2017
by
Matthieu Schaller
Browse files
Updated the cell sanitization tool to give more sensible values to the extreme h outliers.
parent
2983cf56
Changes
8
Hide whitespace changes
Inline
Side-by-side
examples/main.c
View file @
f8be5419
...
...
@@ -418,6 +418,8 @@ int main(int argc, char *argv[]) {
parser_get_param_string
(
params
,
"InitialConditions:file_name"
,
ICfileName
);
const
int
replicate
=
parser_get_opt_param_int
(
params
,
"InitialConditions:replicate"
,
1
);
const
int
clean_h_values
=
parser_get_opt_param_int
(
params
,
"InitialConditions:cleanup_h"
,
0
);
if
(
myrank
==
0
)
message
(
"Reading ICs from file '%s'"
,
ICfileName
);
fflush
(
stdout
);
...
...
@@ -619,7 +621,7 @@ int main(int argc, char *argv[]) {
#endif
/* Initialise the particles */
engine_init_particles
(
&
e
,
flag_entropy_ICs
);
engine_init_particles
(
&
e
,
flag_entropy_ICs
,
clean_h_values
);
/* Write the state of the system before starting time integration. */
engine_dump_snapshot
(
&
e
);
...
...
examples/parameter_example.yml
View file @
f8be5419
...
...
@@ -60,6 +60,7 @@ Gravity:
# Parameters related to the initial conditions
InitialConditions
:
file_name
:
SedovBlast/sedov.hdf5
# The file to read
cleanup_h
:
0
# (Optional) Clean the values of h that are read in. Set to 1 to activate.
h_scaling
:
1.
# (Optional) A scaling factor to apply to all smoothing lengths in the ICs.
shift_x
:
0.
# (Optional) A shift to apply to all particles read from the ICs (in internal units).
shift_y
:
0.
...
...
src/cell.c
View file @
f8be5419
...
...
@@ -941,53 +941,52 @@ void cell_split(struct cell *c, ptrdiff_t parts_offset, ptrdiff_t sparts_offset,
* @brief Sanitizes the smoothing length values of cells by setting large
* outliers to more sensible values.
*
*
We compute the mean and standard deviation of the smoothing lengths in
*
logarithmic space and limit values to mean + 4 sigma
.
*
Each cell with <1000 part will be processed. We limit h to be the size of
*
the cell and replace 0s with a good estimate
.
*
* @param c The cell.
* @param treated Has the cell already been sanitized at this level ?
*/
void
cell_sanitize
(
struct
cell
*
c
)
{
void
cell_sanitize
(
struct
cell
*
c
,
int
treated
)
{
const
int
count
=
c
->
count
;
struct
part
*
parts
=
c
->
parts
;
float
h_max
=
0
.
f
;
/* First collect some statistics */
float
h_mean
=
0
.
f
,
h_mean2
=
0
.
f
;
float
h_min
=
FLT_MAX
,
h_max
=
0
.
f
;
for
(
int
i
=
0
;
i
<
count
;
++
i
)
{
/* Treat cells will <1000 particles */
if
(
count
<
1000
&&
!
treated
)
{
const
float
h
=
logf
(
parts
[
i
].
h
);
h_mean
+=
h
;
h_mean2
+=
h
*
h
;
h_max
=
max
(
h_max
,
h
);
h_min
=
min
(
h_min
,
h
);
}
h_mean
/=
count
;
h_mean2
/=
count
;
const
float
h_var
=
h_mean2
-
h_mean
*
h_mean
;
const
float
h_std
=
(
h_var
>
0
.
f
)
?
sqrtf
(
h_var
)
:
0
.
1
f
*
h_mean
;
/* Choose a cut */
const
float
h_limit
=
expf
(
h_mean
+
4
.
f
*
h_std
);
/* Be verbose this is not innocuous */
message
(
"Cell properties: h_min= %f h_max= %f geometric mean= %f."
,
expf
(
h_min
),
expf
(
h_max
),
expf
(
h_mean
));
/* Get an upper bound on h */
const
float
upper_h_max
=
c
->
dmin
/
(
1
.
2
f
*
kernel_gamma
);
if
(
c
->
h_max
>
h_limit
)
{
/* Apply it */
for
(
int
i
=
0
;
i
<
count
;
++
i
)
{
if
(
parts
[
i
].
h
==
0
.
f
||
parts
[
i
].
h
>
upper_h_max
)
parts
[
i
].
h
=
upper_h_max
;
}
}
message
(
"Smoothing lengths will be limited to (mean + 4sigma)= %f."
,
h_lim
it
)
;
/* Recurse and gather the new h_max values */
if
(
c
->
spl
it
)
{
/* Apply the cut */
for
(
int
i
=
0
;
i
<
count
;
++
i
)
parts
->
h
=
min
(
parts
[
i
].
h
,
h_limit
);
for
(
int
k
=
0
;
k
<
8
;
++
k
)
{
if
(
c
->
progeny
[
k
]
!=
NULL
)
{
c
->
h_max
=
h_limit
;
/* Recurse */
cell_sanitize
(
c
->
progeny
[
k
],
(
count
<
1000
));
/* And collect */
h_max
=
max
(
h_max
,
c
->
progeny
[
k
]
->
h_max
);
}
}
}
else
{
message
(
"Smoothing lengths will not be limited."
);
/* Get the new value of h_max */
for
(
int
i
=
0
;
i
<
count
;
++
i
)
h_max
=
max
(
h_max
,
parts
[
i
].
h
);
}
/* Record the change */
c
->
h_max
=
h_max
;
}
/**
...
...
src/cell.h
View file @
f8be5419
...
...
@@ -344,7 +344,7 @@ struct cell {
void
cell_split
(
struct
cell
*
c
,
ptrdiff_t
parts_offset
,
ptrdiff_t
sparts_offset
,
struct
cell_buff
*
buff
,
struct
cell_buff
*
sbuff
,
struct
cell_buff
*
gbuff
);
void
cell_sanitize
(
struct
cell
*
c
);
void
cell_sanitize
(
struct
cell
*
c
,
int
treated
);
int
cell_locktree
(
struct
cell
*
c
);
void
cell_unlocktree
(
struct
cell
*
c
);
int
cell_glocktree
(
struct
cell
*
c
);
...
...
src/engine.c
View file @
f8be5419
...
...
@@ -2790,8 +2790,10 @@ void engine_print_task_counts(struct engine *e) {
* @brief Rebuild the space and tasks.
*
* @param e The #engine.
* @param clean_h_values Are we cleaning up the values of h before building
* the tasks ?
*/
void
engine_rebuild
(
struct
engine
*
e
)
{
void
engine_rebuild
(
struct
engine
*
e
,
int
clean_h_values
)
{
const
ticks
tic
=
getticks
();
...
...
@@ -2802,7 +2804,7 @@ void engine_rebuild(struct engine *e) {
space_rebuild
(
e
->
s
,
e
->
verbose
);
/* Initial cleaning up session ? */
if
(
e
->
s
->
sanitized
==
0
)
space_sanitize
(
e
->
s
);
if
(
clean_h_values
)
space_sanitize
(
e
->
s
);
/* If in parallel, exchange the cell structure. */
#ifdef WITH_MPI
...
...
@@ -2856,7 +2858,7 @@ void engine_prepare(struct engine *e) {
if
(
e
->
forcerepart
)
engine_repartition
(
e
);
/* Do we need rebuilding ? */
if
(
e
->
forcerebuild
)
engine_rebuild
(
e
);
if
(
e
->
forcerebuild
)
engine_rebuild
(
e
,
0
);
/* Unskip active tasks and check for rebuild */
engine_unskip
(
e
);
...
...
@@ -3218,9 +3220,12 @@ void engine_launch(struct engine *e, int nr_runners) {
*
* @param e The #engine
* @param flag_entropy_ICs Did the 'Internal Energy' of the particles actually
*contain entropy ?
* contain entropy ?
* @param clean_h_values Are we cleaning up the values of h before building
* the tasks ?
*/
void
engine_init_particles
(
struct
engine
*
e
,
int
flag_entropy_ICs
)
{
void
engine_init_particles
(
struct
engine
*
e
,
int
flag_entropy_ICs
,
int
clean_h_values
)
{
struct
space
*
s
=
e
->
s
;
...
...
@@ -3237,7 +3242,7 @@ void engine_init_particles(struct engine *e, int flag_entropy_ICs) {
}
/* Construct all cells and tasks to start everything */
engine_rebuild
(
e
);
engine_rebuild
(
e
,
clean_h_values
);
/* No time integration. We just want the density and ghosts */
engine_skip_force_and_kick
(
e
);
...
...
src/engine.h
View file @
f8be5419
...
...
@@ -272,7 +272,8 @@ void engine_init(struct engine *e, struct space *s,
struct
sourceterms
*
sourceterms
);
void
engine_launch
(
struct
engine
*
e
,
int
nr_runners
);
void
engine_prepare
(
struct
engine
*
e
);
void
engine_init_particles
(
struct
engine
*
e
,
int
flag_entropy_ICs
);
void
engine_init_particles
(
struct
engine
*
e
,
int
flag_entropy_ICs
,
int
clean_h_values
);
void
engine_step
(
struct
engine
*
e
);
void
engine_maketasks
(
struct
engine
*
e
);
void
engine_split
(
struct
engine
*
e
,
struct
partition
*
initial_partition
);
...
...
@@ -281,7 +282,7 @@ void engine_exchange_strays(struct engine *e, size_t offset_parts,
int
*
ind_gpart
,
size_t
*
Ngpart
,
size_t
offset_sparts
,
int
*
ind_spart
,
size_t
*
Nspart
);
void
engine_rebuild
(
struct
engine
*
e
);
void
engine_rebuild
(
struct
engine
*
e
,
int
clean_h_values
);
void
engine_repartition
(
struct
engine
*
e
);
void
engine_repartition_trigger
(
struct
engine
*
e
);
void
engine_makeproxies
(
struct
engine
*
e
);
...
...
src/space.c
View file @
f8be5419
...
...
@@ -962,28 +962,33 @@ void space_split(struct space *s, struct cell *cells, int nr_cells,
}
/**
* @brief Runs through the top-level cells and checks whether tasks associated
* with them can be split. If not, try to sanitize the cells.
* @brief #threadpool mapper function to sanitize the cells
*
* @param s The #space to act upon.
* @param map_data Pointers towards the top-level cells.
* @param num_cells The number of top-level cells.
* @param extra_data Unused parameters.
*/
void
space_sanitize
(
struct
space
*
s
)
{
s
->
sanitized
=
1
;
void
space_sanitize
_mapper
(
void
*
map_data
,
int
num_cells
,
void
*
extra_data
)
{
/* Unpack the inputs. */
s
truct
cell
*
cells_top
=
(
struct
cell
*
)
map_data
;
for
(
int
k
=
0
;
k
<
s
->
nr_cells
;
k
++
)
{
for
(
int
ind
=
0
;
ind
<
num_cells
;
ind
++
)
{
struct
cell
*
c
=
&
cells_top
[
ind
];
cell_sanitize
(
c
,
0
);
}
}
struct
cell
*
c
=
&
s
->
cells_top
[
k
];
const
double
min_width
=
c
->
dmin
;
/**
* @brief Runs through the top-level cells and sanitize their h values
*
* @param s The #space to act upon.
*/
void
space_sanitize
(
struct
space
*
s
)
{
/* Do we have a problem ? */
if
(
c
->
h_max
*
kernel_gamma
*
space_stretch
>
min_width
*
0
.
5
&&
c
->
count
>
space_maxcount
)
{
if
(
s
->
e
->
nodeID
==
0
)
message
(
"Cleaning up unreasonable values of h"
);
/* Ok, clean-up the mess */
cell_sanitize
(
c
);
}
}
threadpool_map
(
&
s
->
e
->
threadpool
,
space_sanitize_mapper
,
s
->
cells_top
,
s
->
nr_cells
,
sizeof
(
struct
cell
),
1
,
NULL
);
}
/**
...
...
@@ -2630,7 +2635,6 @@ void space_init(struct space *s, const struct swift_params *params,
s
->
dim
[
0
]
=
dim
[
0
];
s
->
dim
[
1
]
=
dim
[
1
];
s
->
dim
[
2
]
=
dim
[
2
];
s
->
sanitized
=
0
;
s
->
periodic
=
periodic
;
s
->
gravity
=
gravity
;
s
->
nr_parts
=
Npart
;
...
...
src/space.h
View file @
f8be5419
...
...
@@ -139,9 +139,6 @@ struct space {
/*! Number of queues in the system. */
int
nr_queues
;
/*! Has this space already been sanitized ? */
int
sanitized
;
/*! The associated engine. */
struct
engine
*
e
;
...
...
Write
Preview
Supports
Markdown
0%
Try again
or
attach a new 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