From 17ec4941117c4ff59d3d27066eb1a99b42305210 Mon Sep 17 00:00:00 2001
From: Robert Izzard <r.izzard@surrey.ac.uk>
Date: Mon, 9 Sep 2019 14:24:19 +0100
Subject: [PATCH] add new subroutine to set the effective and nucsyn
metallicities in a more self-consistent way
---
src/binary_c_debug.h | 2 +-
src/binary_c_macros.h | 5 ++++
src/evolution/init_common.c | 4 +--
src/evolution/initialize_parameters.c | 35 +++++++++++---------------
src/evolution/start_of_evolution.c | 1 +
src/nucsyn/nucsyn_parameters.h | 2 +-
src/perl/scripts2/yield_vs_time.pl | 14 +++++------
src/setup/argument_setting_functions.c | 13 ----------
src/setup/cmd_line_args.h | 19 ++++----------
src/setup/cmd_line_macros.h | 2 --
src/setup/parse_arguments.c | 2 ++
src/setup/set_cmd_line_macro_pairs.c | 6 ++++-
src/setup/set_metallicities.c | 25 ++++++++++++++++++
src/setup/setup_prototypes.h | 2 ++
tbse | 19 ++++++--------
15 files changed, 78 insertions(+), 73 deletions(-)
create mode 100644 src/setup/set_metallicities.c
diff --git a/src/binary_c_debug.h b/src/binary_c_debug.h
index fda8f2428..71a1068ed 100644
--- a/src/binary_c_debug.h
+++ b/src/binary_c_debug.h
@@ -87,7 +87,7 @@
*/
-//#define Debug_stop_expression ( ... )
+//#define Debug_stop_expression ( stardata->common.nucsyn_metallicity > 0 )
#undef Debug_stop_expression
/*
diff --git a/src/binary_c_macros.h b/src/binary_c_macros.h
index a244e3d77..79d20e9b1 100644
--- a/src/binary_c_macros.h
+++ b/src/binary_c_macros.h
@@ -1024,4 +1024,9 @@
#define ALPHA_ARGS_SIZE \
(NUMBER_OF_LETTERS_IN_THE_ALPHABET * NUMBER_OF_ALPHABETICAL_CASES)
+/*
+ * Default value for effective_metallicity and nucsyn_metallicity
+ */
+#define DEFAULT_TO_METALLICITY (-1.0)
+
#endif /* BINARY_MACROS_H */
diff --git a/src/evolution/init_common.c b/src/evolution/init_common.c
index 78041cff3..47b0c66a2 100644
--- a/src/evolution/init_common.c
+++ b/src/evolution/init_common.c
@@ -36,9 +36,9 @@ void init_common(struct stardata_t * RESTRICT const stardata)
sizeof(Boolean)*(NLOCKFLAGS+1)*NUMBER_OF_STARS);
#ifdef NUCSYN
- common->nucsyn_metallicity=-1.0;
+ common->nucsyn_metallicity = DEFAULT_TO_METALLICITY;
#endif
- common->effective_metallicity = - 1.0;
+ common->effective_metallicity = DEFAULT_TO_METALLICITY;
//common->random_seed=random_seed();
if(common->random_seed==0) common->random_seed=random_seed();
diff --git a/src/evolution/initialize_parameters.c b/src/evolution/initialize_parameters.c
index ba506c652..e60428d26 100644
--- a/src/evolution/initialize_parameters.c
+++ b/src/evolution/initialize_parameters.c
@@ -10,6 +10,7 @@ void initialize_parameters(struct stardata_t * RESTRICT const stardata)
struct star_t * star1 = &stardata->star[0];
struct star_t * star2 = &stardata->star[1];
+
#ifdef PRE_MAIN_SEQUENCE
double aperi;
#endif
@@ -38,9 +39,9 @@ void initialize_parameters(struct stardata_t * RESTRICT const stardata)
* Set up parameters for the BINT interpolation
* scheme
*/
- star->bint.XHc = 0.760 - 3.0 * stardata->common.metallicity;
+ star->bint.XHc = 0.760 - 3.0 * stardata->common.nucsyn_metallicity;
//0.7381; // Asplund 2009 (solar)
- star->bint.XHec = 0.240 + 2.0 * stardata->common.metallicity;
+ star->bint.XHec = 0.240 + 2.0 * stardata->common.nucsyn_metallicity;
//0.2485; // Asplund 2009 (solar)
}
#endif
@@ -60,24 +61,17 @@ void initialize_parameters(struct stardata_t * RESTRICT const stardata)
}
#endif
- /* set the initial abunds: NB these require the metallicity to be set, which means this must be
- * called after set_up_variables, otherwise you'll get a segfault (which is good! for once!)
+ /*
+ * set the initial abundance mixture using the nucsyn_metallicity
+ */
+ nucsyn_initial_abundances(star1->Xenv,
+ common->nucsyn_metallicity,
+ stardata,
+ stardata->preferences->initial_abundance_mix);
+
+ /*
+ * copy the initial abundances where they are required
*/
- if(common->nucsyn_metallicity>-TINY)
- {
- // use nucsyn_metallicity, not the metallicity of the models
- nucsyn_initial_abundances(star1->Xenv,
- common->nucsyn_metallicity,stardata,
- stardata->preferences->initial_abundance_mix);
- }
- else
- {
- nucsyn_initial_abundances(star1->Xenv,
- common->metallicity,stardata,
- stardata->preferences->initial_abundance_mix);
- }
-
- // copy the initial abundances to lots of places
Starloop(k)
{
SETstar(k);
@@ -99,7 +93,8 @@ void initialize_parameters(struct stardata_t * RESTRICT const stardata)
* overwrite the initial abundance array to be
* the true initial abundances
*/
- Copy_abundances(star1->Xenv,stardata->preferences->the_initial_abundances);
+ Copy_abundances(star1->Xenv,
+ stardata->preferences->the_initial_abundances);
#ifdef NANCHECKS
nancheck("Xsolar in initialize_parameters",common->Xsolar,ISOTOPE_ARRAY_SIZE);
diff --git a/src/evolution/start_of_evolution.c b/src/evolution/start_of_evolution.c
index 94b6d8535..bb7416a07 100644
--- a/src/evolution/start_of_evolution.c
+++ b/src/evolution/start_of_evolution.c
@@ -13,6 +13,7 @@ void start_of_evolution(struct stardata_t *stardata)
/*************************************/
/** set up variables and parameters **/
/*************************************/
+ set_metallicities(stardata);
#ifdef BSE
/* set metallicity parameters */
Dprint("Setting up local variables, calling set_metallicity_parameters\n");
diff --git a/src/nucsyn/nucsyn_parameters.h b/src/nucsyn/nucsyn_parameters.h
index 0b67abc82..c53527665 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/perl/scripts2/yield_vs_time.pl b/src/perl/scripts2/yield_vs_time.pl
index f91673cd0..61dd8ff85 100755
--- a/src/perl/scripts2/yield_vs_time.pl
+++ b/src/perl/scripts2/yield_vs_time.pl
@@ -96,7 +96,7 @@ sub all_output
print color('reset');
}
-sub parse_bse
+sub parser
{
# parser function : takes output from binary_c and
# converts it to useful statistics stored in $h.
@@ -140,6 +140,10 @@ sub parse_bse
{
process_yield_line($h,$1//0,$linearray);
}
+ else
+ {
+ print "$header @$linearray\n";
+ }
# other output is ignored
}
}
@@ -170,7 +174,7 @@ sub defaults
nmod => 0,
# use which parser function?
- parse_bse_function_pointer => \&main::parse_bse,
+ parse_bse_function_pointer => \&main::parser,
return_array_refs => 1,
threads_stack_size => 50, # MBytes
@@ -743,8 +747,7 @@ sub defaults
############################################################
# Nucleosynthesis
initial_abundance_mix => 'NUCSYN_INIT_ABUND_MIX_AG89',
- #gce_metallicity => $population->{_bse_options}{metallicity},
- third_dup => 'True',
+ third_dup => 'True',
delta_mcmin => 0.0,
lambda_min => 0.0,
no_thermohaline_mixing => 'False',
@@ -1360,9 +1363,6 @@ sub local_setup
mkdirhier($working_dir);
$population->parse_args();
- # default GCE (nucsyn) metallicity to the actual metallicity
- $population->{_bse_options}{gce_metallicity} //= $population->{_bse_options}{metallicity};
-
if($population->{_grid_options}{condor})
{
print "Running on Condor : DIR $population->{_grid_options}{condor_dir} on $population->{_grid_options}{condor_njobs} jobs\n";
diff --git a/src/setup/argument_setting_functions.c b/src/setup/argument_setting_functions.c
index eae5d13b6..7b1c7e28d 100644
--- a/src/setup/argument_setting_functions.c
+++ b/src/setup/argument_setting_functions.c
@@ -301,19 +301,6 @@ double set_metallicity(ARG_SUBROUTINE_DECLARATION)
(*c)++;
stardata->common.metallicity = atof(argv[*c]);
}
- /*
- * Set other metallicities if they haven't been set manually
- */
-#ifdef NUCSYN
- if(stardata->common.nucsyn_metallicity<-0.5)
- {
- stardata->common.nucsyn_metallicity=stardata->common.metallicity;
- }
-#endif//NUCSYN
- if(stardata->common.effective_metallicity<-0.5)
- {
- stardata->common.effective_metallicity = stardata->common.metallicity;
- }
Clamp(stardata->common.metallicity,1e-4,0.03);
return stardata->common.metallicity;
}
diff --git a/src/setup/cmd_line_args.h b/src/setup/cmd_line_args.h
index 680933dbe..57e03b802 100644
--- a/src/setup/cmd_line_args.h
+++ b/src/setup/cmd_line_args.h
@@ -364,7 +364,7 @@ struct cmd_line_arg_t
{ \
ARG_SECTION_STARS, \
"metallicity", \
- "This sets the metallicity of the stars, i.e. the amount (by mass) of matter which is not hydrogen or helium. This metallicity sets the value used in both the stellar evolution and the nucleosynthesis algorithms (if NUCSYN is defined) unless gce_metallicity is also used (see gce_metallicity).", \
+ "This sets the metallicity of the stars, i.e. the amount (by mass) of matter which is not hydrogen or helium. If you are using the BSE algorithm, this must be 1e-4 <= metallicity <= 0.03. See also nucsyn_metallicity and effective_metallicity.", \
ARG_SUBROUTINE_RETURN_FLOAT, \
"0.02", \
&set_metallicity, \
@@ -373,10 +373,10 @@ struct cmd_line_arg_t
{ \
ARG_SECTION_STARS, \
"effective_metallicity", \
- "This sets effective metallicity of stars. If not set, is just the same as metallicity. This \"effective\" metallicity is used in, say, mass-loss routines.", \
+ "This sets effective metallicity of stars as used in routines like the Schneider wind loss. If not set, or set to DEFAULT_TO_METALLICITY (==-1, the default), this is just the same as metallicity. The main difference between effective_metallicity and metallicity is the range of validity: 0 <= effective_metallicity <= 1, while metallicity's range of validity is limited by the stellar evolution algorithm (so, for BSE, is 1e-4 <= metallicity <= 0.03).", \
ARG_FLOAT, \
- "0.02", \
- EFFECTIVE_METALLICITY_VAR, \
+ "DEFAULT_TO_METALLICITY", \
+ &(stardata->common.effective_metallicity), \
1.0 \
}, \
{ \
@@ -3097,21 +3097,12 @@ struct cmd_line_arg_t
NUCSYN_REACTION_RATE_MULTIPLIERS_VAR, \
1.0 \
} , \
-{ \
- ARG_SECTION_NUCSYN, \
- "gce_metallicity", \
- "Synonym for nucsyn_metallicity, see that.", \
- ARG_FLOAT, \
- "0.02", \
- NUCSYN_METALLICITY_VAR, \
- 1.0 \
- }, \
{ \
ARG_SECTION_NUCSYN, \
"nucsyn_metallicity", \
"This sets the metallicity of the nucleosynthesis algorithms, i.e. the amount (by mass) of matter which is not hydrogen or helium. Usually you'd just set this with the metallicity parameter, but if you want the nucleosynthesis to be outside the range of the stellar evolution algorithm (e.g. Z=0 or Z=0.04) then you need to use gce_metallicity. That said, it's also outside the range of some of the nucleosynthesis algorithms as well, so you have been warned!", \
ARG_FLOAT, \
- "0.02", \
+ "DEFAULT_TO_METALLICITY", \
NUCSYN_METALLICITY_VAR, \
1.0 \
}, \
diff --git a/src/setup/cmd_line_macros.h b/src/setup/cmd_line_macros.h
index 058f89134..39fdbcaee 100644
--- a/src/setup/cmd_line_macros.h
+++ b/src/setup/cmd_line_macros.h
@@ -120,8 +120,6 @@
#define VROT_2_VAR NULL
#endif
-#define EFFECTIVE_METALLICITY_VAR &(stardata->common.effective_metallicity)
-
#ifdef WD_KICKS
#define WD_SIGMA_VAR &(stardata->preferences->sn_kick_dispersion[SN_WDKICK])
#define WD_KICK_DIRECTION_VAR &(stardata->preferences->wd_kick_direction)
diff --git a/src/setup/parse_arguments.c b/src/setup/parse_arguments.c
index ce1eff8e8..40c32771f 100644
--- a/src/setup/parse_arguments.c
+++ b/src/setup/parse_arguments.c
@@ -547,4 +547,6 @@ static void derived_arguments(struct stardata_t * const stardata)
stardata->preferences->yields_logdt)));
}
#endif // NUCSYN && NUCSYN_YIELDS
+
+
}
diff --git a/src/setup/set_cmd_line_macro_pairs.c b/src/setup/set_cmd_line_macro_pairs.c
index 32aef1fbb..ed9440adf 100644
--- a/src/setup/set_cmd_line_macro_pairs.c
+++ b/src/setup/set_cmd_line_macro_pairs.c
@@ -738,8 +738,12 @@ void set_cmd_line_macro_pairs(struct stardata_t * RESTRICT const stardata,
NOVA_RETENTION_ALGORITHM_CLAEYS2014);
Add_value_macro_pair("nova_retention_method",
NOVA_RETENTION_ALGORITHM_HILLMAN2015);
+ Add_value_macro_pair("effective_metallicity",
+ DEFAULT_TO_METALLICITY);
+ Add_value_macro_pair("nucsyn_metallicity",
+ DEFAULT_TO_METALLICITY);
-
+
/*
* ARGUMENT MACROS
*
diff --git a/src/setup/set_metallicities.c b/src/setup/set_metallicities.c
new file mode 100644
index 000000000..4c87562fa
--- /dev/null
+++ b/src/setup/set_metallicities.c
@@ -0,0 +1,25 @@
+/*
+ * This function initializes the metallicties used in various parts
+ * of binary_c
+ */
+#include "../binary_c.h"
+
+void Pure_function set_metallicities(struct stardata_t * RESTRICT const stardata)
+{
+ /*
+ * If the effective or nucsyn metallicity are < 0,
+ * set it to the metallicity
+ */
+ if(stardata->common.effective_metallicity < -TINY)
+ {
+ stardata->common.effective_metallicity =
+ stardata->common.metallicity;
+ }
+#ifdef NUCSYN
+ if(stardata->common.nucsyn_metallicity < -TINY)
+ {
+ stardata->common.nucsyn_metallicity =
+ stardata->common.metallicity;
+ }
+#endif // NUCSYN
+}
diff --git a/src/setup/setup_prototypes.h b/src/setup/setup_prototypes.h
index d482a74ff..c5248c1af 100644
--- a/src/setup/setup_prototypes.h
+++ b/src/setup/setup_prototypes.h
@@ -91,5 +91,7 @@ void new_system(struct stardata_t ** stardata,
void set_cmd_line_macro_pairs(struct stardata_t * RESTRICT const stardata,
struct cmd_line_arg_t * const cmd_line_args2,
unsigned int arg_count);
+void Pure_function set_metallicities(struct stardata_t * RESTRICT const stardata);
+
#endif /* SETUP_PROTOYPES_H */
diff --git a/tbse b/tbse
index a881c8caf..54e4f8554 100755
--- a/tbse
+++ b/tbse
@@ -122,8 +122,9 @@ METALLICITY=0.0002
# algorithm, and may be used by other algorithms
# such as mass-loss.
#
-# By default it's the same as the actual metallicity.
-EFFECTIVE_METALLICITY=$METALLICITY
+# By default it's the same as the actual metallicity,
+# which is the what happens when you set it to -1.
+EFFECTIVE_METALLICITY=-1
# Stellar structure algorithm
# 0 = STELLAR_STRUCTURE_ALGORITHM_MODIFIED_BSE = modified BSE (default)
@@ -1264,9 +1265,10 @@ MAX_STELLAR_ANGMOM_CHANGE=0.05
# Chemistry and nucleosynthesis
################################
-# metallicity in nucsyn code (defaults to $METALLICITY, but
-# can be <1e-4 or >0.03)
-NUCSYN_METALLICITY=$METALLICITY
+# metallicity in nucsyn code
+#
+# If -1, use the same as the normal metallicity
+NUCSYN_METALLICITY=-1
## new: vary initial abundances individually
#INIT_ABUNDS="--init_abund 2 0.9"
@@ -1286,10 +1288,6 @@ NUCSYN_METALLICITY=$METALLICITY
# NUCSYN_INIT_ABUND_MIX_LODDERS2010 = 8 = Lodders 2010
NUCSYN_INITIAL_ABUNDANCE_MIXTURE=NUCSYN_INIT_ABUND_MIX_AG89
-# envelope "metallicity" (used in nucleosynthesis routines only)
-# Usually the same as $Z (for 1e-4<=Z<=0.03)
-GCE_METALLICITY=$METALLICITY
-
############################################################
# Red giant transition to white dwarf algorithm
############################################################
@@ -1810,7 +1808,6 @@ $NUCMULTS \
--lambda_enthalpy $LAMBDA_ENTHALPY \
--escape_velocity 1000000000 \
--escape_fraction 0 \
---gce_metallicity $GCE_METALLICITY \
--minimum_timestep $MINIMUM_TIMESTEP \
--maximum_timestep $MAXIMUM_TIMESTEP \
--comenv_ns_accretion_fraction $CE_NS_ACC_FRAC \
@@ -1837,8 +1834,6 @@ $NUCMULTS \
--yields_logtimes $YIELDS_LOGTIMES \
--yields_startlogtime $YIELDS_STARTLOGTIME \
$HE_IG_EXTRA \
---nucsyn_metallicity $NUCSYN_METALLICITY \
---effective_metallicity $EFFECTIVE_METALLICITY \
--rotationally_enhanced_mass_loss $ROTATIONALLY_ENHANCED_MASS_LOSS \
--rotationally_enhanced_exponent $ROTATIONALLY_ENHANCED_EXPONENT \
--RLOF_timestep_modulator $RLOF_TIMESTEP_MODULATOR \
--
GitLab