diff --git a/examples/ExternalGravity/makeIC.py b/examples/ExternalGravity/makeIC.py index e82d970eba30f11ce509367fb9f1ab9dc82af50d..bba789c51a648c3295fde8db7f2a18a8e092007d 100644 --- a/examples/ExternalGravity/makeIC.py +++ b/examples/ExternalGravity/makeIC.py @@ -134,7 +134,7 @@ ds[()] = u u = numpy.zeros(1) -ids = 1 + numpy.linspace(0, numPart, numPart, endpoint=False, dtype='L') +ids = 1 + numpy.linspace(0, numPart, numPart, endpoint=False) ds = grp1.create_dataset('ParticleIDs', (numPart, ), 'L') ds[()] = ids diff --git a/examples/main.c b/examples/main.c index 38db83a334f6cbdc04afed45a12446dd82e84f1d..70a42cf2227fa2cfeae742d8ad657c5f08a642e0 100644 --- a/examples/main.c +++ b/examples/main.c @@ -300,6 +300,7 @@ int main(int argc, char *argv[]) { /* Initialise the external potential properties */ struct external_potential potential; if (with_external_gravity) initPotentialProperties(params, &us, &potential); + if (with_external_gravity && myrank == 0) printPotentialProperties(&potential); /* Read particles and space information from (GADGET) ICs */ char ICfileName[200] = ""; diff --git a/src/const.h b/src/const.h index 9e16f9d49bbee77486bc82911747c9d1ed1e0ab6..dd3e34016c1953efffca67e9c7ee1b42f455dec1 100644 --- a/src/const.h +++ b/src/const.h @@ -20,8 +20,6 @@ #ifndef SWIFT_CONST_H #define SWIFT_CONST_H -#include "physical_constants_cgs.h" - /* Hydrodynamical constants. */ #define const_hydro_gamma (5.0f / 3.0f) @@ -73,7 +71,6 @@ //#define DEFAULT_SPH /* Gravity properties */ -/* valid choices EXTERNAL_POTENTIAL_POINTMASS */ #define EXTERNAL_POTENTIAL_POINTMASS #endif /* SWIFT_CONST_H */ diff --git a/src/potentials.c b/src/potentials.c index 6cb2ab6188944f0889bdbfbb4a314b2084c70b06..f14053397e44929ead6a810754ca4a8996e1e665 100644 --- a/src/potentials.c +++ b/src/potentials.c @@ -28,12 +28,16 @@ * @brief Initialises the external potential properties in the internal system * of units. * + * @param parameter_file The parsed parameter file * @param us The current internal system of units * @param potential The external potential properties to initialize */ void initPotentialProperties(const struct swift_params* parameter_file, struct UnitSystem* us, struct external_potential* potential) { + +#ifdef EXTERNAL_POTENTIAL_POINTMASS + potential->point_mass.x = parser_get_param_double(parameter_file, "PointMass:position_x"); potential->point_mass.y = @@ -42,20 +46,20 @@ void initPotentialProperties(const struct swift_params* parameter_file, parser_get_param_double(parameter_file, "PointMass:position_z"); potential->point_mass.mass = parser_get_param_double(parameter_file, "PointMass:mass"); + +#endif /* EXTERNAL_POTENTIAL_POINTMASS */ +} + +/** + * @brief Prints the properties of the external potential to stdout. + * + * @param potential The external potential properties. + */ +void printPotentialProperties(const struct external_potential* potential) { + +#ifdef EXTERNAL_POTENTIAL_POINTMASS message("Point mass properties are (x,y,z) = (%e, %e, %e), M = %e", potential->point_mass.x, potential->point_mass.y, potential->point_mass.z, potential->point_mass.mass); - - /* potential->point_mass.x = */ - /* 50000 * PARSEC_IN_CGS / units_conversion_factor(us, UNIT_CONV_LENGTH); - */ - /* potential->point_mass.y = */ - /* 50000 * PARSEC_IN_CGS / units_conversion_factor(us, UNIT_CONV_LENGTH); - */ - /* potential->point_mass.z = */ - /* 50000 * PARSEC_IN_CGS / units_conversion_factor(us, UNIT_CONV_LENGTH); - */ - /* potential->point_mass.mass = */ - /* 1e10 * SOLAR_MASS_IN_CGS / units_conversion_factor(us, UNIT_CONV_MASS); - */ +#endif /* EXTERNAL_POTENTIAL_POINTMASS */ } diff --git a/src/potentials.h b/src/potentials.h index e489781a62890bf880aee8003df3109938278c11..3bbaab701470033943532f3c1805655d8431de4f 100644 --- a/src/potentials.h +++ b/src/potentials.h @@ -24,36 +24,34 @@ /* Config parameters. */ #include "../config.h" +/* Some standard headers. */ +#include <math.h> + /* Local includes. */ +#include "const.h" #include "error.h" -#include "physical_constants_cgs.h" +#include "part.h" #include "physical_constants.h" #include "units.h" -/* External Potential Constants */ +/* External Potential Properties */ struct external_potential { + +#ifdef EXTERNAL_POTENTIAL_POINTMASS struct { double x, y, z; double mass; } point_mass; +#endif }; /* Include exteral pointmass potential */ #ifdef EXTERNAL_POTENTIAL_POINTMASS -/* #define const_unit_length_in_cgs (1e3 * PARSEC_IN_CGS) */ -/* #define const_unit_mass_in_cgs (SOLAR_MASS_IN_CGS) */ -/* #define External_Potential_X (50000 * PARSEC_IN_CGS / - * const_unit_length_in_cgs) */ -/* #define External_Potential_Y (50000 * PARSEC_IN_CGS / - * const_unit_length_in_cgs) */ -/* #define External_Potential_Z (50000 * PARSEC_IN_CGS / - * const_unit_length_in_cgs) */ -/* #define External_Potential_Mass \ */ -/* (1e10 * SOLAR_MASS_IN_CGS / const_unit_mass_in_cgs) */ /** * @brief Computes the time-step due to the acceleration from a point mass * + * @param potential The properties of the externa potential. * @param phys_cont The physical constants in internal units. * @param gp Pointer to the g-particle data. */ @@ -63,7 +61,7 @@ __attribute__((always_inline)) const struct phys_const* const phys_const, const struct gpart* const g) { - const double G_newton = phys_const->newton_gravity; + const float G_newton = phys_const->newton_gravity; const float dx = g->x[0] - potential->point_mass.x; const float dy = g->x[1] - potential->point_mass.y; const float dz = g->x[2] - potential->point_mass.z; @@ -88,6 +86,7 @@ __attribute__((always_inline)) * @brief Computes the gravitational acceleration of a particle due to a point * mass * + * @param potential The proerties of the external potential. * @param phys_cont The physical constants in internal units. * @param gp Pointer to the g-particle data. */ @@ -95,12 +94,7 @@ __attribute__((always_inline)) INLINE static void external_gravity_pointmass( const struct external_potential* potential, const struct phys_const* const phys_const, struct gpart* g) { - /* message(" (x,y,z) = (%e, %e, %e), M= %e", potential->point_mass.x, */ - /* potential->point_mass.y, potential->point_mass.z, */ - /* potential->point_mass.mass); */ - /* exit(-1); */ - - const double G_newton = phys_const->newton_gravity; + const float G_newton = phys_const->newton_gravity; const float dx = g->x[0] - potential->point_mass.x; const float dy = g->x[1] - potential->point_mass.y; const float dz = g->x[2] - potential->point_mass.z; @@ -115,15 +109,12 @@ __attribute__((always_inline)) INLINE static void external_gravity_pointmass( } #endif /* EXTERNAL_POTENTIAL_POINTMASS */ -/** - * @brief Initialises the external potential properties in the internal system - * of units. - * - * @param us The current internal system of units - * @param potential The external potential properties to initialize - */ +/* Now, some generic functions, defined in the source file */ + void initPotentialProperties(const struct swift_params* parameter_file, struct UnitSystem* us, struct external_potential* potential); +void printPotentialProperties(const struct external_potential* potential); + #endif /* SWIFT_POTENTIALS_H */