From 8d3c309dfbe5902704408d0015a04a3cf718782e Mon Sep 17 00:00:00 2001
From: Robert Izzard <r.izzard@surrey.ac.uk>
Date: Sun, 1 Aug 2021 18:07:59 +0100
Subject: [PATCH] fix issue where the timestep wasn't being limited in a system
with only massless remnants when fixed_timesteps should be triggered
fixed issue where the code got in an infinite loop after the above failure
fix TPAGB issue where massive cores were checked vs MCh just before the TPAGB: we now explicitly check that the star is a TPAGB star to do this
---
meson.build | 4 ++--
src/binary_c_debug.h | 4 ++--
src/binary_c_error_codes.h | 8 ++++++++
src/binary_c_parameters.h | 2 +-
src/evolution/set_next_timestep.c | 11 ++++++++++
src/exit_binaryc.c | 9 ++++-----
src/logging/logfile_string.c | 2 +-
src/logging/open_log_files.c | 6 +++---
src/memory/free_memory.c | 6 ++++--
src/memory/free_persistent_data.c | 12 ++++++++++-
src/stellar_structure/stellar_structure_AGB.c | 9 ++++++---
.../stellar_timescales_high_mass_GB.c | 2 +-
src/timestep/stellar_timestep.c | 8 +++-----
src/timestep/timestep_fixed_timesteps.c | 16 ++++++++++-----
src/timestep/timestep_fixed_trigger.c | 12 +++++------
src/timestep/timestep_hard_limits.c | 20 +++----------------
src/timestep/timestep_limits.c | 3 +--
src/timestep/timestep_prototypes.h | 7 +++++--
18 files changed, 83 insertions(+), 58 deletions(-)
diff --git a/meson.build b/meson.build
index 5ddbb9173..9ae9278ec 100644
--- a/meson.build
+++ b/meson.build
@@ -455,9 +455,9 @@ else
# Note: we disable some of them
'-fno-associative-math',
'-fno-math-errno',
- #'-ffinite-math-only', # not good
+ '-ffinite-math-only', # not good
'-fno-rounding-math',
- #'-fno-signaling-nans', # should be off : nan-signalling is required
+ '-fno-signaling-nans', # should be off : nan-signalling is required
'-fcx-limited-range',
'-fexcess-precision=fast',
# never allow unsafe math optimizations
diff --git a/src/binary_c_debug.h b/src/binary_c_debug.h
index b74b8954d..b11e707f1 100644
--- a/src/binary_c_debug.h
+++ b/src/binary_c_debug.h
@@ -78,7 +78,7 @@
* the same, but stardata is defined everywhere as a global
* variable.
*/
-#define Debug_expression (TRUE)
+#define Debug_expression ()
/*
* If you define Debug_stop_expression, and it is at any time TRUE,
@@ -106,7 +106,7 @@
* this is not shown.
*/
-#define Debug_show_expression " %d %d ",stardata->model.model_number,stardata->star[0].stellar_type
+#define Debug_show_expression " "
//#undef Debug_show_expression
/*
diff --git a/src/binary_c_error_codes.h b/src/binary_c_error_codes.h
index d7c00f526..d33e78735 100644
--- a/src/binary_c_error_codes.h
+++ b/src/binary_c_error_codes.h
@@ -41,4 +41,12 @@ static char * binary_c_error_code_strings[] MAYBE_UNUSED = { BINARY_C_ERROR_CODE
#undef X
+#define Binary_c_error_code_is_bad(CODE) ( \
+ (CODE)!=BINARY_C_NORMAL_EXIT && \
+ (CODE)!=BINARY_C_QUIET_EXIT && \
+ (CODE)!=BINARY_C_SPECIAL_EXIT && \
+ (CODE)!=BINARY_C_NORMAL_BATCHMODE_EXIT \
+ )
+
+
#endif // BINARY_C_ERROR_CODES_H
diff --git a/src/binary_c_parameters.h b/src/binary_c_parameters.h
index dcf45b9e9..04660e8e8 100644
--- a/src/binary_c_parameters.h
+++ b/src/binary_c_parameters.h
@@ -64,7 +64,7 @@
*
* Post-MS evolution has not yet been implemented.
*/
-#define MINT
+//#define MINT
/************************************************************
diff --git a/src/evolution/set_next_timestep.c b/src/evolution/set_next_timestep.c
index e1ca855c0..3d28b4403 100644
--- a/src/evolution/set_next_timestep.c
+++ b/src/evolution/set_next_timestep.c
@@ -37,6 +37,17 @@ void set_next_timestep(struct stardata_t * RESTRICT const stardata)
stardata->model.dtm);
}
+ /*
+ * Apply nonstellar timestep
+ */
+ double dtm = time_remaining;
+ nonstellar_timestep(stardata,&dtm);
+ stardata->model.dtm =
+ Is_zero(stardata->model.dtm) ?
+ dtm :
+ Min(stardata->model.dtm,
+ dtm);
+ Dprint("nonstellar timestep %g\n",stardata->model.dtm);
/*
* Save the timestep on the first
diff --git a/src/exit_binaryc.c b/src/exit_binaryc.c
index 307029470..5a3a3fe24 100644
--- a/src/exit_binaryc.c
+++ b/src/exit_binaryc.c
@@ -4,7 +4,6 @@
*
* Note the unusual header file list: please do not change it! :)
*/
-
#include <stdlib.h>
#include <stdarg.h>
@@ -94,6 +93,7 @@ void Gnu_format_args(6,7)
#endif //BATCHMODE
}
+ Dprint("calling exit_binary_c_handler\n");
va_list args,dummy MAYBE_UNUSED;
va_start(args,format);
va_start(dummy,format);
@@ -264,10 +264,7 @@ static void No_return Gnu_format_args(8,0)
* Non-batchmode or batchmode is off : actually exit
*/
Dprint("binary_c : we really want to exit (not a batchmode long jump)\n");
- if((binary_c_error_code!=BINARY_C_NORMAL_EXIT &&
- binary_c_error_code!=BINARY_C_QUIET_EXIT &&
- binary_c_error_code!=BINARY_C_SPECIAL_EXIT &&
- binary_c_error_code!=BINARY_C_NORMAL_BATCHMODE_EXIT) || (DEBUG))
+ if(Binary_c_error_code_is_bad(binary_c_error_code) || (DEBUG))
{
/*
* If stdout is not a terminal, e.g. it's a file,
@@ -321,7 +318,9 @@ static void No_return Gnu_format_args(8,0)
/*
* Free all allocated memory, including the preferences, store and tmpstore
*/
+ printf("try to free memory\n");fflush(NULL);
if(stardata!=NULL) free_memory(&stardata,TRUE,TRUE,TRUE,TRUE,TRUE);
+ printf("now try to exit\n");fflush(NULL);
#undef exit
diff --git a/src/logging/logfile_string.c b/src/logging/logfile_string.c
index 45d0c3777..aba948b6b 100644
--- a/src/logging/logfile_string.c
+++ b/src/logging/logfile_string.c
@@ -128,7 +128,7 @@ char * logfile_string(const double time,
#else
snprintf(logstring,
LOGSTRINGLENGTH,
- "%s%11.4lf%s%9.3lf%s%s%9.3lf%s%s%3s%s%s%3s%s%s%13.5g%s%s%8.3g%c%s%s%6.2lf%s%s%s%s%8.3lf%s%s%s%8.3lf%s %s%s%s",
+ "%s%11.4lf%s%9.3lf%s%s%9.3lf%s%s%3s%s%s%3s%s%s%13.5g%s%s%9.3g%c%s%s%6.2lf%s%s%s%s%8.3lf%s%s%s%8.3lf%s %s%s%s",
File_log_colour(stardata->store->colours[ORANGE]),
time,
File_log_colour(stardata->store->colours[RED]),
diff --git a/src/logging/open_log_files.c b/src/logging/open_log_files.c
index 3ad788844..896e4ddb5 100644
--- a/src/logging/open_log_files.c
+++ b/src/logging/open_log_files.c
@@ -45,7 +45,7 @@ void open_log_files(FILE *(*log_fp),
{
/* redirect output to stdout */
*log_fp=stdout;
- fprintf(*log_fp,"%s TIME M1 M2 K1 K2 SEP%s PER%s ECC%s R1/ROL1 R2/ROL2 TYPE",
+ fprintf(*log_fp,"%s TIME M1 M2 K1 K2 SEP%s PER%s ECC%s R1/ROL1 R2/ROL2 TYPE",
File_log_colour(stardata->store->colours[WHITE]),
logspace,
logspace,
@@ -64,7 +64,7 @@ void open_log_files(FILE *(*log_fp),
{
/* redirect output to stderr */
*log_fp=stderr;
- fprintf(*log_fp," TIME M1%s M2%s K1%s K2%s SEP%s PER%s ECC%s R1/ROL1%s R2/ROL2%s TYPE%s",
+ fprintf(*log_fp," TIME M1%s M2%s K1%s K2%s SEP%s PER%s ECC%s R1/ROL1%s R2/ROL2%s TYPE%s",
File_log_colour(stardata->store->colours[WHITE]),
logspace,
logspace,
@@ -105,7 +105,7 @@ void open_log_files(FILE *(*log_fp),
else
{
fprintf(*log_fp,
- "%s TIME M1%s M2%s K1%s K2%s SEP%s PER%s ECC%s R1/ROL1%s R2/ROL2%s TYPE RANDOM_SEED=%ld RANDOM_COUNT=%ld%s",
+ "%s TIME M1%s M2%s K1%s K2%s SEP%s PER%s ECC%s R1/ROL1%s R2/ROL2%s TYPE RANDOM_SEED=%ld RANDOM_COUNT=%ld%s",
File_log_colour(stardata->store->colours[WHITE]),
logspace,
logspace,
diff --git a/src/memory/free_memory.c b/src/memory/free_memory.c
index 1610becb4..9fcc16d4d 100644
--- a/src/memory/free_memory.c
+++ b/src/memory/free_memory.c
@@ -46,12 +46,14 @@ void free_memory(struct stardata_t ** RESTRICT const sp,
* Free persistent data including calls to routines
* to output the data if required.
*/
- if(free_persistent == TRUE)
+ if(stardata != NULL && free_persistent == TRUE)
{
FMDebug("persistent_data %p\n",(void*)stardata->persistent_data);
- free_persistent_data(*sp);
+ free_persistent_data(stardata);
}
+ FMDebug("done that\n");
+
#ifdef MEMMAP
FMDebug("Call munmap\n");
munmap(stardata->biglog,BIGLOG_SIZE);
diff --git a/src/memory/free_persistent_data.c b/src/memory/free_persistent_data.c
index 8bd9e6b5a..b278e2a32 100644
--- a/src/memory/free_persistent_data.c
+++ b/src/memory/free_persistent_data.c
@@ -13,6 +13,7 @@ void free_persistent_data(struct stardata_t * const stardata)
(void*)stardata,
(void*)(stardata?stardata->persistent_data:NULL));
+
if(stardata != NULL &&
stardata->persistent_data != NULL)
{
@@ -24,14 +25,20 @@ void free_persistent_data(struct stardata_t * const stardata)
#ifdef STELLAR_POPULATIONS_ENSEMBLE
Dprint("free ensemble hash at %p ?\n",
(void*)stardata->persistent_data->ensemble_hash);
+
if(stardata->persistent_data->ensemble_hash != NULL)
{
Dprint("call ensemble log\n");
+
stardata->model.ensemble_hash = stardata->persistent_data->ensemble_hash;
- ensemble_log(stardata,TRUE);
+ if(!Binary_c_error_code_is_bad(stardata->error_code))
+ {
+ ensemble_log(stardata,TRUE);
+ }
RCHash_free(stardata->persistent_data->ensemble_hash);
Dprint("freed: hash is now %p\n",
(void*)stardata->persistent_data->ensemble_hash);
+
stardata->persistent_data->ensemble_hash = NULL;
stardata->model.ensemble_hash = NULL;
}
@@ -40,12 +47,15 @@ void free_persistent_data(struct stardata_t * const stardata)
Dprint("freeing stardata=%p -> persistent_data = %p\n",
(void*)stardata,
(void*)stardata->persistent_data);
+
Safe_free(stardata->persistent_data);
+
}
else
{
Dprint("free persistent data called on stardata = %p stardata->persistent_data = %p\n",
(void*)stardata,
(void*)(stardata->persistent_data));
+
}
}
diff --git a/src/stellar_structure/stellar_structure_AGB.c b/src/stellar_structure/stellar_structure_AGB.c
index caf4e338e..3867874e4 100644
--- a/src/stellar_structure/stellar_structure_AGB.c
+++ b/src/stellar_structure/stellar_structure_AGB.c
@@ -152,8 +152,6 @@ Stellar_type stellar_structure_AGB(struct star_t * const oldstar,
force_supernova = FALSE;
}
/*
- * RGI
- *
* Only allow stars with CO core mass less than the
* maximum_CO_mass_for_thermal_pulses to have thermal pulses, or
* those that are already thermally pulsing.
@@ -161,10 +159,15 @@ Stellar_type stellar_structure_AGB(struct star_t * const oldstar,
* This is because more massive stars should explode
* soon after the AGB begins, and probably don't ever
* have thermal pulses.
+ *
+ * We also check that TPAGB stars have core masses
+ * < MCh: if not, treat them as EAGB i.e. allow
+ * them to explode
*/
else if((newstar->stellar_type == TPAGB ||
(Is_not_zero(newstar->core_mass[CORE_CO]) &&
- newstar->core_mass[CORE_He] < MCh &&
+ (newstar->stellar_type == EAGB ||
+ newstar->core_mass[CORE_He] < MCh) &&
newstar->core_mass[CORE_CO] < maximum_degenerate_CO_mcbagb)))
{
/* TPAGB star */
diff --git a/src/stellar_timescales/stellar_timescales_high_mass_GB.c b/src/stellar_timescales/stellar_timescales_high_mass_GB.c
index 23b6eee81..3b4542a94 100644
--- a/src/stellar_timescales/stellar_timescales_high_mass_GB.c
+++ b/src/stellar_timescales/stellar_timescales_high_mass_GB.c
@@ -2,7 +2,7 @@
No_empty_translation_unit_warning;
#ifdef BSE
-void stellar_timescales_high_mass_GB(struct stardata_t * const stardata,
+void stellar_timescales_high_mass_GB(struct stardata_t * const stardata MAYBE_UNUSED,
const double mass,
const double mp2,
double * const timescales,
diff --git a/src/timestep/stellar_timestep.c b/src/timestep/stellar_timestep.c
index 0efba17e7..d7998d83b 100644
--- a/src/timestep/stellar_timestep.c
+++ b/src/timestep/stellar_timestep.c
@@ -2,7 +2,7 @@
No_empty_translation_unit_warning;
#include "timestep.h"
-#undef MINT
+
void stellar_timestep(Timestep_prototype_args)
{
/*
@@ -133,15 +133,13 @@ void stellar_timestep(Timestep_prototype_args)
timestep_limits(Timestep_call_args);
timestep_logging(Timestep_call_args);
- timestep_fixed_timesteps(Timestep_call_args);
timestep_modulation(Timestep_call_args);
timestep_hard_limits(Timestep_call_args);
+ //nonstellar_timestep(stardata,dt);
star->stellar_timestep = *dt;
-
- if(stardata->preferences->timestep_logging == TRUE &&
- stardata->star[0].stellar_type==HG)
+ if(stardata->preferences->timestep_logging == TRUE)
{
fprintf(stdout,
"%sstellar_timestep at t=%23.15g model %d reject_same %s reject_shorten %s out star %d stellar_type %d (m=%g mc=%g L=%g R=%g RL=%g (1-R/RL=%g) minit=%g menv=%g SN=%d J*=%g; a=%g e=%g Jorb=%g, ω*=%g ωorb=%g tcirc=%g y), age=%20.12g, tm=%20.12g, tn=%20.12g >>> %sdt=%20.12g Myr%s <<< time_remaining=%g : set by %s%s%s (%d) : logtnow=%g logtnext=%g (Fequal? %s Times_are_equal? %s : diff %g) triggerd? %d %s\n",
diff --git a/src/timestep/timestep_fixed_timesteps.c b/src/timestep/timestep_fixed_timesteps.c
index d78a862be..4316a23f7 100644
--- a/src/timestep/timestep_fixed_timesteps.c
+++ b/src/timestep/timestep_fixed_timesteps.c
@@ -2,7 +2,9 @@
No_empty_translation_unit_warning;
-void timestep_fixed_timesteps(Timestep_prototype_args)
+void timestep_fixed_timesteps(struct stardata_t * const stardata,
+ struct star_t * const star,
+ double * const dtm)
{
/*
* Limit the timestep according to the various
@@ -11,7 +13,7 @@ void timestep_fixed_timesteps(Timestep_prototype_args)
int i;
for(i=0;i<FIXED_TIMESTEP_NUMBER;i++)
{
- struct binary_c_fixed_timestep_t * t = stardata->model.fixed_timesteps + i;
+ struct binary_c_fixed_timestep_t * const t = stardata->model.fixed_timesteps + i;
if(t->enabled == TRUE)
{
@@ -36,14 +38,18 @@ void timestep_fixed_timesteps(Timestep_prototype_args)
dtlim = (t->logarithmic ? exp10(t->next) : t->next) - stardata->model.time;
}
- Dprint("Set dtlim = %30.20g : next time will be %30.20g : Times_are_equal? %s\n",
+ Dprint("At t = %g : Set dtlim = %30.20g : next time will be %30.20g : Times_are_equal? %s\n",
+ T,
dtlim,
- log10(stardata->model.time + dtlim),
+ t->logarithmic ? log10(stardata->model.time + dtlim) : (stardata->model.time + dtlim),
Yesno(Times_are_equal(t->next,
log10(stardata->model.time + dtlim)))
);
- Limit_timestep(*dt,dtlim,star,DT_LIMIT_FIXED_TIMESTEP);
+ Limit_timestep(*dtm,
+ dtlim,
+ star,
+ DT_LIMIT_FIXED_TIMESTEP);
}
}
}
diff --git a/src/timestep/timestep_fixed_trigger.c b/src/timestep/timestep_fixed_trigger.c
index 5b55e9af1..e4bb103b1 100644
--- a/src/timestep/timestep_fixed_trigger.c
+++ b/src/timestep/timestep_fixed_trigger.c
@@ -7,15 +7,15 @@ Boolean timestep_fixed_trigger(struct stardata_t * RESTRICT const stardata,
const int i)
{
/*
- * Check if fixed timestep i should be triggered.
+ * Check if fixed timestep i should be triggered.
* Return TRUE if it should be, FALSE otherwise.
*
* If the timestep is disabled, return FALSE.
- *
+ *
* If we're beyond the end time, or before the start time,
* this can never be true.
*
- * If we're at max_evolution_time and t->final is TRUE, then return TRUE.
+ * If we're at max_evolution_time and t->final is TRUE, then return TRUE.
*/
Boolean ret;
struct binary_c_fixed_timestep_t * t = stardata->model.fixed_timesteps + i;
@@ -30,10 +30,10 @@ Boolean timestep_fixed_trigger(struct stardata_t * RESTRICT const stardata,
else
{
const double T = t->logarithmic ? log10(stardata->model.time) : stardata->model.time;
-
+
if(!In_range(T,t->begin,t->end))
{
- /*
+ /*
* Out of time range of the trigger so cannot trigger
*/
ret = FALSE;
@@ -50,7 +50,7 @@ Boolean timestep_fixed_trigger(struct stardata_t * RESTRICT const stardata,
else
{
/*
- * Check if the time equals, or is more than, the next
+ * Check if the time equals, or is more than, the next
* trigger time
*/
ret = Boolean_(Times_more_or_equal(T,t->next));
diff --git a/src/timestep/timestep_hard_limits.c b/src/timestep/timestep_hard_limits.c
index e1d249e99..cb70d55e0 100644
--- a/src/timestep/timestep_hard_limits.c
+++ b/src/timestep/timestep_hard_limits.c
@@ -47,12 +47,11 @@ void timestep_hard_limits(Timestep_prototype_args)
DT_LIMIT_MAXIMUM_TIMESTEP);
}
- Star_number k;
- Nuclear_burning_Starloop(k)
+ if(NUCLEAR_BURNING(star->stellar_type))
{
Limit_timestep(*dt,
stardata->preferences->maximum_nuclear_burning_timestep,
- &stardata->star[k],
+ star,
DT_LIMIT_MAXIMUM_TIMESTEP);
}
@@ -74,8 +73,7 @@ void timestep_hard_limits(Timestep_prototype_args)
/*
* If there has been an event, force a small timestep
*/
- if(0 &&
- (star->SN_type != SN_NONE ||
+ if((star->SN_type != SN_NONE ||
stardata->previous_stardata->star[star->starnum].SN_type != SN_NONE))
{
Limit_timestep(*dt,
@@ -109,16 +107,4 @@ void timestep_hard_limits(Timestep_prototype_args)
star,
DT_LIMIT_TIME_REMAINING);
}
-
- /*
- * Hard wired minimum : do not apply
- * if the timestep is limited by a fixed timestep
- */
- if(stardata->model.fixed_timestep_triggered == FALSE &&
- star->dtlimiter != DT_LIMIT_FIXED_TIMESTEP)
- {
- *dt = Max(*dt,
- MINIMUM_STELLAR_TIMESTEP);
- }
-
}
diff --git a/src/timestep/timestep_limits.c b/src/timestep/timestep_limits.c
index 3fa81b4e5..dee831040 100644
--- a/src/timestep/timestep_limits.c
+++ b/src/timestep/timestep_limits.c
@@ -648,8 +648,7 @@ void timestep_limits(Timestep_prototype_args)
* if dt_zoomfac < 1, don't let the timestep be >
* 1.2 * dt_zoomfac * (this or previous timestep)
*/
- if(0 &&
- Use_timestep(DT_LIMIT_ZOOMFAC) &&
+ if(Use_timestep(DT_LIMIT_ZOOMFAC) &&
stardata->previous_stardata &&
stardata->model.dt_zoomfac - 1.0 < -TINY)
{
diff --git a/src/timestep/timestep_prototypes.h b/src/timestep/timestep_prototypes.h
index 01db4b7c5..e51867ce0 100644
--- a/src/timestep/timestep_prototypes.h
+++ b/src/timestep/timestep_prototypes.h
@@ -23,7 +23,8 @@ void timestep_logging(Timestep_prototype_args);
void timestep_modulation(Timestep_prototype_args);
void timestep_hard_limits(Timestep_prototype_args);
void timestep_RLOF(Timestep_prototype_args);
-
+void nonstellar_timestep(struct stardata_t * const stardata,
+ double * const dt);
#ifdef DISCS
double Pure_function timestep_disc(struct stardata_t * const stardata,
@@ -31,7 +32,9 @@ double Pure_function timestep_disc(struct stardata_t * const stardata,
#endif // DISCS
void setup_fixed_timesteps(struct stardata_t * RESTRICT const stardata);
-void timestep_fixed_timesteps(Timestep_prototype_args);
+void timestep_fixed_timesteps(struct stardata_t * const stardata,
+ struct star_t * const star,
+ double * const dt);
void timestep_increment_fixed_timesteps(struct stardata_t * RESTRICT const stardata);
Boolean timestep_fixed_trigger(struct stardata_t * RESTRICT stardata,
const int i);
--
GitLab