diff --git a/src/binary_c_event_structures.h b/src/binary_c_event_structures.h
index 9aba3ebf861ebc0d5ef1dc355e4691ca3afc0f89..c08af9c30f1f3e99f5a73e66cf5fbc8f4a1ff655 100644
--- a/src/binary_c_event_structures.h
+++ b/src/binary_c_event_structures.h
@@ -60,5 +60,12 @@ struct binary_c_unstable_RLOF_event_t {
int instability;
};
+struct binary_c_nova_event_t {
+ struct star_t * donor;
+ struct star_t * accretor;
+ double accretion_rate;
+ double steady_burn_rate;
+ double f;
+};
#endif // BINARY_C_EVENT_STRUCUTRES_H
diff --git a/src/binary_c_parameters.h b/src/binary_c_parameters.h
index c3722cb76759759541d892165545edd3bf64c564..253a2b10932f68e81c648a653ebd7e7b6dcff9cb 100644
--- a/src/binary_c_parameters.h
+++ b/src/binary_c_parameters.h
@@ -726,7 +726,8 @@
* rather than time-averaged derivatives.
* (Experimental!)
*/
-//#define INDIVIDUAL_NOVAE
+#define INDIVIDUAL_NOVAE
+#define INDIVIDUAL_NOVAE_EVENTS
/* deprecated: accretion from nova explosions */
//#define NOVA_BACK_ACCRETION
diff --git a/src/binary_star_functions/limit_accretion_rates.c b/src/binary_star_functions/limit_accretion_rates.c
index ead3808741db789ee737c11c89a3aed9607db013..322c2311f3c9a1fe0dfcd7bd56de921b9f2a17d2 100644
--- a/src/binary_star_functions/limit_accretion_rates.c
+++ b/src/binary_star_functions/limit_accretion_rates.c
@@ -362,11 +362,53 @@ void limit_accretion_rates(struct stardata_t * RESTRICT const stardata)
/*
* Increase the surface hydrogen layer mass.
*/
-
accretor->derivative[DERIVATIVE_STELLAR_H_LAYER_MASS] = accretion_rate;
- Dprint("dm NOVA %g\n",
+ printf("dm NOVA %g\n",
accretor->dm_novaH);
+
+#ifdef INDIVIDUAL_NOVAE_EVENTS
+ printf("cf DM %g vs dMCRIT %g\n",
+ accretor->dm_novaH,
+ dMcrit);
+ /*
+ * Check if we have sufficient mass for explosion
+ */
+ if(accretor->dm_novaH > dMcrit)
+ {
+ /*
+ * Trigger a nova event
+ */
+ struct binary_c_nova_event_t * new_nova =
+ Malloc(sizeof(struct binary_c_nova_event_t));
+ new_nova->accretor = accretor;
+ new_nova->donor = donor;
+ new_nova->accretion_rate = accretion_rate;
+ new_nova->steady_burn_rate = steady_burn_rate;
+ new_nova->f = f;
+ printf("novae event\n");
+ if(Add_new_event(stardata,
+ BINARY_C_EVENT_NOVA,
+ &nova_event_handler,
+ &nova_erase_event_handler,
+ new_nova,
+ NORMAL_EVENT)
+ != BINARY_C_EVENT_DENIED
+ )
+ {
+ /*
+ * Nova event added
+ */
+ }
+ else
+ {
+ /*
+ * Failed to allocate new event data
+ */
+ Safe_free(new_nova);
+ }
+ }
+#else
/*
* Check if we have sufficient mass for explosion
@@ -427,6 +469,8 @@ void limit_accretion_rates(struct stardata_t * RESTRICT const stardata)
accretor->novastate = NOVA_STATE_TRIGGER;
}
}
+#endif // INDIVIDUAL_NOVAE_EVENTS
+
else
{
/*
diff --git a/src/librinterpolate/rinterpolate.c b/src/librinterpolate/rinterpolate.c
index 1b93a2214e5cf645d14c9740406f8a1621ac5159..504d5e5796c53dd9ef5e02cf317fc1998e97d354 100644
--- a/src/librinterpolate/rinterpolate.c
+++ b/src/librinterpolate/rinterpolate.c
@@ -192,7 +192,9 @@ struct rinterpolate_data_t * rinterpolate(
*/
if(unlikely(rinterpolate_data==NULL))
{
+ printf("rinterpolate call alloc dataspace\n");
rinterpolate_alloc_dataspace(&rinterpolate_data);
+ printf("rinterpolate post call alloc dataspace %p\n",rinterpolate_data);
table_id = -1;
}
else
@@ -345,4 +347,4 @@ struct rinterpolate_data_t * rinterpolate(
#endif // __HAVE_LIBRINTERPOLATE__
-
\ No newline at end of file
+
diff --git a/src/logging/log_every_timestep.c b/src/logging/log_every_timestep.c
index 5bbbd3517068dd1d4e1ae7ebb48077d13458a2a8..3a34272029604375de69b80f679f2af0f0deab45 100644
--- a/src/logging/log_every_timestep.c
+++ b/src/logging/log_every_timestep.c
@@ -2191,15 +2191,16 @@ void log_every_timestep(struct stardata_t * RESTRICT const stardata)
#ifdef INDIVIDUAL_NOVAE
if(WHITE_DWARF(stardata->star[0].stellar_type))
{
- printf("NOVAE %g %g %g %g %g %g %g\n",
+ printf("NOVAE %20.12e %20.12e %20.12e %20.12e %20.12e %20.12e %20.12e %20.12e %20.12e\n",
stardata->model.time,
stardata->star[0].mass,
dt,
stardata->star[1].radius / stardata->star[1].roche_radius,
- stardata->star[1].mass,
+ stardata->star[1].mass,//5
stardata->star[1].radius,
- Mdot_net(&stardata->star[1])
-
+ Mdot_net(&stardata->star[1]),
+ stardata->star[0].dm_novaH, //8
+ stardata->star[0].derivative[DERIVATIVE_STELLAR_MASS_NOVA]
);
}
#endif//INDIVIDUAL_NOVAE
@@ -2629,7 +2630,7 @@ void log_every_timestep(struct stardata_t * RESTRICT const stardata)
unresolved_stellar_magnitudes(stardata,
stardata->store,
&stardata->tmpstore->unresolved_magnitudes);
- double * unresolved_magnitudes =
+ double * const unresolved_magnitudes =
stardata->tmpstore->unresolved_magnitudes;
/*
diff --git a/src/memory/build_tmpstore_contents.c b/src/memory/build_tmpstore_contents.c
new file mode 100644
index 0000000000000000000000000000000000000000..fa8e2e844cf62e45b941a80743defc98a488ce05
--- /dev/null
+++ b/src/memory/build_tmpstore_contents.c
@@ -0,0 +1,8 @@
+#include "../binary_c.h"
+
+#include "../binary_c_collision_matrix.h"
+
+void build_tmpstore_contents(struct tmpstore_t * RESTRICT const tmpstore)
+{
+ rinterpolate_alloc_dataspace(&tmpstore->rinterpolate_data);
+}
diff --git a/src/memory/free_tmpstore.c b/src/memory/free_tmpstore.c
index e48c06a31c88831dc9853b4b6a2e21ab19172151..5ac8a97cfec8876f25cbb7c1e2237b447bfcc470 100644
--- a/src/memory/free_tmpstore.c
+++ b/src/memory/free_tmpstore.c
@@ -74,11 +74,12 @@ void free_tmpstore(struct tmpstore_t * RESTRICT const tmpstore)
Delete_data_table(tmpstore->Peters_grav_wave_merger_time_table);
- if(tmpstore->rinterpolate_data)
+ if(tmpstore->rinterpolate_data != NULL)
{
#ifdef CHECK_RINTERPOLATE
printf("FREE %p\n",tmpstore->rinterpolate_data);
#endif
+ printf("rinterpolate free %p\n",tmpstore->rinterpolate_data);
rinterpolate_free_data(tmpstore->rinterpolate_data);
Safe_free(tmpstore->rinterpolate_data);
}
diff --git a/src/memory/main_allocations.c b/src/memory/main_allocations.c
index c071defada6be99b4eda5d5e071b0a9fe8cb0604..91bd67a634424226d0a0df903dd950f51855e1cf 100644
--- a/src/memory/main_allocations.c
+++ b/src/memory/main_allocations.c
@@ -83,8 +83,9 @@ void main_allocations(struct stardata_t ** const new_stardata,
initialize_pointers(*new_stardata,(*new_stardata)->preferences,*store);
/*
- * Build the data store
+ * Build the data stores
*/
+ build_tmpstore_contents((*new_stardata)->tmpstore);
build_store_contents(*store);
/*
diff --git a/src/memory/memory_prototypes.h b/src/memory/memory_prototypes.h
index 4f7461eca483d66dea729fc7470de464657d21fa..0ae753c5fed419308d95a3c6ee8b2150de84aa4f 100644
--- a/src/memory/memory_prototypes.h
+++ b/src/memory/memory_prototypes.h
@@ -19,6 +19,7 @@ void initialize_pointers(struct stardata_t * const stardata,
struct store_t * const store);
void build_store_contents(struct store_t * RESTRICT const store);
+void build_tmpstore_contents(struct tmpstore_t * RESTRICT const tmpstore);
void free_store_contents(struct store_t * RESTRICT const store);
void free_tmpstore(struct tmpstore_t * RESTRICT const tmpstore);
void *check_aligned_memcpy_heap_source(void * dest,
diff --git a/src/novae/nova_erase_event_handler.c b/src/novae/nova_erase_event_handler.c
new file mode 100644
index 0000000000000000000000000000000000000000..687ca3f2940b0c936ac31ec9495bf6059f34b543
--- /dev/null
+++ b/src/novae/nova_erase_event_handler.c
@@ -0,0 +1,18 @@
+#include "../binary_c.h"
+
+/*
+ * Event handling function for novae: this deals with
+ * an erase action
+ */
+
+Event_handler_function nova_erase_event_handler(void * eventp,
+ struct stardata_t * stardata,
+ void * data)
+{
+ printf("nova erase event handler\n");
+
+
+
+ return NULL;
+}
+
diff --git a/src/novae/nova_event_handler.c b/src/novae/nova_event_handler.c
new file mode 100644
index 0000000000000000000000000000000000000000..9545549a17453d4a3a98969ddb2953dee8b8e9d6
--- /dev/null
+++ b/src/novae/nova_event_handler.c
@@ -0,0 +1,45 @@
+#include "../binary_c.h"
+
+#ifdef INDIVIDUAL_NOVAE
+
+/*
+ * Event handling function for novae
+ */
+
+Event_handler_function nova_event_handler(void * const eventp,
+ struct stardata_t * const stardata,
+ void * data)
+{
+ //struct binary_c_event_t * const event = eventp;
+
+ struct binary_c_nova_event_t * const nova_data =
+ (struct binary_c_nova_event_t *) data;
+
+ struct star_t * const donor = nova_data->donor;
+ struct star_t * const accretor = nova_data->accretor;
+ const double f = nova_data->f;
+
+ printf("caught nova event at %20.12g : accretor = %d with H-layer mass %g, donor = %d : f = %g\n",
+ stardata->model.time,
+ accretor->starnum,
+ accretor->dm_novaH,
+ donor->starnum,
+ f
+ );
+
+ const double dm_lost = (1.0-f) * accretor->dm_novaH;
+ const double dm_loss_rate = - dm_lost / stardata->model.dt;
+ accretor->derivative[DERIVATIVE_STELLAR_MASS_NOVA] +=
+ dm_loss_rate;
+
+ /*
+ * No hydrogen is left on the surface of the white dwarf
+ */
+ accretor->dm_novaH = 0.0;
+ accretor->mass -= dm_lost;
+ accretor->novastate = NOVA_STATE_POST;
+
+ return NULL;
+}
+
+#endif // INDIVIDUAL_NOVAE
diff --git a/src/novae/novae_prototypes.h b/src/novae/novae_prototypes.h
index bb3cbcb07b1d6176fdef58dc492eac4fefc1ab97..6c4142de1fc0bf458fdec5d79feff5d949b63482 100644
--- a/src/novae/novae_prototypes.h
+++ b/src/novae/novae_prototypes.h
@@ -20,6 +20,14 @@ double Pure_function dwarf_nova_critical_mass_transfer_rate(struct stardata_t *
struct star_t * const donor,
struct star_t * const accretor);
+#ifdef INDIVIDUAL_NOVAE
+Event_handler_function nova_event_handler(void * const eventp,
+ struct stardata_t * const stardata,
+ void * data);
+Event_handler_function nova_erase_event_handler(void * eventp,
+ struct stardata_t * stardata,
+ void * data);
+#endif // INDIVIDUAL_NOVAE
#endif // NOVAE_PROTOTYPES_H
diff --git a/src/nucsyn/nucsyn_ensemble_log.c b/src/nucsyn/nucsyn_ensemble_log.c
index afb142377ef0d5f1ab15a13957b53f5d78b5b63c..18fbcaa421d664171a152b2b4cfe2f65d6a18682 100644
--- a/src/nucsyn/nucsyn_ensemble_log.c
+++ b/src/nucsyn/nucsyn_ensemble_log.c
@@ -712,10 +712,11 @@ void nucsyn_ensemble_log(struct stardata_t * RESTRICT const stardata)
* We can use the stellar_magnitudes wrapper function to
* call the subordinate magnitude functions.
*/
- double ** magnitudes_p = NULL;
- stellar_magnitudes(stardata,star,magnitudes_p);
- double * magnitudes = * magnitudes_p;
+ double ** stellar_mags =
+ stardata->tmpstore->stellar_magnitudes;
+ stellar_magnitudes(stardata,star,&stellar_mags[star->starnum]);
+
Stellar_magnitude_index m;
for(m=0; m<NUMBER_OF_STELLAR_MAGNITUDES; m++)
{
@@ -725,9 +726,8 @@ void nucsyn_ensemble_log(struct stardata_t * RESTRICT const stardata)
* Conversion is done on output.
*/
int n = NUCSYN_ENSEMBLE_MAGNITUDE_U+m;
- ensemble_weight[n] = luminosity_from_magnitude(magnitudes[m]);
+ ensemble_weight[n] = luminosity_from_magnitude(stellar_mags[star->starnum][m]);
}
- Safe_free(magnitudes);
}
Dprint("Spectral types");
/*
diff --git a/src/nucsyn/nucsyn_parameters.h b/src/nucsyn/nucsyn_parameters.h
index 7bd43b964c6d431ff9af38233b90149ea961b683..57128d3cfc526c1a78de83c962b743a6fbe5f4e2 100644
--- a/src/nucsyn/nucsyn_parameters.h
+++ b/src/nucsyn/nucsyn_parameters.h
@@ -740,7 +740,7 @@ undances with s-process
/* define NUCSYN_GCE to modify yield output for gce.pl
* NB this just activates other defines, and turns some off.
*/
-//#define NUCSYN_GCE
+#define NUCSYN_GCE
#ifdef NUCSYN_GCE
diff --git a/src/stellar_colours/eldridge2015_magnitudes.c b/src/stellar_colours/eldridge2015_magnitudes.c
index 16528a3e98c0d1cd8ef745a62bb25410591fdfbf..7debf216ae7c8ae5990016dd8bcb8f6284dc7b68 100644
--- a/src/stellar_colours/eldridge2015_magnitudes.c
+++ b/src/stellar_colours/eldridge2015_magnitudes.c
@@ -42,12 +42,14 @@ void eldridge2015_colours(struct stardata_t * const stardata,
char *c[STELLAR_COLOUR_STRING_LENGTH] = STELLAR_COLOUR_STRINGS;
#endif
+ printf("rinterpolate alloc %p\n",stardata->tmpstore->rinterpolate_data);
rinterpolate(basel,
stardata->tmpstore->rinterpolate_data,
ELDRIDGE2015_METADATA_BASEL,
param,
mags,
0);
+ printf("rinterpolate post alloc %p\n",stardata->tmpstore->rinterpolate_data);
/* set magnitudes */
diff --git a/src/stellar_colours/stellar_magnitudes.c b/src/stellar_colours/stellar_magnitudes.c
index 5018e09916be09991ff082d270063324169ba25b..d4d21bc70d8e7d7bcbcd34b2237a26064e878f83 100644
--- a/src/stellar_colours/stellar_magnitudes.c
+++ b/src/stellar_colours/stellar_magnitudes.c
@@ -29,7 +29,7 @@ void stellar_magnitudes(struct stardata_t * const stardata,
*magnitudes_p == NULL)
{
/*
- *
+ * Allocate space
*/
*magnitudes_p = Calloc(NUMBER_OF_STELLAR_MAGNITUDES,
sizeof(double));
diff --git a/src/stellar_colours/unresolved_stellar_magnitudes.c b/src/stellar_colours/unresolved_stellar_magnitudes.c
index 71f8ad949f8cd10699dd3b283a9db27202f90bb4..f75084fa706151e43b2e93f6294c7f9ae535f20c 100644
--- a/src/stellar_colours/unresolved_stellar_magnitudes.c
+++ b/src/stellar_colours/unresolved_stellar_magnitudes.c
@@ -6,11 +6,11 @@ void unresolved_stellar_magnitudes(struct stardata_t * const stardata,
{
/*
* Calculate unresolved binary magnitudes
+ *
+ * magnitudes_p should not be NULL
*/
- if(magnitudes_p == NULL ||
- *magnitudes_p == NULL)
+ if(*magnitudes_p == NULL)
{
-
*magnitudes_p = Calloc(NUMBER_OF_STELLAR_MAGNITUDES,
sizeof(double));
}
diff --git a/tbse b/tbse
index a04192eb2a728650c6037f4719bc89bee8865701..2ef5d1c97aaca5edf29ec5efb715c29abe78b4e6 100755
--- a/tbse
+++ b/tbse
@@ -101,7 +101,7 @@ M2=1.5
# If the period (days) is given, use it.
# If the period is zero, use the separation (in Rsun) to calculate the period.
# (this requires that the masses are already set)
-ORBITAL_PERIOD=2.0
+ORBITAL_PERIOD=300
ORBITAL_SEPARATION=0.0
# Orbital eccentricity should be in the range 0.0-1.0.
@@ -827,7 +827,7 @@ NOVA_RETENTION_METHOD=NOVA_RETENTION_ALGORITHM_CONSTANT
NOVA_RETENTION_FRACTION=1e-3
# Explode individual novae?
-INDIVIDUAL_NOVAE=False
+INDIVIDUAL_NOVAE=True
# Reverse accretion beta (Shara+ 1986) :
# accrete beta * the ejected mass back onto the companion.