POSYDON-code / POSYDON

POSYDON is a next-generation single and binary-star population synthesis code incorporating full stellar structure and evolution modeling with the use of MESA.
BSD 3-Clause "New" or "Revised" License
30 stars 19 forks source link

Evolving a binary from an arbitrary combination #354

Open maxbriel opened 4 months ago

maxbriel commented 4 months ago

In V1, it was possible to evolve a binary from an arbitrary state, but in V2 you can only start a binary on a grid. For example, the system in the tutorial NS+H-rich_core-h-burning can start on the CO-HMS_RLO grid, but not in the detached step.

This is caused by some values not properly being initialised. I think it should be possible for someone to start a binary with a CO + a fresh star (ZAMS). This is not currently possible.

maxbriel commented 1 month ago

Although I didn't look into solving this issue, it's somewhat related to #387. Parameters are not initialised properly, despite the data being available. This post explains a workaround + the exact problem:

For the single star matching to work, we require a numeric value for the values: 'mass', 'center_h1', 'log_R', 'he_core_mass' If such a value is not given when initialising the binary, these values remain None. As such the matching crashes by failing on the comparison between None + float with the latter being grid values.

The workaround is to provide these values at the binary initialisation. This can be an actual value, or np.nan Below is an example to initialise a NS+HMS binary at solar metallicity. In that case center_h1 and he_core_mass are known, and only log_R has to be set to np.nan.

The traceback --------------------------------------------------------------------------- TypeError Traceback (most recent call last) Cell In[44], line 1 ----> 1 BINARY.evolve() File ~/Documents/POSYDON/posydon/binary_evol/binarystar.py:210, in BinaryStar.evolve(self) 206 while (self.event != 'END' and self.event != 'FAILED' 207 and self.event not in self.properties.end_events 208 and self.state not in self.properties.end_states): 209 signal.alarm(MAXIMUM_STEP_TIME) --> 210 self.run_step() 212 n_steps += 1 213 if max_n_steps is not None: File ~/Documents/POSYDON/posydon/binary_evol/binarystar.py:239, in BinaryStar.run_step(self) 234 raise ValueError( 235 "Next step name '{}' does not correspond to a function in " 236 "SimulationProperties.".format(next_step_name)) 238 self.properties.pre_step(self, next_step_name) --> 239 next_step(self) 240 finally: 241 self.append_state() File ~/Documents/POSYDON/posydon/binary_evol/DT/step_detached.py:1159, in detached_step.__call__(self, binary) 1156 return interp1d, m0, t0 1158 # get the matched data of two stars, respectively -> 1159 interp1d_sec, m0, t0 = get_star_data( 1160 binary, secondary, primary, secondary.htrack, co=False) 1162 primary_not_normal = (primary.co) or (self.non_existent_companion in [1,2]) 1163 primary_normal = (not primary.co) and self.non_existent_companion == 0 File ~/Documents/POSYDON/posydon/binary_evol/DT/step_detached.py:1083, in detached_step.__call__..get_star_data(binary, star1, star2, htrack, co, copy_prev_m0, copy_prev_t0) 1081 else: 1082 t_before_matching = time.time() -> 1083 m0, t0, htrack = self.match_to_single_star(star1, htrack) 1084 t_after_matching = time.time() 1085 if self.verbose or self.verbose == 1: File ~/Documents/POSYDON/posydon/binary_evol/DT/step_detached.py:731, in detached_step.match_to_single_star(self, star, htrack) 729 print('posydon_attributes', posydon_attributes) 730 print(htrack) --> 731 x0 = get_root0(MESA_labels, posydon_attributes, htrack, rs=rs) 733 def sq_diff_function(x): 734 return self.square_difference( 735 x, htrack=htrack, mesa_labels=MESA_labels, 736 posydon_attributes=posydon_attributes, 737 colscalers=colscalers, scales=scales) File ~/Documents/POSYDON/posydon/binary_evol/DT/step_detached.py:612, in detached_step.get_root0(self, keys, x, htrack, rs) 610 idx = np.argmax(np.asanyarray(keys)[:, None] == self.root_keys, axis=1) 611 X = self.rootm[:, :, idx] --> 612 d = np.linalg.norm((X - x[None, None, :]) / rs[None, None, :], axis=-1) 613 idx = np.unravel_index(d.argmin(), X.shape[:-1]) 614 t = self.rootm[idx][np.argmax("age" == self.root_keys)] TypeError: unsupported operand type(s) for -: 'float' and 'NoneType'