Commit 2a8356f1 authored by David Hendriks's avatar David Hendriks
Browse files

updated example population

parent 8ac42d8f
Loading
Loading
Loading
Loading
+64 −26
Original line number Diff line number Diff line
@@ -2,10 +2,8 @@ import os

# import json
# import time
# import sys

from binarycpython.utils.grid import Population
from binarycpython.utils.functions import get_help_all, get_help, create_hdf5
from binarycpython.utils.functions import get_help_all, get_help, create_hdf5, output_lines
from binarycpython.utils.custom_logging_functions import temp_dir

#########################################################
@@ -13,15 +11,9 @@ from binarycpython.utils.custom_logging_functions import temp_dir
# The use of help(<function>) is a good way to inspect what parameters are there to use
#########################################################

## Quick script to get some output
def output_lines(output):
    """
    Function that outputs the lines that were recieved from the binary_c run.
    """
    return output.splitlines()


def parse_function(self, output):
    # EXAMPLE PARSE_FUNCTION

    # extract info from the population instance

    # Get some information from the
@@ -36,14 +28,20 @@ def parse_function(self, output):
    # Create filename
    outfilename = os.path.join(data_dir, base_filename)

    parameters = ["time", "mass", "zams_mass", "probability", "radius", "stellar_type"]

    # Go over the output.
    for el in output_lines(output):
        headerline = el.split()[0]

        # CHeck the header and act accordingly
        if headerline == "MY_STELLAR_DATA":
            parameters = ["time", "mass", "zams_mass", "probability", "radius"]
            values = el.split()[1:]
            print(values)

            if not len(parameters)==len(values):
                print("Amount of column names isnt equal to amount of columns")
                raise ValueError

            if not os.path.exists(outfilename):
                with open(outfilename, "w") as f:
@@ -52,7 +50,6 @@ def parse_function(self, output):
            with open(outfilename, "a") as f:
                f.write(seperator.join(values) + "\n")


# Create population object
example_pop = Population()

@@ -70,9 +67,14 @@ example_pop.set(
    orbital_period=45000000080,  # bse_options
    max_evolution_time=15000,  # bse_options
    eccentricity=0.02,  # bse_options

    # Set companion to low mass
    M_2=0.08, # Since in the example we run a single system, we should set the companion mass here. If we donm't do this, the code will complain. 

    # grid_options
    amt_cores=1,  # grid_options
    amt_cores=2,  # grid_options
    verbose=1,  # verbosity. Not fully configured correctly yet but having it value of 1 prints alot of stuff

    # Custom options # TODO: need to be set in grid_options probably
    data_dir=os.path.join(
        temp_dir(), "example_python_population_result"
@@ -80,7 +82,6 @@ example_pop.set(
    base_filename="example_pop.dat",  # custom_options
)


# Creating a parsing function
example_pop.set(
    parse_function=parse_function,  # Setting the parse function thats used in the evolve_population
@@ -97,20 +98,21 @@ example_pop.set(
# )


# Log the moment when the star turns into a hertzsprung-gap
# Log the moment when the star turns into neutron
example_pop.set(
    C_logging_code="""
if(stardata->star[0].stellar_type >= 2)    
if(stardata->star[0].stellar_type >= 13)    
{
    if (stardata->model.time < stardata->model.max_evolution_time)
    {
        Printf("MY_STELLAR_DATA %30.12e %g %g %g %g\\n",
        Printf("MY_STELLAR_DATA %30.12e %g %g %g %g %d\\n",
            // 
            stardata->model.time, // 1
            stardata->star[0].mass, // 2
            stardata->star[0].pms_mass, //4
            stardata->common.zero_age.mass[0], // 4
            stardata->model.probability, // 5
            stardata->star[0].radius // 6
            stardata->star[0].radius, // 6
            stardata->star[0].stellar_type // 7
      );
    };
    /* Kill the simulation to save time */
@@ -120,7 +122,11 @@ if(stardata->star[0].stellar_type >= 2)
)

# Add grid variables
resolution = {"M_1": 100}
resolution = {
    "M_1": 20, 
    'q': 20, 
    'per': 40
}

# Mass
example_pop.add_grid_variable(
@@ -136,6 +142,38 @@ example_pop.add_grid_variable(
    condition="",  # Impose a condition on this grid variable. Mostly for a check for yourself
)

# # Mass ratio
# test_pop.add_grid_variable(
#     name="q",
#     longname="Mass ratio",
#     valuerange=["0.1/M_1", 1],
#     resolution="{}".format(resolution['q']),
#     spacingfunc="const(0.1/M_1, 1, {})".format(resolution['q']),
#     probdist="flatsections(q, [{'min': 0.1/M_1, 'max': 1.0, 'height': 1}])",
#     dphasevol="dq",
#     precode="M_2 = q * M_1",
#     parameter_name="M_2",
#     condition="",  # Impose a condition on this grid variable. Mostly for a check for yourself    
# )

# # 
# test_pop.add_grid_variable(
#    name="log10per", # in days
#    longname="log10(Orbital_Period)", 
#    valuerange=[0.15, 5.5],
#    resolution="{}".format(resolution["per"]),
#    spacingfunc="const(0.15, 5.5, {})".format(resolution["per"]),
#    precode="""orbital_period = 10** log10per
# sep = calc_sep_from_period(M_1, M_2, orbital_period)
# sep_min = calc_sep_from_period(M_1, M_2, 10**0.15)
# sep_max = calc_sep_from_period(M_1, M_2, 10**5.5)""", 
#    probdist="sana12(M_1, M_2, sep, orbital_period, sep_min, sep_max, math.log10(10**0.15), math.log10(10**5.5), -0.55)",
#    parameter_name="orbital_period",
#    dphasevol="dlog10per",
# )



# Exporting of all the settings can be done with .export_all_info()
# on default it exports everything, but can be supressed by turning it off:
#   population settings (bse_options, grid_options, custom_options), turn off with include_population
@@ -152,12 +190,12 @@ example_pop.export_all_info()

## Executing a single system
## This uses the M_1 orbital period etc set with the set function
output = example_pop.evolve_single()
print(output)
# output = example_pop.evolve_single()
# print(output)

## Executing a population
## This uses the values generated by the grid_variables
# example_pop.evolve_population_mp_chunks() # TODO: update this function call
example_pop.evolve() # TODO: update this function call

# Wrapping up the results to an hdf5 file can be done by using the create_hdf5
# (<directory containing data and settings>) This function takes the settings file