|
1 | | -from matplotlib import pyplot |
2 | | -from amuse.lab import * |
3 | | -from prepare_figure import * |
4 | | -from distinct_colours import get_distinct |
5 | | -import numpy |
| 1 | +import numpy as np |
| 2 | +import matplotlib.pyplot as plt |
| 3 | +from amuse.units import units |
| 4 | +from amuse.support.console import set_printing_strategy |
| 5 | +from amuse.datamodel import Particles, Particle |
| 6 | +from amuse.community.seba import Seba |
| 7 | + |
| 8 | +set_printing_strategy( |
| 9 | + "custom", # nbody_converter = converter, |
| 10 | + preferred_units=[units.MSun, units.au, units.Myr], |
| 11 | + precision=4, |
| 12 | + prefix="", |
| 13 | + separator=" [", |
| 14 | + suffix="]", |
| 15 | +) |
6 | 16 |
|
7 | | -set_printing_strategy("custom", #nbody_converter = converter, |
8 | | - preferred_units = [units.MSun, units.AU, units.Myr], |
9 | | - precision = 4, prefix = "", |
10 | | - separator = " [", suffix = "]") |
11 | 17 |
|
12 | 18 | def maximum_stellar_radius(): |
13 | | - M = (10**numpy.arange(-0.2, 2, 0.05)) | units.MSun |
14 | | - R = [] | units.AU |
| 19 | + M = (10 ** np.arange(-0.2, 2, 0.05)) | units.MSun |
| 20 | + R = [] | units.au |
15 | 21 | for mi in M: |
16 | | - stellar = SeBa() |
| 22 | + stellar = Seba() |
17 | 23 | s = Particle() |
18 | 24 | s.mass = mi |
19 | 25 | stellar.particles.add_particle(s) |
20 | | - rmax = 0|units.AU |
21 | | - while stellar.particles[0].stellar_type<9|units.stellar_type: |
| 26 | + rmax = 0 | units.au |
| 27 | + while stellar.particles[0].stellar_type < 9 | units.stellar_type: |
22 | 28 | stellar.evolve_model() |
23 | 29 | rmax = max(rmax, stellar.particles[0].radius) |
24 | | - print("Time=", mi, stellar.particles[0].mass, \ |
25 | | - stellar.model_time, stellar.particles[0].radius) |
| 30 | + print( |
| 31 | + "Time=", |
| 32 | + mi, |
| 33 | + stellar.particles[0].mass, |
| 34 | + stellar.model_time, |
| 35 | + stellar.particles[0].radius, |
| 36 | + ) |
26 | 37 | R.append(rmax) |
27 | 38 | return M, R |
28 | | - |
29 | | -names = [r"$\tau$ CMa", "V* CQ Dra", r"$\xi$ Tau", "V* V1334 Cyg", |
30 | | - "V* DL Vir", "V* d Ser", "HD 97131", "KIC002856960"] |
| 39 | + |
| 40 | + |
| 41 | +names = [ |
| 42 | + r"$\tau$ CMa", |
| 43 | + "V* CQ Dra", |
| 44 | + r"$\xi$ Tau", |
| 45 | + "V* V1334 Cyg", |
| 46 | + "V* DL Vir", |
| 47 | + "V* d Ser", |
| 48 | + "HD 97131", |
| 49 | + "KIC002856960", |
| 50 | +] |
31 | 51 | Mprim = [50.0, 6.3, 5.5, 4.4, 3.68, 2.5, 1.5, 0.76] | units.MSun |
32 | | -Rprim = [0.95, 2.64, 0.42, 3.99, 2.34, 0.62, 0.26, 0.12]| units.AU |
| 52 | +Rprim = [0.95, 2.64, 0.42, 3.99, 2.34, 0.62, 0.26, 0.12] | units.au |
33 | 53 |
|
34 | | -colors = get_distinct(4) |
35 | | -figure = pyplot.figure(figsize=(16, 12)) |
36 | | -ax = pyplot.gca() |
37 | | -ax.minorticks_on() # switch on the minor ticks |
| 54 | +figure = plt.figure() |
| 55 | +ax = figure.add_subplot(1, 1, 1) |
| 56 | +ax.minorticks_on() # switch on the minor ticks |
38 | 57 | ax.locator_params(nbins=3) |
39 | 58 |
|
40 | 59 | M, Rmax = maximum_stellar_radius() |
41 | 60 |
|
42 | 61 | x_label = "m [$M_\odot$]" |
43 | | -y_label = "$r [AU]$" |
44 | | -pyplot.xlabel(x_label) |
45 | | -pyplot.ylabel(y_label) |
46 | | -pyplot.semilogx() |
47 | | -pyplot.semilogy() |
48 | | -pyplot.xlim(0.6, 100.0) |
49 | | -pyplot.ylim(0.07, 20.0) |
50 | | -pyplot.scatter(Mprim.value_in(units.MSun), Rprim.value_in(units.AU), |
51 | | - c=colors[0], lw=0, s=300) |
52 | | -pyplot.plot(M.value_in(units.MSun), Rmax.value_in(units.AU), |
53 | | - c=colors[1], lw=4) |
| 62 | +y_label = "$r [au]$" |
| 63 | +plt.xlabel(x_label) |
| 64 | +plt.ylabel(y_label) |
| 65 | +plt.semilogx() |
| 66 | +plt.semilogy() |
| 67 | +plt.xlim(0.6, 100.0) |
| 68 | +plt.ylim(0.07, 20.0) |
| 69 | +plt.scatter( |
| 70 | + Mprim.value_in(units.MSun), Rprim.value_in(units.au), lw=0, s=300 |
| 71 | +) |
| 72 | +plt.plot(M.value_in(units.MSun), Rmax.value_in(units.au), lw=4) |
54 | 73 | for i in range(len(names)): |
55 | | - pyplot.text(1.1*Mprim[i].value_in(units.MSun), |
56 | | - 0.9*Rprim[i].value_in(units.AU), names[i], fontsize=20) |
| 74 | + plt.text( |
| 75 | + 1.1 * Mprim[i].value_in(units.MSun), |
| 76 | + 0.9 * Rprim[i].value_in(units.au), |
| 77 | + names[i], |
| 78 | + fontsize=20, |
| 79 | + ) |
57 | 80 |
|
58 | | -save_file = 'fig_RLOF_triples_R.png' |
59 | | -pyplot.savefig(save_file) |
60 | | -print('\nSaved figure in file', save_file,'\n') |
61 | | -pyplot.show() |
| 81 | +save_file = "fig_RLOF_triples_R.png" |
| 82 | +plt.savefig(save_file) |
| 83 | +print("\nSaved figure in file", save_file, "\n") |
| 84 | +plt.show() |
0 commit comments