Commit 21a4773f authored by Matthieu Schaller's avatar Matthieu Schaller
Browse files

Re-factored the whole external potential mechanism to match what is done for...

Re-factored the whole external potential mechanism to match what is done for cooling. Renamed 'disk' to 'disc'.
parent c5417916
...@@ -763,6 +763,7 @@ INPUT = @top_srcdir@ @top_srcdir@/src @top_srcdir@/tests @top_ ...@@ -763,6 +763,7 @@ INPUT = @top_srcdir@ @top_srcdir@/src @top_srcdir@/tests @top_
INPUT += @top_srcdir@/src/hydro/Minimal INPUT += @top_srcdir@/src/hydro/Minimal
INPUT += @top_srcdir@/src/gravity/Default INPUT += @top_srcdir@/src/gravity/Default
INPUT += @top_srcdir@/src/riemann INPUT += @top_srcdir@/src/riemann
INPUT += @top_srcdir@/src/potential/point_mass
INPUT += @top_srcdir@/src/cooling/const_du INPUT += @top_srcdir@/src/cooling/const_du
# This tag can be used to specify the character encoding of the source files # This tag can be used to specify the character encoding of the source files
......
...@@ -19,7 +19,7 @@ Statistics: ...@@ -19,7 +19,7 @@ Statistics:
# Parameters governing the snapshots # Parameters governing the snapshots
Snapshots: Snapshots:
basename: Disk-Patch # Common part of the name of output files basename: Disc-Patch # Common part of the name of output files
time_first: 0. # Time of the first output (in internal units) time_first: 0. # Time of the first output (in internal units)
delta_time: 8. # Time difference between consecutive outputs (in internal units) delta_time: 8. # Time difference between consecutive outputs (in internal units)
...@@ -33,11 +33,11 @@ SPH: ...@@ -33,11 +33,11 @@ SPH:
# Parameters related to the initial conditions # Parameters related to the initial conditions
InitialConditions: InitialConditions:
file_name: Disk-Patch.hdf5 # The file to read file_name: Disc-Patch.hdf5 # The file to read
# External potential parameters # External potential parameters
Disk-PatchPotential: DiscPatchPotential:
surface_density: 10. surface_density: 10.
scale_height: 100. scale_height: 100.
z_disk: 300. z_disc: 300.
timestep_mult: 0.03 timestep_mult: 0.03
...@@ -26,7 +26,7 @@ import random ...@@ -26,7 +26,7 @@ import random
# Generates N particles in a box of [0:BoxSize,0:BoxSize,-2scale_height:2scale_height] # Generates N particles in a box of [0:BoxSize,0:BoxSize,-2scale_height:2scale_height]
# see Creasey, Theuns & Bower, 2013, for the equations: # see Creasey, Theuns & Bower, 2013, for the equations:
# disk parameters are: surface density sigma # disc parameters are: surface density sigma
# scale height b # scale height b
# density: rho(z) = (sigma/2b) sech^2(z/b) # density: rho(z) = (sigma/2b) sech^2(z/b)
# isothermal velocity dispersion = <v_z^2? = b pi G sigma # isothermal velocity dispersion = <v_z^2? = b pi G sigma
...@@ -79,7 +79,7 @@ N = int(sys.argv[1]) # Number of particles ...@@ -79,7 +79,7 @@ N = int(sys.argv[1]) # Number of particles
rho = 2. # Density rho = 2. # Density
P = 1. # Pressure P = 1. # Pressure
gamma = 5./3. # Gas adiabatic index gamma = 5./3. # Gas adiabatic index
fileName = "Disk-Patch.hdf5" fileName = "Disc-Patch.hdf5"
#--------------------------------------------------- #---------------------------------------------------
......
#!/bin/bash #!/bin/bash
# Generate the initial conditions if they are not present. # Generate the initial conditions if they are not present.
if [ ! -e Isothermal.hdf5 ] if [ ! -e Disc-Patch.hdf5 ]
then then
echo "Generating initial conditions for the disk-patch example..." echo "Generating initial conditions for the disc-patch example..."
python makeIC.py 1000 python makeIC.py 1000
fi fi
../../swift -g -t 2 disk-patch.yml ../../swift -g -t 2 disc-patch.yml
...@@ -8,7 +8,7 @@ iplot = 1 ; if iplot = 1, make plot of E/Lz conservation, else, simply compare f ...@@ -8,7 +8,7 @@ iplot = 1 ; if iplot = 1, make plot of E/Lz conservation, else, simply compare f
@physunits @physunits
indir = './' indir = './'
basefile = 'Disk-Patch_' basefile = 'Disc-Patch_'
; set properties of potential ; set properties of potential
uL = phys.pc ; unit of length uL = phys.pc ; unit of length
......
Generates and evolves a disk-patch, where gas is in hydrostatic Generates and evolves a disc-patch, where gas is in hydrostatic
equilibrium with an imposed external gravitational force, using the equilibrium with an imposed external gravitational force, using the
equations from Creasey, Theuns & Bower, 2013, MNRAS, Volume 429, equations from Creasey, Theuns & Bower, 2013, MNRAS, Volume 429,
Issue 3, p.1922-1948. Issue 3, p.1922-1948.
...@@ -10,11 +10,11 @@ To generate ICs ready for a scientific run: ...@@ -10,11 +10,11 @@ To generate ICs ready for a scientific run:
2) Generate pre-ICs by running the 'makeIC.py' script. 2) Generate pre-ICs by running the 'makeIC.py' script.
3) Run SWIFT with an isothermal EoS, no cooling nor feedback, and the 3) Run SWIFT with an isothermal EoS, no cooling nor feedback, and the
disk-patch potential switched on and using the parameters from disc-patch potential switched on and using the parameters from
'disk-patch-icc.yml' 'disc-patch-icc.yml'
4) The ICs are then ready to be run for a science problem. Rename the last 4) The ICs are then ready to be run for a science problem. Rename the last
output to 'Disk-Patch-dynamic.hdf5'. These are now the ICs for the actual test. output to 'Disc-Patch-dynamic.hdf5'. These are now the ICs for the actual test.
When running SWIFT with the parameters from 'disk-patch.yml' and an When running SWIFT with the parameters from 'disc-patch.yml' and an
ideal gas EoS on these ICs the disk should stay in equilibrium. ideal gas EoS on these ICs the disc should stay in equilibrium.
...@@ -19,7 +19,7 @@ Statistics: ...@@ -19,7 +19,7 @@ Statistics:
# Parameters governing the snapshots # Parameters governing the snapshots
Snapshots: Snapshots:
basename: Disk-Patch # Common part of the name of output files basename: Disc-Patch # Common part of the name of output files
time_first: 0. # Time of the first output (in internal units) time_first: 0. # Time of the first output (in internal units)
delta_time: 12. # Time difference between consecutive outputs (in internal units) delta_time: 12. # Time difference between consecutive outputs (in internal units)
...@@ -33,12 +33,12 @@ SPH: ...@@ -33,12 +33,12 @@ SPH:
# Parameters related to the initial conditions # Parameters related to the initial conditions
InitialConditions: InitialConditions:
file_name: Disk-Patch.hdf5 # The file to read file_name: Disc-Patch.hdf5 # The file to read
# External potential parameters # External potential parameters
Disk-PatchPotential: DiscPatchPotential:
surface_density: 10. surface_density: 10.
scale_height: 100. scale_height: 100.
z_disk: 200. z_disc: 200.
timestep_mult: 0.03 timestep_mult: 0.03
growth_time: 5. growth_time: 5.
...@@ -19,7 +19,7 @@ Statistics: ...@@ -19,7 +19,7 @@ Statistics:
# Parameters governing the snapshots # Parameters governing the snapshots
Snapshots: Snapshots:
basename: Disk-Patch-dynamic # Common part of the name of output files basename: Disc-Patch-dynamic # Common part of the name of output files
time_first: 968. # Time of the first output (in internal units) time_first: 968. # Time of the first output (in internal units)
delta_time: 24. # Time difference between consecutive outputs (in internal units) delta_time: 24. # Time difference between consecutive outputs (in internal units)
...@@ -33,11 +33,11 @@ SPH: ...@@ -33,11 +33,11 @@ SPH:
# Parameters related to the initial conditions # Parameters related to the initial conditions
InitialConditions: InitialConditions:
file_name: Disk-Patch-dynamic.hdf5 # The file to read file_name: Disc-Patch-dynamic.hdf5 # The file to read
# External potential parameters # External potential parameters
Disk-PatchPotential: DiscPatchPotential:
surface_density: 10. surface_density: 10.
scale_height: 100. scale_height: 100.
z_disk: 200. z_disc: 200.
timestep_mult: 0.03 timestep_mult: 0.03
...@@ -25,9 +25,9 @@ import math ...@@ -25,9 +25,9 @@ import math
import random import random
import matplotlib.pyplot as plt import matplotlib.pyplot as plt
# Generates a disk-patch in hydrostatic equilibrium # Generates a disc-patch in hydrostatic equilibrium
# see Creasey, Theuns & Bower, 2013, for the equations: # see Creasey, Theuns & Bower, 2013, for the equations:
# disk parameters are: surface density sigma # disc parameters are: surface density sigma
# scale height b # scale height b
# density: rho(z) = (sigma/2b) sech^2(z/b) # density: rho(z) = (sigma/2b) sech^2(z/b)
# isothermal velocity dispersion = <v_z^2? = b pi G sigma # isothermal velocity dispersion = <v_z^2? = b pi G sigma
...@@ -79,7 +79,7 @@ Radius = 100. # maximum radius of particles [kpc] ...@@ -79,7 +79,7 @@ Radius = 100. # maximum radius of particles [kpc]
G = const_G G = const_G
# File # File
fileName = "Disk-Patch.hdf5" fileName = "Disc-Patch.hdf5"
#--------------------------------------------------- #---------------------------------------------------
mass = 1 mass = 1
...@@ -145,7 +145,7 @@ mass = 0.*h + pmass ...@@ -145,7 +145,7 @@ mass = 0.*h + pmass
entropy_flag = 0 entropy_flag = 0
vel = 0 + 0 * pos vel = 0 + 0 * pos
# move centre of disk to middle of box # move centre of disc to middle of box
pos[:,:] += boxSize/2 pos[:,:] += boxSize/2
......
...@@ -8,7 +8,7 @@ iplot = 1 ; if iplot = 1, make plot of E/Lz conservation, else, simply compare f ...@@ -8,7 +8,7 @@ iplot = 1 ; if iplot = 1, make plot of E/Lz conservation, else, simply compare f
@physunits @physunits
indir = './' indir = './'
basefile = 'Disk-Patch_' basefile = 'Disc-Patch_'
; set properties of potential ; set properties of potential
uL = phys.pc ; unit of length uL = phys.pc ; unit of length
......
...@@ -38,7 +38,7 @@ InitialConditions: ...@@ -38,7 +38,7 @@ InitialConditions:
shift_z: 50. shift_z: 50.
# External potential parameters # External potential parameters
PointMass: PointMassPotential:
position_x: 50. # location of external point mass in internal units position_x: 50. # location of external point mass in internal units
position_y: 50. position_y: 50.
position_z: 50. position_z: 50.
......
...@@ -66,7 +66,7 @@ DomainDecomposition: ...@@ -66,7 +66,7 @@ DomainDecomposition:
# Parameters related to external potentials -------------------------------------------- # Parameters related to external potentials --------------------------------------------
# Point mass external potentials # Point mass external potentials
PointMass: PointMassPotential:
position_x: 50. # location of external point mass (internal units) position_x: 50. # location of external point mass (internal units)
position_y: 50. position_y: 50.
position_z: 50. position_z: 50.
...@@ -82,12 +82,12 @@ IsothermalPotential: ...@@ -82,12 +82,12 @@ IsothermalPotential:
timestep_mult: 0.03 # Dimensionless pre-factor for the time-step condition timestep_mult: 0.03 # Dimensionless pre-factor for the time-step condition
# Disk-patch potential parameters # Disk-patch potential parameters
Disk-PatchPotential: DiscPatchPotential:
surface_density: 10. # Surface density of the disk (internal units) surface_density: 10. # Surface density of the disc (internal units)
scale_height: 100. # Scale height of the disk (internal units) scale_height: 100. # Scale height of the disc (internal units)
z_disk: 200. # Disk height (internal units) z_disc: 200. # Disc height (internal units)
timestep_mult: 0.03 # Dimensionless pre-factor for the time-step condition timestep_mult: 0.03 # Dimensionless pre-factor for the time-step condition
growth_time: 5. # (Optional) Time for the disk to grow to its final size (multiple of the dynamical time) growth_time: 5. # (Optional) Time for the disc to grow to its final size (multiple of the dynamical time)
# Parameters related to cooling function ---------------------------------------------- # Parameters related to cooling function ----------------------------------------------
......
...@@ -42,7 +42,7 @@ endif ...@@ -42,7 +42,7 @@ endif
include_HEADERS = space.h runner.h queue.h task.h lock.h cell.h part.h const.h \ include_HEADERS = space.h runner.h queue.h task.h lock.h cell.h part.h const.h \
engine.h swift.h serial_io.h timers.h debug.h scheduler.h proxy.h parallel_io.h \ engine.h swift.h serial_io.h timers.h debug.h scheduler.h proxy.h parallel_io.h \
common_io.h single_io.h multipole.h map.h tools.h partition.h clocks.h parser.h \ common_io.h single_io.h multipole.h map.h tools.h partition.h clocks.h parser.h \
physical_constants.h physical_constants_cgs.h potentials.h version.h \ physical_constants.h physical_constants_cgs.h potential.h version.h \
hydro_properties.h threadpool.h cooling.h cooling_struct.h hydro_properties.h threadpool.h cooling.h cooling_struct.h
...@@ -51,7 +51,7 @@ AM_SOURCES = space.c runner.c queue.c task.c cell.c engine.c \ ...@@ -51,7 +51,7 @@ AM_SOURCES = space.c runner.c queue.c task.c cell.c engine.c \
serial_io.c timers.c debug.c scheduler.c proxy.c parallel_io.c \ serial_io.c timers.c debug.c scheduler.c proxy.c parallel_io.c \
units.c common_io.c single_io.c multipole.c version.c map.c \ units.c common_io.c single_io.c multipole.c version.c map.c \
kernel_hydro.c tools.c part.c partition.c clocks.c parser.c \ kernel_hydro.c tools.c part.c partition.c clocks.c parser.c \
physical_constants.c potentials.c hydro_properties.c \ physical_constants.c potential.c hydro_properties.c \
runner_doiact_fft.c threadpool.c cooling.c runner_doiact_fft.c threadpool.c cooling.c
# Include files for distribution, not installation. # Include files for distribution, not installation.
......
...@@ -90,9 +90,10 @@ ...@@ -90,9 +90,10 @@
#define const_gravity_eta 0.025f #define const_gravity_eta 0.025f
/* External gravity properties */ /* External gravity properties */
#define EXTERNAL_POTENTIAL_POINTMASS #define EXTERNAL_POTENTIAL_NONE
//#define EXTERNAL_POTENTIAL_POINTMASS
//#define EXTERNAL_POTENTIAL_ISOTHERMALPOTENTIAL //#define EXTERNAL_POTENTIAL_ISOTHERMALPOTENTIAL
//#define EXTERNAL_POTENTIAL_DISK_PATCH //#define EXTERNAL_POTENTIAL_DISC_PATCH
/* Cooling properties */ /* Cooling properties */
#define COOLING_NONE #define COOLING_NONE
......
...@@ -41,7 +41,7 @@ ...@@ -41,7 +41,7 @@
#include "cooling_struct.h" #include "cooling_struct.h"
#include "parser.h" #include "parser.h"
#include "partition.h" #include "partition.h"
#include "potentials.h" #include "potential.h"
#include "runner.h" #include "runner.h"
#include "scheduler.h" #include "scheduler.h"
#include "space.h" #include "space.h"
......
...@@ -22,37 +22,7 @@ ...@@ -22,37 +22,7 @@
#include <float.h> #include <float.h>
#include "minmax.h" #include "minmax.h"
#include "potentials.h" #include "physical_constants.h"
/**
* @brief Computes the gravity time-step of a given particle due to an external
*potential.
*
* This function only branches towards the potential chosen by the user.
*
* @param potential The properties of the external potential.
* @param phys_const The physical constants in internal units.
* @param gp Pointer to the g-particle data.
*/
__attribute__((always_inline)) INLINE static float
gravity_compute_timestep_external(const struct external_potential* potential,
const struct phys_const* const phys_const,
const struct gpart* const gp) {
float dt = FLT_MAX;
#ifdef EXTERNAL_POTENTIAL_POINTMASS
dt = min(dt, external_gravity_pointmass_timestep(potential, phys_const, gp));
#endif
#ifdef EXTERNAL_POTENTIAL_ISOTHERMALPOTENTIAL
dt = min(dt, external_gravity_isothermalpotential_timestep(potential,
phys_const, gp));
#endif
#ifdef EXTERNAL_POTENTIAL_DISK_PATCH
dt = min(dt, external_gravity_disk_patch_timestep(potential, phys_const, gp));
#endif
return dt;
}
/** /**
* @brief Computes the gravity time-step of a given particle due to self-gravity * @brief Computes the gravity time-step of a given particle due to self-gravity
...@@ -125,31 +95,6 @@ __attribute__((always_inline)) INLINE static void gravity_end_force( ...@@ -125,31 +95,6 @@ __attribute__((always_inline)) INLINE static void gravity_end_force(
gp->a_grav[2] *= const_G; gp->a_grav[2] *= const_G;
} }
/**
* @brief Computes the gravitational acceleration induced by external potentials
*
* This function only branches towards the potential chosen by the user.
*
* @param time The current time in internal units.
* @param potential The properties of the external potential.
* @param phys_const The physical constants in internal units.
* @param gp The particle to act upon.
*/
__attribute__((always_inline)) INLINE static void external_gravity(
double time, const struct external_potential* potential,
const struct phys_const* const phys_const, struct gpart* gp) {
#ifdef EXTERNAL_POTENTIAL_POINTMASS
external_gravity_pointmass(potential, phys_const, gp);
#endif
#ifdef EXTERNAL_POTENTIAL_ISOTHERMALPOTENTIAL
external_gravity_isothermalpotential(potential, phys_const, gp);
#endif
#ifdef EXTERNAL_POTENTIAL_DISK_PATCH
external_gravity_disk_patch_potential(time, potential, phys_const, gp);
#endif
}
/** /**
* @brief Kick the additional variables * @brief Kick the additional variables
* *
......
/*******************************************************************************
* This file is part of SWIFT.
* Copyright (c) 2016 Tom Theuns (tom.theuns@durham.ac.uk)
* Matthieu Schaller (matthieu.schaller@durham.ac.uk)
*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU Lesser General Public License as published
* by the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* 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 License for more details.
*
* You should have received a copy of the GNU Lesser General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*
******************************************************************************/
/* Config parameters. */
#include "../config.h"
/* This object's header. */
#include "potential.h"
/**
* @brief Initialises the external potential properties in the internal system
* of units.
*
* @param parameter_file The parsed parameter file
* @param phys_const Physical constants in internal units
* @param us The current internal system of units
* @param potential The external potential properties to initialize
*/
void potential_init(const struct swift_params* parameter_file,
const struct phys_const* phys_const,
const struct UnitSystem* us,
struct external_potential* potential) {
potential_init_backend(parameter_file, phys_const, us, potential);
}
/**
* @brief Prints the properties of the external potential to stdout.
*
* @param potential The external potential properties.
*/
void potential_print(const struct external_potential* potential) {
potential_print_backend(potential);
}
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment