#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>
#include <time.h>
#include <mpi.h>
#include "allvars.h"
#include "proto.h"
Go to the source code of this file.
Defines | |
| #define | NTAB 1000 |
| #define | NEAREST(x) (((x)>boxhalf)?((x)-boxsize):(((x)<-boxhalf)?((x)+boxsize):(x))) |
| #define | EN 64 |
Functions | |
| int | force_treebuild (int npart) |
| int | force_treebuild_single (int npart) |
| void | force_create_empty_nodes (int no, int topnode, int bits, int x, int y, int z, int *nodecount, int *nextfree) |
| void | force_insert_pseudo_particles (void) |
| void | force_update_node_recursive (int no, int sib, int father) |
| void | force_update_pseudoparticles (void) |
| void | force_exchange_pseudodata (void) |
| void | force_treeupdate_pseudos (void) |
| void | force_flag_localnodes (void) |
| void | force_update_len (void) |
| void | force_update_node_len_local (void) |
| void | force_update_node_len_toptree (void) |
| void | force_update_hmax (void) |
| void | force_update_node_hmax_local (void) |
| void | force_update_node_hmax_toptree (void) |
| int | force_treeevaluate (int target, int mode, double *ewaldcountsum) |
| int | force_treeevaluate_shortrange (int target, int mode) |
| int | force_treeevaluate_ewald_correction (int target, int mode, double pos_x, double pos_y, double pos_z, double aold) |
| void | force_treeevaluate_potential (int target, int mode) |
| void | force_treeevaluate_potential_shortrange (int target, int mode) |
| void | force_treeallocate (int maxnodes, int maxpart) |
| void | force_treefree (void) |
| int | force_treeevaluate_direct (int target, int mode) |
| void | dump_particles (void) |
| void | ewald_init (void) |
| void | ewald_corr (double dx, double dy, double dz, double *fper) |
| double | ewald_pot_corr (double dx, double dy, double dz) |
| double | ewald_psi (double x[3]) |
| void | ewald_force (int iii, int jjj, int kkk, double x[3], double force[3]) |
This file contains the computation of the gravitational force by means of a tree. The type of tree implemented is a geometrical oct-tree, starting from a cube encompassing all particles. This cube is automatically found in the domain decomposition, which also splits up the global "top-level" tree along node boundaries, moving the particles of different parts of the tree to separate processors. Tree nodes can be dynamically updated in drift/kick operations to avoid having to reconstruct the tree every timestep.
Definition in file forcetree.c.
|
|
Size of 3D lock-up table for Ewald correction force Definition at line 45 of file forcetree.c. Referenced by ewald_pot_corr(), and force_treeevaluate_ewald_correction(). |
|
|
Macro that maps a distance to the nearest periodic neighbour Definition at line 43 of file forcetree.c. Referenced by force_treeevaluate_ewald_correction(). |
|
|
length of lock-up table for short-range force kernel in TreePM algorithm Definition at line 31 of file forcetree.c. |
|
|
This function dumps some of the basic particle data to a file. In case the tree construction fails, it is called just before the run terminates with an error message. Examination of the generated file may then give clues to what caused the problem. Definition at line 2763 of file forcetree.c. Referenced by force_create_empty_nodes(), and force_treebuild_single(). |
|
||||||||||||||||||||
|
This function looks up the correction force due to the infinite number of periodic particle/node images. We here use trilinear interpolation to get it from the precomputed tables, which contain one octant around the target particle at the origin. The other octants are obtained from it by exploiting the symmetry properties. Definition at line 2951 of file forcetree.c. |
|
||||||||||||||||||||||||
|
This function computes the force correction term (difference between full force of infinite lattice and nearest image) by Ewald summation. Definition at line 3134 of file forcetree.c. |
|
|
This function initializes tables with the correction force and the correction potential due to the periodic images of a point mass located at the origin. These corrections are obtained by Ewald summation. (See Hernquist, Bouchet, Suto, ApJS, 1991, 75, 231) The correction fields are used to obtain the full periodic force if periodic boundaries combined with the pure tree algorithm are used. For the TreePM algorithm, the Ewald correction is not used. The correction fields are stored on disk once they are computed. If a corresponding file is found, they are loaded from disk to speed up the initialization. The Ewald summation is done in parallel, i.e. the processors share the work to compute the tables if needed. Definition at line 2802 of file forcetree.c. References ThisTask. Referenced by begrun(). |
|
||||||||||||||||
|
This function looks up the correction potential due to the infinite number of periodic particle/node images. We here use tri-linear interpolation to get it from the precomputed table, which contains one octant around the target particle at the origin. The other octants are obtained from it by exploiting symmetry properties. Definition at line 3040 of file forcetree.c. References EN. |
|
|
This function computes the potential correction term by means of Ewald summation. Definition at line 3093 of file forcetree.c. |
|
||||||||||||||||||||||||||||||||||||
|
This function recursively creates a set of empty tree nodes which corresponds to the top-level tree for the domain grid. This is done to ensure that this top-level tree is always "complete" so that we can easily associate the pseudo-particles of other CPUs with tree-nodes at a given level in the tree, even when the particle population is so sparse that some of these nodes are actually empty. Definition at line 294 of file forcetree.c. References NODE::center, topnode_data::Daughter, DomainNodeIndex, dump_particles(), endrun(), topnode_data::Leaf, NODE::len, MaxNodes, Nodes, peano_hilbert_key(), NODE::suns, ThisTask, TopNodes, and NODE::u. Referenced by force_treebuild_single(). |
|
|
This function communicates the values of the multipole moments of the top-level tree-nodes of the domain grid. This data can then be used to update the pseudo-particles on each CPU accordingly. Definition at line 687 of file forcetree.c. References NODE::bitflags, NODE::d, DomainMoment, DomainNodeIndex, Extnodes, DomainNODE::mass, NODE::mass, Nodes, DomainNODE::s, NODE::s, NODE::u, DomainNODE::vs, and extNODE::vs. |
|
|
This function flags nodes in the top-level tree that are dependent on local particle data. Definition at line 819 of file forcetree.c. References NODE::bitflags, NODE::d, DomainNodeIndex, NODE::father, Nodes, and NODE::u. Referenced by force_treebuild(). |
|
|
this function inserts pseudo-particles which will represent the mass distribution of the other CPUs. Initially, the mass of the pseudo-particles is set to zero, and their coordinate is set to the center of the domain-cell they correspond to. These quantities will be updated later on. Definition at line 348 of file forcetree.c. References NODE::center, DomainMoment, DomainNodeIndex, DomainNODE::mass, Nodes, and DomainNODE::s. Referenced by force_treebuild_single(). |
|
||||||||||||
|
This function allocates the memory used for storage of the tree and of auxiliary arrays needed for tree-walk and link-lists. Usually, maxnodes approximately equal to 0.7*maxpart is sufficient to store the tree for up to maxpart particles. < Maximum number of nodes in the top-level tree used for domain decomposition < Maximum number of nodes in the top-level tree used for domain decomposition Definition at line 2551 of file forcetree.c. References endrun(), MaxNodes, and Nodes_base. Referenced by init(), pmforce_periodic(), pmpotential_periodic(), and restart(). |
|
|
This function is a driver routine for constructing the gravitational oct-tree, which is done by calling a small number of other functions. Definition at line 61 of file forcetree.c. References All, force_flag_localnodes(), force_treebuild_single(), force_update_pseudoparticles(), Numnodestree, global_data_all_processes::Time, and TimeOfLastTreeConstruction. Referenced by compute_potential(), gravity_tree(), and ngb_treebuild(). |
|
|
Constructs the gravitational oct-tree. The index convention for accessing tree nodes is the following: the indices 0...NumPart-1 reference single particles, the indices All.MaxPart.... All.MaxPart+nodes-1 reference tree nodes. `Nodes_base' points to the first tree node, while `nodes' is shifted such that nodes[All.MaxPart] gives the first tree node. Finally, node indices with values 'All.MaxPart + MaxNodes' and larger indicate "pseudo particles", i.e. multipole moments of top-level nodes that lie on different CPUs. If such a node needs to be opened, the corresponding particle must be exported to that CPU. The 'Extnodes' structure parallels that of 'Nodes'. Its information is only needed for the SPH part of the computation. (The data is split onto these two structures as a tuning measure. If it is merged into 'Nodes' a somewhat bigger size of the nodes also for gravity would result, which would reduce cache utilization slightly. < Bits per dimension available for Peano-Hilbert order. Note: If peanokey is defined as type int, the allowed maximum is 10. If 64-bit integers are used, the maximum is 21 Definition at line 93 of file forcetree.c. References All, BITS_PER_DIMENSION, NODE::center, NODE::d, topnode_data::Daughter, DomainCenter, DomainCorner, DomainFac, DomainNodeIndex, dump_particles(), endrun(), force_create_empty_nodes(), force_insert_pseudo_particles(), force_update_node_recursive(), global_data_all_processes::ForceSoftening, get_random_number(), particle_data::GravCost, topnode_data::Leaf, NODE::len, MaxNodes, global_data_all_processes::MaxPart, Nextnode, NODE::nextnode, Nodes, P, peano_hilbert_key(), peanokey, particle_data::Pos, topnode_data::Size, topnode_data::StartKey, NODE::suns, ThisTask, TopNodes, particle_data::Type, and NODE::u. Referenced by force_treebuild(). |
|
||||||||||||||||
|
This routine computes the gravitational force for a given local particle, or for a particle in the communication buffer. Depending on the value of TypeOfOpeningCriterion, either the geometrical BH cell-opening criterion, or the `relative' opening criterion is used. Definition at line 1111 of file forcetree.c. References All, global_data_all_processes::ErrTolForceAcc, sph_particle_data::Hsml, particle_data::OldAcc, P, particle_data::Pos, SphP, and particle_data::Type. Referenced by gravity_tree(). |
|
||||||||||||
|
This function does the force computation with direct summation for the specified particle in the communication buffer. This can be useful for debugging purposes, in particular for explicit checks of the force accuracy. Definition at line 2632 of file forcetree.c. References All, P, particle_data::Pos, and particle_data::Type. Referenced by gravity_forcetest(). |
|
||||||||||||||||||||||||||||
|
This function computes the Ewald correction, and is needed if periodic boundary conditions together with a pure tree algorithm are used. Note that the ordinary tree walk does not carry out this correction directly as it was done in Gadget-1.1. Instead, the tree is walked a second time. This is actually faster because the "Ewald-Treewalk" can use a different opening criterion than the normal tree walk. In particular, the Ewald correction is negligible for particles that are very close, but it is large for particles that are far away (this is quite different for the normal direct force). So we can here use a different opening criterion. Sufficient accuracy is usually obtained if the node length has dropped to a certain fraction ~< 0.25 of the BoxLength. However, we may only short-cut the interaction list of the normal full Ewald tree walk if we are sure that the whole node and all daughter nodes "lie on the same side" of the periodic boundary, i.e. that the real tree walk would not find a daughter node or particle that was mapped to a different nearest neighbour position when the tree walk would be further refined. Definition at line 1734 of file forcetree.c. References gravdata_in::Acc, All, NODE::bitflags, global_data_all_processes::BoxSize, NODE::center, NODE::d, DomainTask, EN, global_data_all_processes::ErrTolTheta, Exportflag, particle_data::GravAccel, particle_data::GravCost, GravDataResult, NODE::len, particle_data::Mass, NODE::mass, global_data_all_processes::MaxPart, NEAREST, Nextnode, NODE::nextnode, gravdata_in::Ninteractions, Nodes, P, particle_data::Pos, NODE::s, NODE::sibling, NODE::u, gravdata_in::u, and gravdata_in::w. |
|
||||||||||||
|
This routine computes the gravitational potential by walking the tree. The same opening criteria is used as for the gravitational force walk. Definition at line 2000 of file forcetree.c. References All, global_data_all_processes::ErrTolForceAcc, sph_particle_data::Hsml, particle_data::OldAcc, P, particle_data::Pos, SphP, and particle_data::Type. Referenced by compute_potential(). |
|
||||||||||||
|
This function computes the short-range potential when the TreePM algorithm is used. This potential is the Newtonian potential, modified by a complementary error function. Definition at line 2247 of file forcetree.c. References All, global_data_all_processes::ErrTolForceAcc, sph_particle_data::Hsml, particle_data::OldAcc, P, particle_data::Pos, SphP, and particle_data::Type. Referenced by compute_potential(). |
|
||||||||||||
|
In the TreePM algorithm, the tree is walked only locally around the target coordinate. Tree nodes that fall outside a box of half side-length Rcut= RCUT*ASMTH*MeshSize can be discarded. The short-range potential is modified by a complementary error function, multiplied with the Newtonian form. The resulting short-range suppression compared to the Newtonian force is tabulated, because looking up from this table is faster than recomputing the corresponding factor, despite the memory-access panelty (which reduces cache performance) incurred by the table. Definition at line 1393 of file forcetree.c. References All, global_data_all_processes::ErrTolForceAcc, sph_particle_data::Hsml, particle_data::OldAcc, P, particle_data::Pos, SphP, and particle_data::Type. Referenced by gravity_tree(). |
|
|
This function frees the memory allocated for the tree, i.e. it frees the space allocated by the function force_treeallocate(). Definition at line 2615 of file forcetree.c. Referenced by pmforce_periodic(), and pmpotential_periodic(). |
|
|
This function updates the top-level tree after the multipole moments of the pseudo-particles have been updated. Definition at line 731 of file forcetree.c. References All, NODE::bitflags, NODE::d, DomainMoment, DomainNodeIndex, Extnodes, NODE::father, FLOAT, global_data_all_processes::ForceSoftening, NODE::mass, DomainNODE::mass, Nodes, NODE::s, DomainNODE::s, NODE::u, extNODE::vs, and DomainNODE::vs. |
|
|
This function updates the hmax-values in tree nodes that hold SPH particles. These values are needed to find all neighbors in the hydro-force computation. Since the Hsml-values are potentially changed in the SPH-denity computation, force_update_hmax() should be carried out just before the hydrodynamical SPH forces are computed, i.e. after density(). Definition at line 999 of file forcetree.c. References DomainHmax, DomainNodeIndex, and Extnodes. Referenced by compute_accelerations(). |
|
|
This function updates the side-length of tree nodes in case the tree is not reconstructed, but only drifted. The grouping of particles to tree nodes is not changed in this case, but some tree nodes may need to be enlarged because particles moved out of their original bounds. Definition at line 870 of file forcetree.c. References DomainNodeIndex, DomainTreeNodeLen, and Nodes. Referenced by move_particles(). |
|
|
This routine updates the hmax-values of local tree nodes. Definition at line 1036 of file forcetree.c. References NODE::d, Extnodes, Father, NODE::father, extNODE::hmax, sph_particle_data::Hsml, Nodes, SphP, and NODE::u. |
|
|
This function recursively sets the hmax-values of the top-level tree. Definition at line 1072 of file forcetree.c. References NODE::d, DomainHmax, DomainNodeIndex, Extnodes, NODE::father, extNODE::hmax, Nodes, and NODE::u. |
|
|
This function recursively enlarges nodes such that they always contain all their daughter nodes and daughter particles. Definition at line 908 of file forcetree.c. References NODE::center, NODE::d, Father, NODE::father, FLOAT, NODE::len, Nodes, P, particle_data::Pos, and NODE::u. |
|
|
This function recursively enlarges nodes of the top-level tree such that they always contain all their daughter nodes. Definition at line 955 of file forcetree.c. References NODE::center, NODE::d, DomainNodeIndex, DomainTreeNodeLen, NODE::father, FLOAT, NODE::len, Nodes, and NODE::u. |
|
||||||||||||||||
|
this routine determines the multipole moments for a given internal node and all its subnodes using a recursive computation. The result is stored in the Nodes[] structure in the sequence of this tree-walk. Note that the bitflags-variable for each node is used to store in the lowest bits some special information: Bit 0 flags whether the node belongs to the top-level tree corresponding to the domain decomposition, while Bit 1 signals whether the top-level node is dependent on local mass. If UNEQUALSOFTENINGS is set, bits 2-4 give the particle type with the maximum softening among the particles in the node, and bit 5 flags whether the node contains any particles with lower softening than that. Definition at line 425 of file forcetree.c. References All, NODE::bitflags, NODE::center, NODE::d, Extnodes, NODE::father, FLOAT, global_data_all_processes::ForceSoftening, extNODE::hmax, sph_particle_data::Hsml, NODE::mass, particle_data::Mass, global_data_all_processes::MaxPart, Nextnode, NODE::nextnode, Nodes, P, particle_data::Pos, NODE::s, NODE::sibling, SphP, NODE::suns, particle_data::Type, NODE::u, particle_data::Vel, and extNODE::vs. Referenced by force_treebuild_single(). |
|
|
This function updates the multipole moments of the pseudo-particles that represent the mass distribution on different CPUs. For that purpose, it first exchanges the necessary data, and then updates the top-level tree accordingly. The detailed implementation of these two tasks is done in separate functions. Definition at line 674 of file forcetree.c. Referenced by force_treebuild(), and move_particles(). |
1.3.9.1