Line data Source code
1 : ! ***********************************************************************
2 : !
3 : ! Copyright (C) 2010-2022 Bill Paxton, Pablo Marchant, Matthias Fabry
4 : ! & The MESA Team
5 : !
6 : ! This program is free software: you can redistribute it and/or modify
7 : ! it under the terms of the GNU Lesser General Public License
8 : ! as published by the Free Software Foundation,
9 : ! either version 3 of the License, or (at your option) any later version.
10 : !
11 : ! This program is distributed in the hope that it will be useful,
12 : ! but WITHOUT ANY WARRANTY; without even the implied warranty of
13 : ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
14 : ! See the GNU Lesser General Public License for more details.
15 : !
16 : ! You should have received a copy of the GNU Lesser General Public License
17 : ! along with this program. If not, see <https://www.gnu.org/licenses/>.
18 : !
19 : ! ***********************************************************************
20 :
21 : module binary_history
22 :
23 : use const_def, only: dp, pi, lsun, msun, rsun, secyer, secday, two_thirds
24 : use math_lib
25 : use binary_def
26 : use binary_private_def
27 : use binary_history_specs
28 :
29 : implicit none
30 :
31 : contains
32 :
33 0 : integer function how_many_binary_history_columns(binary_id)
34 : integer, intent(in) :: binary_id
35 : integer :: numcols, ierr
36 : type (binary_info), pointer :: b
37 :
38 : ierr = 0
39 0 : call binary_ptr(binary_id, b, ierr)
40 0 : if (ierr /= 0) then
41 0 : write(*, *) 'failed in binary_ptr'
42 0 : numcols = 0
43 0 : return
44 : end if
45 :
46 0 : if (.not. associated(b% binary_history_column_spec)) then
47 : numcols = 0
48 : else
49 0 : numcols = size(b% binary_history_column_spec, dim = 1)
50 : end if
51 :
52 0 : how_many_binary_history_columns = numcols
53 : end function how_many_binary_history_columns
54 :
55 :
56 0 : subroutine data_for_binary_history_columns(&
57 0 : binary_id, n, names, vals, ierr)
58 : integer, intent(in) :: binary_id, n
59 : character (len = 80) :: names(n)
60 : real(dp) :: vals(n)
61 : integer, intent(out) :: ierr
62 :
63 : type (binary_info), pointer :: b
64 : integer :: c, int_val, j
65 : logical :: is_int_val
66 0 : real(dp) :: val
67 :
68 : ierr = 0
69 0 : call binary_ptr(binary_id, b, ierr)
70 0 : if (ierr /= 0) then
71 0 : write(*, *) 'failed in binary_ptr'
72 0 : return
73 : end if
74 :
75 0 : do j = 1, n
76 0 : c = b% binary_history_column_spec(j)
77 0 : names(j) = trim(binary_history_column_name(c))
78 : call binary_history_getval(&
79 0 : b, c, val, int_val, is_int_val, ierr)
80 0 : if (ierr /= 0) then
81 0 : write(*, *) "Unknown binary_history_columns.list column"
82 0 : return
83 : end if
84 0 : if (is_int_val) then
85 0 : vals(j) = int_val
86 : else
87 0 : vals(j) = val
88 : end if
89 : end do
90 : end subroutine data_for_binary_history_columns
91 :
92 0 : subroutine write_binary_history_info(b, ierr)
93 : type (binary_info), pointer :: b
94 : integer, intent(out) :: ierr
95 : logical, parameter :: write_flag = .true.
96 :
97 0 : call do_binary_history_info(b, write_flag, ierr)
98 0 : end subroutine write_binary_history_info
99 :
100 0 : subroutine do_get_data_for_binary_history_columns(b, ierr)
101 : type (binary_info), pointer :: b
102 : integer, intent(out) :: ierr
103 : logical, parameter :: write_flag = .false.
104 :
105 0 : call do_binary_history_info(b, write_flag, ierr)
106 0 : end subroutine do_get_data_for_binary_history_columns
107 :
108 0 : subroutine do_binary_history_info(b, write_flag, ierr)
109 : use utils_lib, only : integer_dict_create_hash, integer_dict_free
110 : type (binary_info), pointer :: b
111 : logical, intent(in) :: write_flag
112 : integer, intent(out) :: ierr
113 :
114 : character (len = strlen) :: fname, dbl_fmt, int_fmt, txt_fmt
115 : integer :: numcols, io, i, col, j, i0, n
116 :
117 : integer :: num_extra_header_items, num_extra_cols
118 :
119 : character (len = maxlen_history_column_name), pointer, dimension(:) :: &
120 0 : extra_header_item_names, extra_col_names
121 : real(dp), pointer, dimension(:) :: &
122 0 : extra_header_item_vals, extra_col_vals
123 :
124 : logical :: binary_history_file_exists
125 : character (len = maxlen_history_column_name), pointer :: &
126 0 : names(:) ! (num_history_columns)
127 0 : real(dp), pointer :: vals(:) ! (num_history_columns)
128 0 : logical, pointer :: is_int(:) ! (num_history_columns)
129 :
130 : include 'formats'
131 :
132 0 : extra_header_item_names => null()
133 0 : extra_header_item_vals => null()
134 :
135 0 : extra_col_names => null()
136 0 : extra_col_vals => null()
137 :
138 0 : dbl_fmt = b% history_dbl_format
139 0 : int_fmt = b% history_int_format
140 0 : txt_fmt = b% history_txt_format
141 :
142 0 : ierr = 0
143 :
144 0 : if (.not. associated(b% binary_history_column_spec)) then
145 0 : numcols = 0
146 : else
147 0 : numcols = size(b% binary_history_column_spec, dim = 1)
148 : end if
149 :
150 0 : num_extra_cols = b% how_many_extra_binary_history_columns(b% binary_id)
151 0 : n = numcols + num_extra_cols
152 :
153 0 : if (n == 0) then
154 0 : write(*, *) 'WARNING: do not have any output specified for binary logs.'
155 0 : return
156 : end if
157 : ! write(*, *) " got num of cols"
158 0 : if (b% number_of_binary_history_columns < 0) then
159 0 : b% number_of_binary_history_columns = n
160 0 : else if (b% number_of_binary_history_columns /= n) then
161 0 : if (associated(b% binary_history_values)) then
162 0 : deallocate(b% binary_history_values)
163 0 : nullify(b% binary_history_values)
164 : end if
165 0 : if (associated(b% binary_history_names)) then
166 0 : deallocate(b% binary_history_names)
167 0 : nullify(b% binary_history_names)
168 : end if
169 0 : if (associated(b% binary_history_value_is_integer)) then
170 0 : deallocate(b% binary_history_value_is_integer)
171 0 : nullify(b% binary_history_value_is_integer)
172 : end if
173 0 : if (associated(b% binary_history_names_dict)) then
174 0 : call integer_dict_free(b% binary_history_names_dict)
175 0 : nullify(b% binary_history_names_dict)
176 : end if
177 0 : b% need_to_set_binary_history_names_etc = .true.
178 0 : b% number_of_binary_history_columns = n
179 : end if
180 :
181 0 : if (.not. associated(b% binary_history_values)) then
182 0 : allocate(b% binary_history_values(n))
183 0 : else if (size(b% binary_history_values, dim = 1) /= n) then
184 0 : ierr = -1
185 0 : write(*, 3) 'bad size b% binary_history_values', &
186 0 : size(b% binary_history_values, dim = 1), n
187 : end if
188 0 : vals => b% binary_history_values
189 :
190 0 : if (.not. associated(b% binary_history_names)) then
191 0 : allocate(b% binary_history_names(n))
192 0 : else if (size(b% binary_history_names, dim = 1) /= n) then
193 0 : ierr = -1
194 0 : write(*, 3) 'bad size b% binary_history_names', &
195 0 : size(b% binary_history_names, dim = 1), n
196 : end if
197 :
198 0 : if (b% need_to_set_binary_history_names_etc) then
199 0 : names => b% binary_history_names
200 : else
201 0 : nullify(names)
202 : end if
203 :
204 0 : if (.not. associated(b% binary_history_value_is_integer)) then
205 0 : allocate(b% binary_history_value_is_integer(n))
206 0 : else if (size(b% binary_history_value_is_integer, dim = 1) /= n) then
207 0 : ierr = -1
208 0 : write(*, 2) 'bad size b% binary_history_value_is_integer', &
209 0 : size(b% binary_history_value_is_integer, dim = 1), n
210 : end if
211 0 : if (b% need_to_set_binary_history_names_etc) then
212 0 : is_int => b% binary_history_value_is_integer
213 : else
214 0 : nullify(is_int)
215 : end if
216 : ! write(*, *) " associated arrays"
217 :
218 0 : nullify(extra_col_names)
219 0 : nullify(extra_col_vals)
220 :
221 0 : if (num_extra_cols > 0) then
222 : allocate(&
223 0 : extra_col_names(num_extra_cols), extra_col_vals(num_extra_cols), stat = ierr)
224 0 : if (ierr /= 0) then
225 0 : call dealloc
226 0 : return
227 : end if
228 0 : extra_col_names(1:num_extra_cols) = 'unknown'
229 0 : extra_col_vals(1:num_extra_cols) = -1d99
230 : call b% data_for_extra_binary_history_columns(&
231 0 : b% binary_id, num_extra_cols, extra_col_names, extra_col_vals, ierr)
232 0 : if (ierr /= 0) then
233 0 : call dealloc
234 0 : return
235 : end if
236 0 : do i = 1, num_extra_cols
237 0 : if(trim(extra_col_names(i))=='unknown') then
238 0 : write(*, *) "Warning empty history name for extra_binary_history_column ", i
239 : end if
240 : end do
241 : end if
242 :
243 : ! write(*, *) " starting write loop ", write_flag
244 0 : i0 = 1
245 0 : if (write_flag .and. (open_close_log .or. b% s_donor% model_number == -100)) then
246 0 : fname = trim(b% log_directory) // '/' // trim(b% history_name)
247 0 : inquire(file = trim(fname), exist = binary_history_file_exists)
248 0 : if ((.not. binary_history_file_exists) .or. b% open_new_history_file) then
249 0 : ierr = 0
250 0 : open(newunit = io, file = trim(fname), action = 'write', iostat = ierr)
251 0 : b% open_new_history_file = .false.
252 : else
253 0 : i0 = 3
254 0 : open(newunit = io, file = trim(fname), action = 'write', position = 'append', iostat = ierr)
255 : end if
256 0 : if (ierr /= 0) then
257 0 : write(*, *) 'failed to open ' // trim(fname)
258 0 : call dealloc
259 0 : return
260 : end if
261 : end if
262 :
263 0 : if (write_flag .and. i0 == 1) then ! write parameters at start of log
264 :
265 0 : num_extra_header_items = b% how_many_extra_binary_history_header_items(b% binary_id)
266 :
267 0 : if (num_extra_header_items > 0) then
268 : allocate(&
269 : extra_header_item_names(num_extra_header_items), &
270 0 : extra_header_item_vals(num_extra_header_items), stat = ierr)
271 0 : if (ierr /= 0) then
272 0 : call dealloc
273 0 : return
274 : end if
275 0 : extra_header_item_names(1:num_extra_header_items) = 'unknown'
276 0 : extra_header_item_vals(1:num_extra_header_items) = -1d99
277 : call b% data_for_extra_binary_history_header_items(&
278 : b% binary_id, num_extra_header_items, &
279 0 : extra_header_item_names, extra_header_item_vals, ierr)
280 0 : if (ierr /= 0) then
281 0 : call dealloc
282 0 : return
283 : end if
284 0 : do i = 1, num_extra_header_items
285 0 : if(trim(extra_header_item_names(i))=='unknown') then
286 0 : write(*, *) "Warning empty history name for extra_binary_history_header ", i
287 : end if
288 : end do
289 : end if
290 :
291 0 : do i = 1, 3
292 0 : col = 0
293 0 : call write_string(io, col, i, 'version_number', version_number)
294 0 : call write_val(io, col, i, 'initial_don_mass', initial_mass(1))
295 0 : call write_val(io, col, i, 'initial_acc_mass', initial_mass(2))
296 : call write_val(io, col, i, 'initial_period_days', &
297 0 : initial_binary_period / (3600 * 24))
298 :
299 0 : call write_string(io, col, i, 'compiler', compiler_name)
300 0 : call write_string(io, col, i, 'build', compiler_version_name)
301 0 : call write_string(io, col, i, 'MESA_SDK_version', mesasdk_version_name)
302 0 : call write_string(io, col, i, 'date', date)
303 :
304 0 : do j = 1, num_extra_header_items
305 : call write_val(io, col, i, &
306 0 : extra_header_item_names(j), extra_header_item_vals(j))
307 : end do
308 :
309 0 : write(io, *)
310 : end do
311 0 : write(io, *)
312 : end if
313 :
314 0 : do i = i0, 3 ! add a row to the log
315 0 : col = 0
316 : ! write(*, *) "doing cols pass", i
317 0 : do j = 1, numcols
318 0 : call do_col(i, j)
319 : end do
320 0 : do j = 1, num_extra_cols
321 0 : call do_extra_col(i, j)
322 : end do
323 0 : if (write_flag) write(io, *)
324 : end do
325 : ! write (*, *) "cols handled"
326 0 : if (open_close_log) close(io)
327 :
328 0 : call dealloc
329 : ! write(*, *) "history written"
330 :
331 0 : b% model_number_of_binary_history_values = b% model_number
332 :
333 0 : if (b% need_to_set_binary_history_names_etc) then
334 : ! write(*, *) " creating hash of the history dict"
335 0 : call integer_dict_create_hash(b% binary_history_names_dict, ierr)
336 0 : if (ierr /= 0) then
337 0 : write(*, *) "hash failed"
338 0 : return
339 : end if
340 : end if
341 :
342 0 : b% need_to_set_binary_history_names_etc = .false.
343 :
344 : contains
345 :
346 :
347 0 : subroutine dealloc
348 0 : if (associated(extra_header_item_names)) deallocate(extra_header_item_names)
349 0 : if (associated(extra_header_item_vals)) deallocate(extra_header_item_vals)
350 0 : if (associated(extra_col_names)) deallocate(extra_col_names)
351 0 : if (associated(extra_col_vals)) deallocate(extra_col_vals)
352 0 : end subroutine dealloc
353 :
354 :
355 0 : subroutine do_extra_col(pass, j)
356 : integer, intent(in) :: pass, j
357 0 : if (pass == 1) then
358 0 : if (write_flag) write(io, fmt = int_fmt, advance = 'no') j + numcols
359 0 : else if (pass == 2) then
360 0 : call do_name(j + numcols, extra_col_names(j))
361 0 : else if (pass == 3) then
362 0 : call do_val(j + numcols, extra_col_vals(j))
363 : end if
364 0 : end subroutine do_extra_col
365 :
366 :
367 0 : subroutine do_name(j, col_name)
368 : use utils_lib, only: integer_dict_define
369 : integer, intent(in) :: j
370 : character (len = *), intent(in) :: col_name
371 0 : if (write_flag) write(io, fmt = txt_fmt, advance = 'no') trim(col_name)
372 0 : if (associated(names)) names(j) = trim(col_name)
373 0 : if (b% need_to_set_binary_history_names_etc) then
374 0 : call integer_dict_define(b% binary_history_names_dict, col_name, j, ierr)
375 0 : if (ierr /= 0) write(*, *) 'failed in dict define ' // trim(col_name)
376 : end if
377 :
378 0 : end subroutine do_name
379 :
380 :
381 0 : subroutine do_col(pass, j)
382 : integer, intent(in) :: pass, j
383 0 : if (pass == 1) then
384 0 : call do_col_pass1
385 0 : else if (pass == 2) then
386 0 : call do_col_pass2(j)
387 0 : else if (pass == 3) then
388 0 : call do_col_pass3(b% binary_history_column_spec(j))
389 : end if
390 0 : end subroutine do_col
391 :
392 :
393 0 : subroutine do_col_pass1 ! write the column number
394 0 : col = col + 1
395 0 : if (write_flag) write(io, fmt = int_fmt, advance = 'no') col
396 0 : end subroutine do_col_pass1
397 :
398 :
399 0 : subroutine do_col_pass2(j) ! get the column name
400 : integer, intent(in) :: j
401 : character (len = 100) :: col_name
402 : integer :: c
403 0 : c = b% binary_history_column_spec(j)
404 0 : col_name = trim(binary_history_column_name(c))
405 0 : call do_name(j, col_name)
406 0 : end subroutine do_col_pass2
407 :
408 :
409 0 : subroutine do_col_pass3(c) ! get the column value
410 : integer, intent(in) :: c
411 : integer :: k, int_val
412 : logical :: is_int_val
413 : real(dp) :: val
414 : int_val = 0; val = 0; is_int_val = .false.
415 : call binary_history_getval(&
416 0 : b, c, val, int_val, is_int_val, ierr)
417 0 : if (ierr /= 0) then
418 0 : write(*, *) 'missing log info for ' // trim(binary_history_column_name(c)), j, k
419 0 : return
420 : end if
421 0 : if (is_int_val) then
422 0 : call do_int_val(j, int_val)
423 : else
424 0 : call do_val(j, val)
425 : end if
426 : end subroutine do_col_pass3
427 :
428 :
429 0 : subroutine do_val(j, val)
430 : use utils_lib, only : is_bad
431 : integer, intent(in) :: j
432 : real(dp), intent(in) :: val
433 0 : if (write_flag) then
434 0 : if (is_bad(val)) then
435 0 : write(io, fmt = dbl_fmt, advance = 'no') -1d99
436 : else
437 0 : write(io, fmt = dbl_fmt, advance = 'no') val
438 : end if
439 : end if
440 0 : if (associated(vals)) vals(j) = val
441 0 : if (associated(is_int)) is_int(j) = .false.
442 0 : end subroutine do_val
443 :
444 :
445 0 : subroutine do_int_val(j, val)
446 : integer, intent(in) :: j
447 : integer, intent(in) :: val
448 0 : if (write_flag) write(io, fmt = int_fmt, advance = 'no') val
449 0 : if (associated(vals)) vals(j) = dble(val)
450 0 : if (associated(is_int)) is_int(j) = .true.
451 0 : end subroutine do_int_val
452 :
453 :
454 : subroutine write_integer(io, col, pass, name, val)
455 : integer, intent(in) :: io, pass
456 : integer, intent(inout) :: col
457 : character (len = *), intent(in) :: name
458 : integer, intent(in) :: val
459 : if (pass == 1) then
460 : col = col + 1
461 : write(io, fmt = int_fmt, advance = 'no') col
462 : else if (pass == 2) then
463 : write(io, fmt = txt_fmt, advance = 'no') trim(name)
464 : else if (pass == 3) then
465 : write(io, fmt = int_fmt, advance = 'no') val
466 : end if
467 : end subroutine write_integer
468 :
469 :
470 0 : subroutine write_val(io, col, pass, name, val) ! for header items only
471 : integer, intent(in) :: io, pass
472 : integer, intent(inout) :: col
473 : character (len = *), intent(in) :: name
474 : real(dp), intent(in) :: val
475 0 : if (pass == 1) then
476 0 : col = col + 1
477 0 : write(io, fmt = int_fmt, advance = 'no') col
478 0 : else if (pass == 2) then
479 0 : write(io, fmt = txt_fmt, advance = 'no') trim(name)
480 0 : else if (pass == 3) then
481 0 : write(io, fmt = dbl_fmt, advance = 'no') val
482 : end if
483 0 : end subroutine write_val
484 :
485 :
486 0 : subroutine write_string(io, col, pass, name, val) !for header items only
487 : integer, intent(in) :: io, pass
488 : integer, intent(inout) :: col
489 : character(len = *), intent(in) :: name, val
490 : character(len = strlen) :: my_val
491 :
492 0 : my_val = '"' // trim(val) // '"'
493 0 : if (pass == 1) then
494 0 : col = col + 1
495 0 : write(io, fmt = int_fmt, advance = 'no') col
496 0 : else if (pass == 2) then
497 0 : write(io, fmt = txt_fmt, advance = 'no') trim(name)
498 0 : else if (pass == 3) then
499 0 : write(io, fmt = txt_fmt, advance = 'no') trim(my_val)
500 : end if
501 0 : end subroutine write_string
502 :
503 :
504 : end subroutine do_binary_history_info
505 :
506 :
507 0 : subroutine binary_history_getval(b, c, val, int_val, is_int_val, ierr)
508 : type (binary_info), pointer :: b
509 : integer, intent(in) :: c
510 : real(dp), intent(out) :: val
511 : integer, intent(out) :: int_val
512 : logical, intent(out) :: is_int_val
513 : integer, intent(out) :: ierr
514 :
515 : include 'formats'
516 :
517 0 : ierr = 0
518 0 : is_int_val = .false.
519 0 : int_val = 0
520 0 : val = 0
521 0 : select case(c)
522 :
523 : case(bh_model_number)
524 0 : int_val = b% model_number
525 0 : is_int_val = .true.
526 : case(bh_age)
527 0 : val = b% binary_age
528 : case(bh_donor_index)
529 0 : int_val = b% d_i
530 0 : is_int_val = .true.
531 : case(bh_period_days)
532 0 : val = b% period / secday
533 : case(bh_period_hr)
534 0 : val = b% period / (60d0 * 60d0)
535 : case(bh_period_minutes)
536 0 : val = b% period / 60d0
537 : case(bh_lg_separation)
538 0 : val = safe_log10(b% separation)
539 : case(bh_binary_separation)
540 0 : val = b% separation / Rsun
541 : case(bh_eccentricity)
542 0 : val = b% eccentricity
543 : case(bh_star_1_radius)
544 0 : val = b% r(1) / Rsun
545 : case(bh_star_2_radius)
546 0 : val = b% r(2) / Rsun
547 : case(bh_rl_1)
548 0 : val = b% rl(1) / Rsun
549 : case(bh_rl_2)
550 0 : val = b% rl(2) / Rsun
551 : case(bh_rl_overflow_1)
552 0 : val = (b% r(1) - b% rl(1)) / Rsun
553 : case(bh_rl_overflow_2)
554 0 : val = (b% r(2) - b% rl(2)) / Rsun
555 : case(bh_rl_relative_overflow_1)
556 0 : val = b% rl_relative_gap(1)
557 : case(bh_rl_relative_overflow_2)
558 0 : val = b% rl_relative_gap(2)
559 : case(bh_P_rot_div_P_orb_1)
560 0 : if (b% point_mass_i /= 1) then
561 0 : val = 2 * pi / b% s1% omega_avg_surf / b% period
562 : else
563 : val = 0.0d0
564 : end if
565 : case(bh_P_rot_div_P_orb_2)
566 0 : if (b% point_mass_i /= 2) then
567 0 : val = 2 * pi / b% s2% omega_avg_surf / b% period
568 : else
569 0 : if (.not. b% model_twins_flag) then
570 : val = 0.0d0
571 : else
572 0 : val = 2 * pi / b% s1% omega_avg_surf / b% period
573 : end if
574 : end if
575 : case(bh_lg_t_sync_1)
576 0 : val = safe_log10(abs(b% t_sync_1) / secyer)
577 : case(bh_lg_t_sync_2)
578 0 : val = safe_log10(abs(b% t_sync_2) / secyer)
579 : case(bh_star_1_mass)
580 0 : val = b% m(1) / Msun
581 : case(bh_lg_star_1_mass)
582 0 : val = safe_log10(b% m(1) / Msun)
583 : case(bh_star_2_mass)
584 0 : val = b% m(2) / Msun
585 : case(bh_lg_star_2_mass)
586 0 : val = safe_log10(b% m(2) / Msun)
587 : case(bh_sum_of_masses)
588 0 : val = (b% m(1) + b% m(2)) / Msun
589 : case(bh_mass_ratio)
590 0 : val = b% m(2) / b% m(1)
591 : case(bh_obs_mass_ratio)
592 0 : val = min(b% m(2) / b% m(1), b% m(1) / b% m(2))
593 : case(bh_lg_mtransfer_rate)
594 0 : val = safe_log10(abs(b% step_mtransfer_rate) / Msun * secyer)
595 : case(bh_lg_mstar_dot_1)
596 0 : val = safe_log10(abs(b% component_mdot(1)) / Msun * secyer)
597 : case(bh_lg_mstar_dot_2)
598 0 : val = safe_log10(abs(b% component_mdot(2)) / Msun * secyer)
599 : case(bh_lg_system_mdot_1)
600 0 : val = safe_log10(abs(b% mdot_system_transfer(1)) / Msun * secyer)
601 : case(bh_lg_system_mdot_2)
602 0 : val = safe_log10(abs(b% mdot_system_transfer(2)) / Msun * secyer)
603 : case(bh_lg_wind_mdot_1)
604 0 : val = safe_log10(abs(b% mdot_system_wind(1)) / Msun * secyer)
605 : case(bh_lg_wind_mdot_2)
606 0 : val = safe_log10(abs(b% mdot_system_wind(2)) / Msun * secyer)
607 : case(bh_star_1_div_star_2_mass)
608 0 : val = b% m(1) / b% m(2)
609 : case(bh_delta_star_1_mass)
610 0 : val = b% m(1) - initial_mass(1)
611 : case(bh_delta_star_2_mass)
612 0 : val = b% m(2) - initial_mass(2)
613 : case(bh_lg_F_irr)
614 0 : val = safe_log10(b% s_donor% irradiation_flux)
615 : case(bh_fixed_xfer_fraction)
616 0 : val = b% fixed_xfer_fraction
617 : case(bh_eff_xfer_fraction)
618 0 : if (b% component_mdot(b% d_i) == 0d0) then
619 0 : val = 1d0
620 : else
621 0 : val = (-b% component_mdot(b% a_i)) / (b% component_mdot(b% d_i))
622 : end if
623 : case(bh_lg_mdot_edd)
624 0 : if (b% limit_retention_by_mdot_edd) then
625 0 : val = safe_log10(b% mdot_edd / Msun * secyer)
626 : else
627 0 : val = safe_log10(0d0)
628 : end if
629 : case(bh_mdot_edd_eta)
630 0 : if (b% limit_retention_by_mdot_edd) then
631 0 : val = b% mdot_edd_eta
632 : else
633 : val = 0d0
634 : end if
635 : case(bh_lg_accretion_luminosity)
636 0 : val = safe_log10(b% accretion_luminosity / Lsun)
637 : case(bh_bh_spin)
638 0 : if (b% point_mass_i /= 0) then
639 : val = sqrt(two_thirds) &
640 : * (b% eq_initial_bh_mass / min(b% m(b% point_mass_i), sqrt(6d0) * b% eq_initial_bh_mass)) &
641 : * (4d0 - sqrt(18d0 * pow2(b% eq_initial_bh_mass / &
642 0 : min(b% m(b% point_mass_i), sqrt(6d0) * b% eq_initial_bh_mass)) - 2d0))
643 : else
644 : val = 0
645 : end if
646 : case(bh_v_orb_1)
647 0 : val = 2.0d0 * pi * b% m(2) / (b% m(1) + b% m(2)) * b% separation / b% period / 1.0d5
648 : case(bh_v_orb_2)
649 0 : val = 2.0d0 * pi * b% m(1) / (b% m(1) + b% m(2)) * b% separation / b% period / 1.0d5
650 : case(bh_J_orb)
651 0 : val = b% angular_momentum_j
652 : case(bh_J_spin_1)
653 0 : if (b% point_mass_i /= 1) then
654 0 : val = b% s1% total_angular_momentum
655 : else
656 : val = 0d0
657 : end if
658 : case(bh_J_spin_2)
659 0 : if (b% point_mass_i /= 2) then
660 0 : val = b% s2% total_angular_momentum
661 : else
662 0 : if (.not. b% model_twins_flag) then
663 : val = 0d0
664 : else
665 0 : val = b% s1% total_angular_momentum
666 : end if
667 : end if
668 : case(bh_J_total)
669 0 : val = b% angular_momentum_j
670 0 : if (b% point_mass_i /= 1) &
671 0 : val = val + b% s1% total_angular_momentum
672 0 : if (b% point_mass_i /= 2) then
673 0 : val = val + b% s2% total_angular_momentum
674 0 : else if (b% model_twins_flag) then
675 0 : val = val + b% s1% total_angular_momentum
676 : end if
677 0 : val = val
678 : case(bh_Jdot)
679 0 : val = b% jdot
680 : case(bh_jdot_mb)
681 0 : val = b% jdot_mb
682 : case(bh_jdot_gr)
683 0 : val = b% jdot_gr
684 : case(bh_jdot_ml)
685 0 : val = b% jdot_ml
686 : case(bh_jdot_ls)
687 0 : val = b% jdot_ls
688 : case(bh_jdot_missing_wind)
689 0 : val = b% jdot_missing_wind
690 : case(bh_extra_jdot)
691 0 : val = b% extra_jdot
692 : case(bh_accretion_mode)
693 0 : int_val = b% accretion_mode
694 0 : is_int_val = .true.
695 : case(bh_acc_am_div_kep_am)
696 0 : val = b% acc_am_div_kep_am
697 : case(bh_edot)
698 0 : val = b% edot
699 : case(bh_edot_tidal)
700 0 : val = b% edot_tidal
701 : case(bh_edot_enhance)
702 0 : val = b% edot_enhance
703 : case(bh_extra_edot)
704 0 : val = b% extra_edot
705 : case(bh_point_mass_index)
706 0 : is_int_val = .true.
707 0 : int_val = b% point_mass_i
708 : case(bh_ignore_rlof_flag)
709 0 : is_int_val = .true.
710 0 : if (b% ignore_rlof_flag) then
711 0 : int_val = 1d0
712 : else
713 : int_val = 0d0
714 : end if
715 : case(bh_model_twins_flag)
716 0 : is_int_val = .true.
717 0 : if (b% model_twins_flag) then
718 0 : int_val = 1d0
719 : else
720 : int_val = 0d0
721 : end if
722 : case(bh_CE_flag)
723 0 : is_int_val = .true.
724 0 : if (b% CE_flag) then
725 0 : int_val = 1d0
726 : else
727 : int_val = 0d0
728 : end if
729 : case(bh_CE_lambda1)
730 0 : val = b% CE_lambda1
731 : case(bh_CE_lambda2)
732 0 : val = b% CE_lambda2
733 : case(bh_CE_Ebind1)
734 0 : val = b% CE_Ebind1
735 : case(bh_CE_Ebind2)
736 0 : val = b% CE_Ebind2
737 : case(bh_CE_num1)
738 0 : is_int_val = .true.
739 0 : int_val = b% CE_num1
740 : case(bh_CE_num2)
741 0 : is_int_val = .true.
742 0 : int_val = b% CE_num2
743 :
744 : case default
745 0 : ierr = -1
746 :
747 : end select
748 :
749 0 : end subroutine binary_history_getval
750 :
751 0 : subroutine get_binary_history_specs(b, num, names, specs, report)
752 :
753 : use utils_lib
754 : use utils_def
755 :
756 : type (binary_info), pointer :: b
757 : integer, intent(in) :: num
758 : character (len = *), intent(in) :: names(:)
759 : integer, intent(out) :: specs(:)
760 : logical, intent(in) :: report
761 :
762 : integer :: i, ierr, n, j, iounit, t
763 : character (len = strlen) :: buffer, string
764 :
765 : include 'formats'
766 0 : ierr = 0
767 0 : if (num <= 0) return
768 0 : iounit = -1
769 0 : specs(1:num) = 0
770 0 : do i = 1, num
771 0 : buffer = names(i)
772 0 : n = len_trim(buffer) + 1
773 0 : buffer(n:n) = ' '
774 0 : j = 0
775 0 : t = token(iounit, n, j, buffer, string)
776 0 : if (t /= name_token) then
777 0 : if (len_trim(names(i)) > 0 .and. report) &
778 0 : write(*, *) 'bad value for name of history item ' // trim(names(i))
779 0 : specs(i) = -1
780 0 : ierr = 0
781 0 : cycle
782 : end if
783 : specs(i) = do1_binary_history_spec(&
784 0 : iounit, t, n, j, string, buffer, report, ierr)
785 0 : if (ierr /= 0) then
786 0 : if (report) write(*, *) 'get_binary_history_specs failed for ' // trim(names(i))
787 0 : specs(i) = -1
788 0 : ierr = 0
789 : end if
790 : end do
791 :
792 : end subroutine get_binary_history_specs
793 :
794 :
795 0 : subroutine get_binary_history_values(b, num, specs, &
796 0 : is_int_value, int_values, values, failed_to_find_value)
797 : ! note: this doesn't handle user-defined extra columns
798 :
799 : use utils_lib
800 : use utils_def
801 :
802 : type (binary_info), pointer :: b
803 : integer, intent(in) :: num
804 : integer, intent(in) :: specs(:)
805 : logical, intent(out) :: is_int_value(:)
806 : integer, intent(out) :: int_values(:)
807 : real(dp), intent(inout) :: values(:)
808 : logical, intent(out) :: failed_to_find_value(:)
809 :
810 : integer :: i, c, ierr
811 :
812 : include 'formats'
813 0 : ierr = 0
814 :
815 0 : do i = 1, num
816 0 : failed_to_find_value(i) = .false.
817 0 : c = specs(i)
818 0 : if (c <= 0) then
819 0 : failed_to_find_value(i) = .true.
820 : else
821 : call binary_history_getval(&
822 0 : b, c, values(i), int_values(i), is_int_value(i), ierr)
823 0 : if (ierr /= 0) then
824 0 : failed_to_find_value(i) = .true.
825 0 : ierr = 0
826 : end if
827 : end if
828 : end do
829 :
830 0 : end subroutine get_binary_history_values
831 :
832 0 : logical function get1_binary_hist_value(b, name, val)
833 : ! includes other_history_columns from run_binary_extras
834 : use utils_lib, only : integer_dict_lookup
835 : type (binary_info), pointer :: b
836 : character (len = *) :: name
837 : real(dp), intent(out) :: val
838 : integer :: i, ierr, num_extra_cols
839 0 : character (len = 80), pointer, dimension(:) :: extra_col_names
840 0 : real(dp), pointer, dimension(:) :: extra_col_vals
841 : include 'formats'
842 :
843 0 : get1_binary_hist_value = .false.
844 0 : call integer_dict_lookup(b% binary_history_names_dict, name, i, ierr)
845 0 : if (ierr /= 0 .or. i <= 0) return ! didn't find it
846 0 : if (associated(b% pg% pgbinary_hist)) then
847 0 : if (associated(b% pg% pgbinary_hist% vals)) then
848 0 : if (size(b% pg% pgbinary_hist% vals, dim = 1) >= i) then
849 0 : val = b% pg% pgbinary_hist% vals(i)
850 0 : get1_binary_hist_value = .true.
851 0 : return
852 : end if
853 : end if
854 : end if
855 :
856 : ! try extras
857 0 : if (associated(b% how_many_extra_binary_history_columns) .and. &
858 : associated(b% data_for_extra_binary_history_columns)) then
859 0 : num_extra_cols = b% how_many_extra_binary_history_columns(b% binary_id)
860 0 : if (num_extra_cols > 0) then
861 : allocate(&
862 : extra_col_names(num_extra_cols), &
863 0 : extra_col_vals(num_extra_cols), stat = ierr)
864 : call b% data_for_extra_binary_history_columns(&
865 0 : b% binary_id, num_extra_cols, extra_col_names, extra_col_vals, ierr)
866 0 : do i = 1, num_extra_cols
867 0 : if (extra_col_names(i) == name) then
868 0 : val = extra_col_vals(i)
869 0 : get1_binary_hist_value = .true.
870 0 : exit
871 : end if
872 : end do
873 0 : deallocate(extra_col_names, extra_col_vals)
874 0 : if (get1_binary_hist_value) return
875 : end if
876 : end if
877 :
878 0 : end function get1_binary_hist_value
879 :
880 : end module binary_history
|
