From 557ab4d258d13bddc21cee78f4097974e79e5833 Mon Sep 17 00:00:00 2001
From: Jacob Kegerreis <jacob.kegerreis@durham.ac.uk>
Date: Thu, 19 Apr 2018 13:06:55 +0100
Subject: [PATCH] Add (placeholder for now) Tillotson equation of state
---
README_planetary | 15 +-
configure.ac | 5 +-
src/equation_of_state.h | 2 +
.../tillotson_iron/equation_of_state.h | 203 ++++++++++++++++++
4 files changed, 221 insertions(+), 4 deletions(-)
create mode 100644 src/equation_of_state/tillotson_iron/equation_of_state.h
diff --git a/README_planetary b/README_planetary
index 4e352ea388..ca809a7e5a 100644
--- a/README_planetary
+++ b/README_planetary
@@ -1,8 +1,17 @@
-SWIFT branch to add equations of state for planetary materials such as iron,
-silicate rocks, ices, and atmospheres; as well as allowing different particles
-to have different materials for the purpose of e.g. giant impact simulations.
+SWIFT 'planetary' branch
Jacob Kegerreis
+
jacob.kegerreis@durham.ac.uk
+# Aims:
+
++ Implement various equations of state, e.g. Tillotson, ANEOS
+
++ Allow multiple SPH materials with different EOS for different particles
+
++ Add examples, e.g. canonical Moon-forming giant impact
+
+# Status:
+
diff --git a/configure.ac b/configure.ac
index e730619ea1..bb7c0eb10c 100644
--- a/configure.ac
+++ b/configure.ac
@@ -890,7 +890,7 @@ esac
# Equation of state
AC_ARG_WITH([equation-of-state],
[AS_HELP_STRING([--with-equation-of-state=<EoS>],
- [equation of state @<:@ideal-gas, isothermal-gas default: ideal-gas@:>@]
+ [equation of state @<:@ideal-gas, isothermal-gas, tillotson-iron default: ideal-gas@:>@]
)],
[with_eos="$withval"],
[with_eos="ideal-gas"]
@@ -902,6 +902,9 @@ case "$with_eos" in
isothermal-gas)
AC_DEFINE([EOS_ISOTHERMAL_GAS], [1], [Isothermal gas equation of state])
;;
+ tillotson-iron)
+ AC_DEFINE([EOS_TILLOTSON_IRON], [1], [Tillotson iron equation of state])
+ ;;
*)
AC_MSG_ERROR([Unknown equation of state: $with_eos])
;;
diff --git a/src/equation_of_state.h b/src/equation_of_state.h
index 195b52514f..1f7bf5d77f 100644
--- a/src/equation_of_state.h
+++ b/src/equation_of_state.h
@@ -34,6 +34,8 @@
#include "./equation_of_state/ideal_gas/equation_of_state.h"
#elif defined(EOS_ISOTHERMAL_GAS)
#include "./equation_of_state/isothermal/equation_of_state.h"
+#elif defined(EOS_TILLOTSON_IRON)
+#include "./equation_of_state/tillotson_iron/equation_of_state.h"
#else
#error "Invalid choice of equation of state"
#endif
diff --git a/src/equation_of_state/tillotson_iron/equation_of_state.h b/src/equation_of_state/tillotson_iron/equation_of_state.h
new file mode 100644
index 0000000000..54a1e15b66
--- /dev/null
+++ b/src/equation_of_state/tillotson_iron/equation_of_state.h
@@ -0,0 +1,203 @@
+/*******************************************************************************
+ * This file is part of SWIFT.
+ * Copyright (c) 2016 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/>.
+ *
+ ******************************************************************************/
+#ifndef SWIFT_TILLOTSON_IRON_EQUATION_OF_STATE_H
+#define SWIFT_TILLOTSON_IRON_EQUATION_OF_STATE_H
+
+/* Some standard headers. */
+#include <math.h>
+
+/* Local headers. */
+#include "adiabatic_index.h"
+#include "common_io.h"
+#include "debug.h"
+#include "inline.h"
+
+extern struct eos_parameters eos;
+/* ------------------------------------------------------------------------- */
+
+struct eos_parameters {};
+
+/**
+ * @brief Returns the internal energy given density and entropy
+ *
+ * Computes \f$u = \frac{S\rho^{\gamma-1} }{\gamma - 1}\f$.
+ *
+ * @param density The density \f$\rho\f$.
+ * @param entropy The entropy \f$S\f$.
+ */
+__attribute__((always_inline)) INLINE static float
+gas_internal_energy_from_entropy(float density, float entropy) {
+
+ return entropy * pow_gamma_minus_one(density) *
+ hydro_one_over_gamma_minus_one;
+}
+
+/**
+ * @brief Returns the pressure given density and entropy
+ *
+ * Computes \f$P = S\rho^\gamma\f$.
+ *
+ * @param density The density \f$\rho\f$.
+ * @param entropy The entropy \f$S\f$.
+ */
+__attribute__((always_inline)) INLINE static float gas_pressure_from_entropy(
+ float density, float entropy) {
+
+ return entropy * pow_gamma(density);
+}
+
+/**
+ * @brief Returns the entropy given density and pressure.
+ *
+ * Computes \f$A = \frac{P}{\rho^\gamma}\f$.
+ *
+ * @param density The density \f$\rho\f$.
+ * @param pressure The pressure \f$P\f$.
+ * @return The entropy \f$A\f$.
+ */
+__attribute__((always_inline)) INLINE static float gas_entropy_from_pressure(
+ float density, float pressure) {
+
+ return pressure * pow_minus_gamma(density);
+}
+
+/**
+ * @brief Returns the sound speed given density and entropy
+ *
+ * Computes \f$c = \sqrt{\gamma S \rho^{\gamma-1}}\f$.
+ *
+ * @param density The density \f$\rho\f$.
+ * @param entropy The entropy \f$S\f$.
+ */
+__attribute__((always_inline)) INLINE static float gas_soundspeed_from_entropy(
+ float density, float entropy) {
+
+ return sqrtf(hydro_gamma * pow_gamma_minus_one(density) * entropy);
+}
+
+/**
+ * @brief Returns the entropy given density and internal energy
+ *
+ * Computes \f$S = \frac{(\gamma - 1)u}{\rho^{\gamma-1}}\f$.
+ *
+ * @param density The density \f$\rho\f$
+ * @param u The internal energy \f$u\f$
+ */
+__attribute__((always_inline)) INLINE static float
+gas_entropy_from_internal_energy(float density, float u) {
+
+ return hydro_gamma_minus_one * u * pow_minus_gamma_minus_one(density);
+}
+
+/**
+ * @brief Returns the pressure given density and internal energy
+ *
+ * Computes \f$P = (\gamma - 1)u\rho\f$.
+ *
+ * @param density The density \f$\rho\f$
+ * @param u The internal energy \f$u\f$
+ */
+__attribute__((always_inline)) INLINE static float
+gas_pressure_from_internal_energy(float density, float u) {
+
+ return hydro_gamma_minus_one * u * density;
+}
+
+/**
+ * @brief Returns the internal energy given density and pressure.
+ *
+ * Computes \f$u = \frac{1}{\gamma - 1}\frac{P}{\rho}\f$.
+ *
+ * @param density The density \f$\rho\f$.
+ * @param pressure The pressure \f$P\f$.
+ * @return The internal energy \f$u\f$.
+ */
+__attribute__((always_inline)) INLINE static float
+gas_internal_energy_from_pressure(float density, float pressure) {
+ return hydro_one_over_gamma_minus_one * pressure / density;
+}
+
+/**
+ * @brief Returns the sound speed given density and internal energy
+ *
+ * Computes \f$c = \sqrt{\gamma (\gamma - 1) u }\f$.
+ *
+ * @param density The density \f$\rho\f$
+ * @param u The internal energy \f$u\f$
+ */
+__attribute__((always_inline)) INLINE static float
+gas_soundspeed_from_internal_energy(float density, float u) {
+
+ return sqrtf(u * hydro_gamma * hydro_gamma_minus_one);
+}
+
+/**
+ * @brief Returns the sound speed given density and pressure
+ *
+ * Computes \f$c = \sqrt{\frac{\gamma P}{\rho} }\f$.
+ *
+ * @param density The density \f$\rho\f$
+ * @param P The pressure \f$P\f$
+ */
+__attribute__((always_inline)) INLINE static float gas_soundspeed_from_pressure(
+ float density, float P) {
+
+ const float density_inv = 1.f / density;
+ return sqrtf(hydro_gamma * P * density_inv);
+}
+
+/**
+ * @brief Initialize the eos parameters
+ *
+ * @param e The #eos_paramters
+ * @param params The parsed parameters
+ */
+__attribute__((always_inline)) INLINE static void eos_init(
+ struct eos_parameters *e, const struct swift_params *params) {}
+
+/**
+ * @brief Print the equation of state
+ *
+ * @param e The #eos_parameters
+ */
+__attribute__((always_inline)) INLINE static void eos_print(
+ const struct eos_parameters *e) {
+
+ message("Equation of state: Ideal gas.");
+
+ message("Adiabatic index gamma: %f.", hydro_gamma);
+}
+
+#if defined(HAVE_HDF5)
+/**
+ * @brief Write equation of state information to the snapshot
+ *
+ * @param h_grpsph The HDF5 group in which to write
+ * @param e The #eos_parameters
+ */
+__attribute__((always_inline)) INLINE static void eos_print_snapshot(
+ hid_t h_grpsph, const struct eos_parameters *e) {
+
+ io_write_attribute_f(h_grpsph, "Adiabatic index", hydro_gamma);
+
+ io_write_attribute_s(h_grpsph, "Equation of state", "Ideal gas");
+}
+#endif
+
+#endif /* SWIFT_TILLOTSON_IRON_EQUATION_OF_STATE_H */
--
GitLab