Skip to content

Magneto-rotational physics

mlpoppyns.simulator.magneto_rotational_physics.initial_period

Initial probability distribution of pulsar periods.

Authors:

    Vanessa Graber (graber@ice.csic.es)
    Michele Ronchi (ronchi@ice.csic.es)

pdf_period_lognormal(mean, sigma, NS_number)

Log-normal distribution for the initial spin periods as suggested in Igoshev et al. (2022). The mean and standard deviation are defined in the configuration file.

Parameters:

Name Type Description Default
mean float

Mean of the Gaussian initial period distribution, in [s].

 required
sigma float

Standard deviation of the initial period distribution, in [s].

 required
NS_number int

Total number of neutron stars created in the simulation.

 required

Returns:

Type Description
ndarray

Initial pulsar period in [s] drawn from a Log-normal distribution.
Source code in mlpoppyns/simulator/magneto_rotational_physics/initial_period.py 
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
def pdf_period_lognormal(
    mean: float, sigma: float, NS_number: int
) -> np.ndarray:
    """
    Log-normal distribution for the initial spin periods as suggested in Igoshev et al. (2022).
    The mean and standard deviation are defined in the configuration file.

    Args:
        mean (float): Mean of the Gaussian initial period distribution, in [s].
        sigma (float): Standard deviation of the initial period distribution, in [s].
        NS_number (int): Total number of neutron stars created in the simulation.

    Returns:
        (np.ndarray): Initial pulsar period in [s] drawn from a Log-normal distribution.
    """

    P_initial = 10 ** np.random.normal(mean, sigma, NS_number)

    return P_initial

pdf_period_normal(mean, sigma, NS_number)

Normal (Gaussian) distribution for the initial spin periods as suggested in Faucher-Giguère & Kaspi (2006) and Gullon et al. (2014). The mean and standard deviation are defined in the configuration file. Note that for physical reasons, we reject negative spin-periods and redraw them again from the Gaussian distribution.

Parameters:

Name Type Description Default
mean float

Mean of the Gaussian initial period distribution, in [s].

 required
sigma float

Standard deviation of the initial period distribution, in [s].

 required
NS_number int

Total number of neutron stars created in the simulation.

 required

Returns:

Type Description
ndarray

Initial pulsar period in [s] drawn from a Gaussian distribution.
Source code in mlpoppyns/simulator/magneto_rotational_physics/initial_period.py 
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
def pdf_period_normal(mean: float, sigma: float, NS_number: int) -> np.ndarray:
    """
    Normal (Gaussian) distribution for the initial spin periods as suggested in
    Faucher-Giguère & Kaspi (2006) and Gullon et al. (2014).
    The mean and standard deviation are defined in the configuration file. Note
    that for physical reasons, we reject negative spin-periods and redraw them
    again from the Gaussian distribution.

    Args:
        mean (float): Mean of the Gaussian initial period distribution, in [s].
        sigma (float): Standard deviation of the initial period distribution, in [s].
        NS_number (int): Total number of neutron stars created in the simulation.

    Returns:
        (np.ndarray): Initial pulsar period in [s] drawn from a Gaussian distribution.
    """

    P_initial = np.zeros(NS_number)

    for i in range(NS_number):
        P_initial[i] = np.random.normal(mean, sigma, 1)
        while P_initial[i] <= 0:
            P_initial[i] = np.random.normal(mean, sigma, 1)

    return P_initial

mlpoppyns.simulator.magneto_rotational_physics.magnetic_field_evolution

This module compute the evolution in time of the pulsar dipolar magnetic field.

Authors:

Vanessa Graber (graber@ice.csic.es)
Michele Ronchi (ronchi@ice.csic.es)

magnetic_field_evolution_analytical(B_initial, t, B_asymptotic, L, sigma, n_e)

Calculating the evolution of the magnetic field strength of a pulsar based on a simplified model (see eq. (17) of Aguilera et al. (2008)) that captures the characteristics of more complicated numerical simulations of pulsar magnetic field evolution, i.e., at early timescales the Hall evolution dominates, while at late times the exponential magnetic field decay due to Ohmic dissipation kicks in. Note that as explained in Aguilera et al. (2008) the Hall timescale corresponds to that of the initial field strength.

Parameters:

Name Type Description Default
B_initial float

Initial magnetic field strength in [G].

 required
t float

Time in [s].

 required
B_asymptotic float

Asymptotic magnetic field strength at late times in [G].

 required
L float

Characteristic length scale on which the magnetic field varies, measured in [cm].

 required
sigma float

Conductivity of the dominating dissipative process, measure in [1/s].

 required
n_e float

Electron density, measured in [g/cm^3].

 required

Returns:

Type Description
float

Magnetic field derivatives for a simulated pulsars in [G/yr].
Source code in mlpoppyns/simulator/magneto_rotational_physics/magnetic_field_evolution.py 
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
@jit(float64(float64, float64, float64, float64, float64, float64))
def magnetic_field_evolution_analytical(
    B_initial: float,
    t: float,
    B_asymptotic: float,
    L: float,
    sigma: float,
    n_e: float,
) -> float:
    """
    Calculating the evolution of the magnetic field strength of a pulsar based on a simplified model (see eq. (17)
    of Aguilera et al. (2008)) that captures the characteristics of more complicated numerical simulations of pulsar
    magnetic field evolution, i.e., at early timescales the Hall evolution dominates, while at late times the
    exponential magnetic field decay due to Ohmic dissipation kicks in. Note that as explained in Aguilera et al.
    (2008) the Hall timescale corresponds to that of the initial field strength.

    Args:
        B_initial (float): Initial magnetic field strength in [G].
        t (float): Time in [s].
        B_asymptotic (float): Asymptotic magnetic field strength at late times in [G].
        L (float): Characteristic length scale on which the magnetic field varies, measured in [cm].
        sigma (float): Conductivity of the dominating dissipative process, measure in [1/s].
        n_e (float): Electron density, measured in [g/cm^3].

    Returns:
        (float): Magnetic field derivatives for a simulated pulsars in [G/yr].
    """

    tau_ohm = timescale_ohmic(L, sigma)
    tau_hall = timescale_Hall(B_initial, L, n_e)

    B = (
        B_initial
        * math.exp(-t / tau_ohm)
        / (1.0 + tau_ohm / tau_hall * (1.0 - math.exp(-t / tau_ohm)))
    )

    # If the magnetic field becomes lower than an asymptotic value derived from the old millisecond pulsar population,
    # fix the magnetic field to that constant asymptotic value unless the initial field is already lower than this
    # asymptotic value.
    if B_initial < B_asymptotic:
        B = B_initial
    elif B < B_asymptotic:
        B = B_asymptotic

    return B

magnetic_field_evolution_analytical_numpy(B_initial, t, B_asymptotic, L, sigma, n_e)

Calculating the evolution of the magnetic field strength of a pulsar based on a simplified model (see eq. (17) of Aguilera et al. (2008)) that captures the characteristics of more complicated numerical simulations of pulsar magnetic field evolution, i.e., at early timescales the Hall evolution dominates, while at late times the exponential magnetic field decay due to Ohmic dissipation kicks in. Note that as explained in Aguilera et al. (2008) the Hall timescale corresponds to that of the initial field strength.

This method is compatible with NumPy arrays and is used if one wants to save the entire evolution output in the magneto_rotational_evolution method.

Parameters:

Name Type Description Default
B_initial float

Initial magnetic field strength in [G].

 required
t ndarray

Time in [s].

 required
B_asymptotic float

Asymptotic magnetic field strength at late times in [G].

 required
L float

Characteristic length scale on which the magnetic field varies, measured in [cm].

 required
sigma float

Conductivity of the dominating dissipative process, measure in [1/s].

 required
n_e float

Electron density, measured in [g/cm^3].

 required

Returns:

Type Description
float

Magnetic field derivatives for a simulated pulsars in [G/yr].
Source code in mlpoppyns/simulator/magneto_rotational_physics/magnetic_field_evolution.py 
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
def magnetic_field_evolution_analytical_numpy(
    B_initial: float,
    t: np.ndarray,
    B_asymptotic: float,
    L: float,
    sigma: float,
    n_e: float,
) -> np.ndarray:
    """
    Calculating the evolution of the magnetic field strength of a pulsar based on a simplified model (see eq. (17)
    of Aguilera et al. (2008)) that captures the characteristics of more complicated numerical simulations of pulsar
    magnetic field evolution, i.e., at early timescales the Hall evolution dominates, while at late times the
    exponential magnetic field decay due to Ohmic dissipation kicks in. Note that as explained in Aguilera et al.
    (2008) the Hall timescale corresponds to that of the initial field strength.

    This method is compatible with NumPy arrays and is used if one wants to save the entire evolution output in the
    magneto_rotational_evolution method.

    Args:
        B_initial (float): Initial magnetic field strength in [G].
        t (np.ndarray): Time in [s].
        B_asymptotic (float): Asymptotic magnetic field strength at late times in [G].
        L (float): Characteristic length scale on which the magnetic field varies, measured in [cm].
        sigma (float): Conductivity of the dominating dissipative process, measure in [1/s].
        n_e (float): Electron density, measured in [g/cm^3].

    Returns:
        (float): Magnetic field derivatives for a simulated pulsars in [G/yr].
    """

    tau_ohm = timescale_ohmic(L, sigma)
    tau_hall = timescale_Hall(B_initial, L, n_e)

    B = (
        B_initial
        * np.exp(-t / tau_ohm)
        / (1.0 + tau_ohm / tau_hall * (1.0 - np.exp(-t / tau_ohm)))
    )

    # If the magnetic field becomes lower than an asymptotic value derived from the old millisecond pulsar population,
    # fix the magnetic field to that constant asymptotic value unless the initial field is already lower than this
    # asymptotic value.
    if B_initial < B_asymptotic:
        B = np.ones(len(t)) * B_initial
    elif B_initial > B_asymptotic:
        B[B < B_asymptotic] = B_asymptotic

    return B

magnetic_field_evolution_fit(B_initial, t, B_asymptotic, a1, a2, A1, A2, b1, b2, tau_late, a_late)

An analytical fit for the magnetic field evolution curves from the magneto-thermal evolution simulations. This method is used when solving the differential equations if one wants to save only the final state in the magneto_rotational_evolution method.

The fit parameters for each magneto-thermal model specified in the config_simulator.py file were adjusted by hand (see the notebook tutorials/analysis_notebooks/magnetic_field_evolution_fit.ipynb for more details).

Parameters:

Name Type Description Default
B_initial float

Initial magnetic field strength in [G].

 required
t float

Time in [s].

 required
B_asymptotic float

Asymptotic magnetic field strength at late times in [G].

 required
a1 float

Power-law index for the first power-law component for the magnetic-field evolution fit.

 required
a2 float

Power-law index for the second power-law component for the magnetic-field evolution fit.

 required
A1 float

Normalization for the timescale parameter of the first power-law component for the magnetic-field evolution fit.

 required
A2 float

Normalization for the timescale parameter of the second power-law component for the magnetic-field evolution fit.

 required
b1 float

Power-law index for the timescale parameter of the first power-law component for the magnetic-field evolution fit.

 required
b2 float

Power-law index for the timescale parameter of the second power-law component for the magnetic-field evolution fit.

 required
tau_late float

Timescale in [yr] when transitioning from the simulated curves to the simple late-time power-law evolution of the magnetic field strength.

 required
a_late float

Power-law index for the late-time evolution of the magnetic field strength.

 required

Returns:

Type Description
float

Magnetic field value in [G] at time t.
Source code in mlpoppyns/simulator/magneto_rotational_physics/magnetic_field_evolution.py 
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
@jit(
    [
        float64(
            float64,
            float64,
            float64,
            float64,
            float64,
            float64,
            float64,
            float64,
            float64,
            float64,
            float64,
        )
    ],
    nopython=True,
)
def magnetic_field_evolution_fit(
    B_initial: float,
    t: float,
    B_asymptotic: float,
    a1: float,
    a2: float,
    A1: float,
    A2: float,
    b1: float,
    b2: float,
    tau_late: float,
    a_late: float,
) -> float:
    """
    An analytical fit for the magnetic field evolution curves from the magneto-thermal evolution simulations.
    This method is used when solving the differential equations if one wants to save only the final state in the
    magneto_rotational_evolution method.

    The fit parameters for each magneto-thermal model specified in the config_simulator.py file were adjusted by hand
    (see the notebook tutorials/analysis_notebooks/magnetic_field_evolution_fit.ipynb for more details).

    Args:
        B_initial (float): Initial magnetic field strength in [G].
        t (float): Time in [s].
        B_asymptotic (float): Asymptotic magnetic field strength at late times in [G].
        a1 (float): Power-law index for the first power-law component for the magnetic-field evolution fit.
        a2 (float): Power-law index for the second power-law component for the magnetic-field evolution fit.
        A1 (float): Normalization for the timescale parameter of the first power-law component for the
            magnetic-field evolution fit.
        A2 (float): Normalization for the timescale parameter of the second power-law component for the
            magnetic-field evolution fit.
        b1 (float): Power-law index for the timescale parameter of the first power-law component for the
            magnetic-field evolution fit.
        b2 (float): Power-law index for the timescale parameter of the second power-law component for the
            magnetic-field evolution fit.
        tau_late (float): Timescale in [yr] when transitioning from the simulated curves to the simple late-time
            power-law evolution of the magnetic field strength.
        a_late (float): Power-law index for the late-time evolution of the magnetic field strength.

    Returns:
        (float): Magnetic field value in [G] at time t.
    """
    # We consider an evolution determined by three main timescales tau1, tau2 and tau_late.
    # At early times the curves are fixed to reproduce the simulated magnetic field evolution from the magneto-thermal
    # code. At late times after a timescale tau_late we assume that the evolution is determined by a simple power-law.

    B: float = B_initial

    # Define the two timescales as a function of the initial B field.
    tau1 = A1 * B_initial**b1
    tau2 = A2 * B_initial**b2

    if tau2 < tau_late:
        B = (
            B_initial
            * (1 + t / tau1) ** a1
            * (1 + t / tau2) ** (a2 - a1)
            * (1 + t / tau_late) ** (a_late - a2)
        )

    elif (tau1 < tau_late) & (tau_late < tau2):
        B = (
            B_initial
            * (1 + t / tau1) ** a1
            * (1 + t / tau_late) ** (a_late - a1)
        )
    elif tau_late < tau1:
        B = B_initial * (1 + t / tau_late) ** a_late

    # If the magnetic field becomes lower than an asymptotic value derived from the old millisecond pulsar population,
    # fix the magnetic field to that constant asymptotic value unless the initial field is already lower than this
    # asymptotic value.
    if B_initial < B_asymptotic:
        B = B_initial
    elif B < B_asymptotic:
        B = B_asymptotic

    return B

magnetic_field_evolution_fit_numpy(B_initial, t, B_asymptotic, a1, a2, A1, A2, b1, b2, tau_late, a_late)

An analytical function for the magnetic field evolution curves from the magneto-thermal evolution simulations. This method is compatible with NumPy arrays and is used if one wants to save the entire evolution output in the magneto_rotational_evolution method.

The fit parameters for each magneto-thermal model specified in the config_simulator.py file were adjusted by hand (see the notebook tutorials/analysis_notebooks/magnetic_field_evolution_fit.ipynb for more details).

Parameters:

Name Type Description Default
B_initial float

Initial magnetic field strength in [G].

 required
t ndarray

Time in [s].

 required
B_asymptotic float

Asymptotic magnetic field strength at late times in [G].

 required
a1 float

Power-law index for the first power-law component for the magnetic-field evolution fit.

 required
a2 float

Power-law index for the second power-law component for the magnetic-field evolution fit.

 required
A1 float

Normalization for the timescale parameter of the first power-law component for the magnetic-field evolution fit.

 required
A2 float

Normalization for the timescale parameter of the second power-law component for the magnetic-field evolution fit.

 required
b1 float

Power-law index for the timescale parameter of the first power-law component for the magnetic-field evolution fit.

 required
b2 float

Power-law index for the timescale parameter of the second power-law component for the magnetic-field evolution fit.

 required
tau_late float

Timescale in [yr] when transitioning from the simulated curves to the simple late-time power-law evolution of the magnetic field strength.

 required
a_late float

Power-law index for the late-time evolution of the magnetic field strength.

 required

Returns:

Type Description
ndarray

Magnetic field evolution in [G] as a function of time t.
Source code in mlpoppyns/simulator/magneto_rotational_physics/magnetic_field_evolution.py 
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
def magnetic_field_evolution_fit_numpy(
    B_initial: float,
    t: np.ndarray,
    B_asymptotic: float,
    a1: float,
    a2: float,
    A1: float,
    A2: float,
    b1: float,
    b2: float,
    tau_late: float,
    a_late: float,
) -> np.ndarray:
    """
    An analytical function for the magnetic field evolution curves from the magneto-thermal evolution simulations.
    This method is compatible with NumPy arrays and is used if one wants to save the entire evolution output in the
    magneto_rotational_evolution method.

    The fit parameters for each magneto-thermal model specified in the config_simulator.py file were adjusted by hand
    (see the notebook tutorials/analysis_notebooks/magnetic_field_evolution_fit.ipynb for more details).

    Args:
        B_initial (float): Initial magnetic field strength in [G].
        t (np.ndarray): Time in [s].
        B_asymptotic (float): Asymptotic magnetic field strength at late times in [G].
        a1 (float): Power-law index for the first power-law component for the magnetic-field evolution fit.
        a2 (float): Power-law index for the second power-law component for the magnetic-field evolution fit.
        A1 (float): Normalization for the timescale parameter of the first power-law component for the
            magnetic-field evolution fit.
        A2 (float): Normalization for the timescale parameter of the second power-law component for the
            magnetic-field evolution fit.
        b1 (float): Power-law index for the timescale parameter of the first power-law component for the
            magnetic-field evolution fit.
        b2 (float): Power-law index for the timescale parameter of the second power-law component for the
            magnetic-field evolution fit.
        tau_late (float): Timescale in [yr] when transitioning from the simulated curves to the simple late-time
            power-law evolution of the magnetic field strength.
        a_late (float): Power-law index for the late-time evolution of the magnetic field strength.

    Returns:
        (np.ndarray): Magnetic field evolution in [G] as a function of time t.
    """
    # We consider an evolution determined by three main timescales tau1, tau2 and tau_late.
    # At early times the curves are fixed to reproduce the simulated magnetic field evolution from the magneto-thermal
    # code. At late times after a timescale tau_late we assume that the evolution is determined by a simple power-law.

    B = np.zeros(len(t))

    # Define the two timescales as a function of the initial B field.
    tau1 = A1 * B_initial**b1
    tau2 = A2 * B_initial**b2

    if tau2 < tau_late:
        B = (
            B_initial
            * (1 + t / tau1) ** a1
            * (1 + t / tau2) ** (a2 - a1)
            * (1 + t / tau_late) ** (a_late - a2)
        )

    elif (tau1 < tau_late) & (tau_late < tau2):
        B = (
            B_initial
            * (1 + t / tau1) ** a1
            * (1 + t / tau_late) ** (a_late - a1)
        )
    elif tau_late < tau1:
        B = B_initial * (1 + t / tau_late) ** a_late

    # If the magnetic field becomes lower than an asymptotic value derived from the old millisecond pulsar population,
    # fix the magnetic field to that constant asymptotic value unless the initial field is already lower than this
    # asymptotic value.
    if B_initial < B_asymptotic:
        B = np.ones(len(t)) * B_initial
    elif B_initial > B_asymptotic:
        B[B < B_asymptotic] = B_asymptotic

    return B

timescale_Hall(B, L, n_e)

Calculating the Hall timescale for a given field strength, characteristic magnetic field length scale and electron density. Note that for our purposes, we neglect the fact that all quantities can vary significantly with depth inside the neutron star, and we simply use effective quantities that reflect the conservative Hall process. For our choices of L and n_e see the configuration file. B will be identified with the initial dipolar magnetic field components at the pulsars' pole.

Parameters:

Name Type Description Default
B float

Local magnetic field strength, measured in [G].

 required
L float

Characteristic length scale on which the magnetic field varies, measured in [cm].

 required
n_e float

Electron density, measured in [g/cm^3].

 required

Returns:

Type Description
float

Hall timescale in [yr].
Source code in mlpoppyns/simulator/magneto_rotational_physics/magnetic_field_evolution.py 
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
@jit(float64(float64, float64, float64))
def timescale_Hall(B: float, L: float, n_e: float) -> float:
    """
    Calculating the Hall timescale for a given field strength, characteristic magnetic field length
    scale and electron density. Note that for our purposes, we neglect the fact that all quantities can
    vary significantly with depth inside the neutron star, and we simply use effective quantities
    that reflect the conservative Hall process. For our choices of L and n_e see the configuration file.
    B will be identified with the initial dipolar magnetic field components at the pulsars' pole.

    Args:
        B (float): Local magnetic field strength, measured in [G].
        L (float): Characteristic length scale on which the magnetic field varies, measured in [cm].
        n_e (float): Electron density, measured in [g/cm^3].

    Returns:
        (float): Hall timescale in [yr].
    """

    tau_Hall = (
        4 * np.pi * const.E * n_e * L**2 / (const.C * B * const.YR_TO_S)
    )

    return tau_Hall

timescale_ohmic(L, sigma)

Calculating the ohmic diffusion timescale for a given conductivity and characteristic magnetic field length scale. Note that for our purposes, we neglect the fact that both quantities can vary significantly with depth inside the neutron star, and we simply use effective quantities that reflect the ohmic diffusion process. For our choices see the configuration file.

Parameters:

Name Type Description Default
L float

Characteristic length scale on which the magnetic field varies, measured in [cm].

 required
sigma float

Conductivity of the dominating dissipative process, measure in [1/s].

 required

Returns:

Type Description
float

Ohmic diffusion timescale in [yr].
Source code in mlpoppyns/simulator/magneto_rotational_physics/magnetic_field_evolution.py 
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
@jit(float64(float64, float64))
def timescale_ohmic(L: float, sigma: float) -> float:
    """
    Calculating the ohmic diffusion timescale for a given conductivity and characteristic magnetic
    field length scale. Note that for our purposes, we neglect the fact that both quantities can
    vary significantly with depth inside the neutron star, and we simply use effective quantities
    that reflect the ohmic diffusion process. For our choices see the configuration file.

    Args:
        L (float): Characteristic length scale on which the magnetic field varies, measured in [cm].
        sigma (float): Conductivity of the dominating dissipative process, measure in [1/s].

    Returns:
        (float): Ohmic diffusion timescale in [yr].
    """

    tau_ohm = 4 * np.pi * sigma * L**2 / (const.C**2 * const.YR_TO_S)

    return tau_ohm

mlpoppyns.simulator.magneto_rotational_physics.magneto_rotational_evolution

Combined evolution of the pulsar period, misalignment angle and magnetic field relying on an analytical approximation for the magnetic field evolution or fits to magneto-thermal simulations.

Authors:

Vanessa Graber (graber@ice.csic.es)
Michele Ronchi (ronchi@ice.csic.es)

combined_derivatives_analytical(t, y, B_initial, B_asymptotic, L, sigma, n_e, NS_mass, NS_radius)

Combining the two derivative functions for the misalignment angle and the spin period (combined into a single two-component vector y) into a single function to allow combined integration.

This function is applicable for the case, where the magnetic field is prescribed analytically.

Parameters:

Name Type Description Default
t float

Unused time variable, required for the integration below.

 required
y ndarray

Two magneto-rotational parameters, i.e., chi in [rad] and P in [s] for a single pulsar at a given time.

 required
B_initial float

Initial magnetic field magnitude for one pulsar, measured in [G].

 required
B_asymptotic float

Asymptotic magnetic field strength at late times in [G].

 required
L float

Characteristic length scale on which the magnetic field varies, measured in [cm].

 required
sigma float

Conductivity of the dominating dissipative process, measure in [1/s].

 required
n_e float

Electron density, measured in [g/cm^3].

 required
NS_mass float

Neutron star mass, measured in [g].

 required
NS_radius float

Neutron star radius, measured in [cm].

 required

Returns:

Type Description
ndarray

Derivative of the two magneto-rotational parameters for one pulsar,
ndarray

quantities are referred to in respective changes per [yr].
Source code in mlpoppyns/simulator/magneto_rotational_physics/magneto_rotational_evolution.py 
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
@jit(
    [
        float64[:](
            float64,
            float64[:],
            float64,
            float64,
            float64,
            float64,
            float64,
            float64,
            float64,
        )
    ],
    nopython=True,
)
def combined_derivatives_analytical(
    t: float,
    y: np.ndarray,
    B_initial: float,
    B_asymptotic: float,
    L: float,
    sigma: float,
    n_e: float,
    NS_mass: float,
    NS_radius: float,
) -> np.ndarray:
    """
    Combining the two derivative functions for the misalignment angle and the spin period (combined into a single
    two-component vector y) into a single function to allow combined integration.

    This function is applicable for the case, where the magnetic field is prescribed analytically.

    Args:
        t (float): Unused time variable, required for the integration below.
        y (np.ndarray): Two magneto-rotational parameters, i.e., chi in [rad]
            and P in [s] for a single pulsar at a given time.
        B_initial (float): Initial magnetic field magnitude for one pulsar, measured in [G].
        B_asymptotic (float): Asymptotic magnetic field strength at late times in [G].
        L (float): Characteristic length scale on which the magnetic field varies, measured in [cm].
        sigma (float): Conductivity of the dominating dissipative process, measure in [1/s].
        n_e (float): Electron density, measured in [g/cm^3].
        NS_mass (float): Neutron star mass, measured in [g].
        NS_radius (float): Neutron star radius, measured in [cm].

    Returns:
        (np.ndarray): Derivative of the two magneto-rotational parameters for one pulsar,
        quantities are referred to in respective changes per [yr].
    """

    # Unpacking the two components of the vector y.
    chi, P = y

    # Specifying the two derivatives.
    dy = np.zeros(len(y), dtype=np.float64)

    B = mfev.magnetic_field_evolution_analytical(
        B_initial, t, B_asymptotic, L, sigma, n_e
    )

    dy[0] = madv.misalignment_angle_derivative(B, chi, P, NS_mass, NS_radius)
    dy[1] = pdv.period_derivative(B, chi, P, NS_mass, NS_radius)

    return dy

combined_derivatives_fit(t, y, B_initial, B_asymptotic, a1, a2, A1, A2, b1, b2, tau_late, a_late, NS_mass, NS_radius)

Combining the two derivative functions for the misalignment angle and the spin period (combined into a single two-component vector y) into a single function to allow combined integration.

This function is applicable for the case, where the magnetic field is prescribed as a numerical fit.

Parameters:

Name Type Description Default
t float

Unused time variable, required for the integration below.

 required
y ndarray

Two magneto-rotational parameters, i.e., chi in [rad] and P in [s] for a single pulsar at a given time.

 required
B_initial float

Initial magnetic field magnitude for one pulsar, measured in [G].

 required
B_asymptotic float

Asymptotic magnetic field strength at late times in [G].

 required
a1 float

Power-law index for the first power-law component for the magnetic-field evolution fit.

 required
a2 float

Power-law index for the second power-law component for the magnetic-field evolution fit.

 required
A1 float

Normalization for the timescale parameter of the first power-law component for the magnetic-field evolution fit.

 required
A2 float

Normalization for the timescale parameter of the second power-law component for the magnetic-field evolution fit.

 required
b1 float

Power-law index for the timescale parameter of the first power-law component for the magnetic-field evolution fit.

 required
b2 float

Power-law index for the timescale parameter of the second power-law component for the magnetic-field evolution fit.

 required
tau_late float

Timescale in [yr] when transitioning from the simulated curves to the simple late-time power-law evolution of the magnetic field strength.

 required
a_late float

Power-law index for the late-time evolution of the magnetic field strength.

 required
NS_mass float

Neutron star mass, measured in [g].

 required
NS_radius float

Neutron star radius, measured in [cm].

 required

Returns:

Type Description
ndarray

Derivative of the two magneto-rotational parameters for one pulsar,
ndarray

quantities are referred to in respective changes per [yr].
Source code in mlpoppyns/simulator/magneto_rotational_physics/magneto_rotational_evolution.py 
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
@jit(
    [
        float64[:](
            float64,
            float64[:],
            float64,
            float64,
            float64,
            float64,
            float64,
            float64,
            float64,
            float64,
            float64,
            float64,
            float64,
            float64,
        )
    ],
    nopython=True,
)
def combined_derivatives_fit(
    t: float,
    y: np.ndarray,
    B_initial: float,
    B_asymptotic: float,
    a1: float,
    a2: float,
    A1: float,
    A2: float,
    b1: float,
    b2: float,
    tau_late: float,
    a_late: float,
    NS_mass: float,
    NS_radius: float,
) -> np.ndarray:
    """
    Combining the two derivative functions for the misalignment angle and the spin period (combined into a single
    two-component vector y) into a single function to allow combined integration.

    This function is applicable for the case, where the magnetic field is prescribed as a numerical fit.

    Args:
        t (float): Unused time variable, required for the integration below.
        y (np.ndarray): Two magneto-rotational parameters, i.e., chi in [rad]
            and P in [s] for a single pulsar at a given time.
        B_initial (float): Initial magnetic field magnitude for one pulsar, measured in [G].
        B_asymptotic (float): Asymptotic magnetic field strength at late times in [G].
        a1 (float): Power-law index for the first power-law component for the magnetic-field evolution fit.
        a2 (float): Power-law index for the second power-law component for the magnetic-field evolution fit.
        A1 (float): Normalization for the timescale parameter of the first power-law component for the
            magnetic-field evolution fit.
        A2 (float): Normalization for the timescale parameter of the second power-law component for the
            magnetic-field evolution fit.
        b1 (float): Power-law index for the timescale parameter of the first power-law component for the
            magnetic-field evolution fit.
        b2 (float): Power-law index for the timescale parameter of the second power-law component for the
            magnetic-field evolution fit.
        tau_late (float): Timescale in [yr] when transitioning from the simulated curves to the simple late-time
            power-law evolution of the magnetic field strength.
        a_late (float): Power-law index for the late-time evolution of the magnetic field strength.
        NS_mass (float): Neutron star mass, measured in [g].
        NS_radius (float): Neutron star radius, measured in [cm].

    Returns:
        (np.ndarray): Derivative of the two magneto-rotational parameters for one pulsar,
        quantities are referred to in respective changes per [yr].
    """

    # Unpacking the two components of the vector y.
    chi, P = y

    # Specifying the two derivatives.
    dy = np.zeros(len(y), dtype=np.float64)

    B = mfev.magnetic_field_evolution_fit(
        B_initial, t, B_asymptotic, a1, a2, A1, A2, b1, b2, tau_late, a_late
    )

    dy[0] = madv.misalignment_angle_derivative(B, chi, P, NS_mass, NS_radius)
    dy[1] = pdv.period_derivative(B, chi, P, NS_mass, NS_radius)

    return dy

evolve_population_magrot(dict_pop_initial_magrot, output_path)

Evolve the magneto-rotational properties of a neutron star population over time based on initial conditions.

Parameters:

Name Type Description Default
dict_pop_initial_magrot dict

Dictionary containing initial magneto-rotational properties of the population.

 required
output_path Path

The path where the evolution data will be saved if enabled in the configuration.

 required

Returns:

Type Description
dict

A dictionary containing the properties of the evolved neutron star population.
Source code in mlpoppyns/simulator/magneto_rotational_physics/magneto_rotational_evolution.py 
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
def evolve_population_magrot(
    dict_pop_initial_magrot: dict,
    output_path: pathlib.Path,
) -> dict:
    """
    Evolve the magneto-rotational properties of a neutron star population over time based on initial conditions.

    Args:
        dict_pop_initial_magrot (dict): Dictionary containing initial magneto-rotational properties of the population.
        output_path (pathlib.Path): The path where the evolution data will be saved if enabled in the configuration.

    Returns:
        (dict): A dictionary containing the properties of the evolved neutron star population.
    """
    age = dict_pop_initial_magrot["age"]
    B_initial = dict_pop_initial_magrot["B"]
    chi_initial = dict_pop_initial_magrot["chi"]
    P_initial = dict_pop_initial_magrot["P"]

    # Determine the evolved magnetic field, misalignment angle and rotation period.
    (
        B_final,
        chi_final,
        P_final,
        magrot_evol_dict,
    ) = magneto_rotational_evolution(
        B_initial,
        chi_initial,
        P_initial,
        age,
    )

    if cfg["save_magrot_evolution"]:
        # Save dictionary containing evolution information to output path in a .json file.
        magrot_evolution_dump_path = pathlib.Path().joinpath(
            output_path, "magrot_evolution.json"
        )

        with open(magrot_evolution_dump_path, "wb") as f:
            f.write(
                orjson.dumps(
                    dict(magrot_evol_dict),
                    option=orjson.OPT_SERIALIZE_NUMPY
                    | orjson.OPT_NON_STR_KEYS
                    | orjson.OPT_SORT_KEYS,
                )
            )

    # Determining the final period derivative.
    P_dot_final = pdv.period_derivative_numpy(
        B_final,
        chi_final,
        P_final,
        cfg["NS_mass"],
        cfg["NS_radius"],
    )

    dictionary_final_pop_magrot = {
        "B_initial": B_initial,
        "B": B_final,
        "chi": chi_final,
        "P": P_final,
        "P_dot": P_dot_final,
    }

    return dictionary_final_pop_magrot

initialize_population_magrot(age)

Initialize the magneto-rotational properties of a neutron star population.

Parameters:

Name Type Description Default
age ndarray

An array of ages in [yr] fot the neutron stars that has to be initialized for the magneto-rotational evolution.

 required

Returns:

Type Description
dict

A dictionary containing the initialized magneto-rotational properties of the neutron star population in the survey sky coverage.
Source code in mlpoppyns/simulator/magneto_rotational_physics/magneto_rotational_evolution.py 
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
def initialize_population_magrot(age: np.ndarray) -> dict:
    """
    Initialize the magneto-rotational properties of a neutron star population.

    Args:
        age (np.ndarray): An array of ages in [yr] fot the neutron stars that has to be initialized for the
            magneto-rotational evolution.

    Returns:
        (dict): A dictionary containing the initialized magneto-rotational properties of the neutron star population
            in the survey sky coverage.
    """

    # Initialize neutron star population properties.
    pop_initial = ipop.InitialNeutronStarPopulation(NS_number=len(age))

    # Computing the initial field strengths, misalignment angles, and periods.
    B_initial = pop_initial.magnetic_field()
    chi_initial = pop_initial.misalignment_angle()
    P_initial = pop_initial.period()
    P_dot_initial = pdv.period_derivative_numpy(
        B_initial, chi_initial, P_initial, cfg["NS_mass"], cfg["NS_radius"]
    )

    dictionary_initial_pop_magrot = {
        "age": age,
        "B": B_initial,
        "chi": chi_initial,
        "P": P_initial,
        "P_dot": P_dot_initial,
    }

    return dictionary_initial_pop_magrot

magneto_rotational_evolution(B_initial, chi_initial, P_initial, t_age)

Evolving the neutron stars' magnetic fields, misalignment angles and periods according to their respective ages forward in time to obtain their current magnetic field strengths, misalignment angles and periods. Note that right now the times at which these three parameters are evaluated (apart from the current time) do not agree for pulsars.

Parameters:

Name Type Description Default
B_initial ndarray

Pulsars' initial magnetic field magnitudes, measured in [G].

 required
chi_initial ndarray

Pulsars' initial misalignment angles, measured in [rad].

 required
P_initial ndarray

Pulsars' initial rotation periods, measured in [s].

 required
t_age ndarray

Array of neutron star ages in [yr].

 required

Returns:

Type Description
Tuple[ndarray, ndarray, ndarray, dict]

Tuple consisting of three arrays defining the neutron stars' final magnetic field strengths in [G], misalignment angles in [rad] and rotation periods in [s] and a dictionary containing the time evolution of these quantities for each neutron star (if the option to save the time evolution is enabled).
Source code in mlpoppyns/simulator/magneto_rotational_physics/magneto_rotational_evolution.py 
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
def magneto_rotational_evolution(
    B_initial: np.ndarray,
    chi_initial: np.ndarray,
    P_initial: np.ndarray,
    t_age: np.ndarray,
) -> Tuple[np.ndarray, np.ndarray, np.ndarray, dict]:
    """
    Evolving the neutron stars' magnetic fields, misalignment angles and periods according to their respective ages
    forward in time to obtain their current magnetic field strengths, misalignment angles and periods. Note that right
    now the times at which these three parameters are evaluated (apart from the current time) do not agree for pulsars.

    Args:
        B_initial (np.ndarray): Pulsars' initial magnetic field magnitudes, measured in [G].
        chi_initial (np.ndarray): Pulsars' initial misalignment angles, measured in [rad].
        P_initial (np.ndarray): Pulsars' initial rotation periods, measured in [s].
        t_age (np.ndarray): Array of neutron star ages in [yr].

    Returns:
        (Tuple[np.ndarray, np.ndarray, np.ndarray, dict]): Tuple consisting of three arrays defining the neutron stars'
            final magnetic field strengths in [G], misalignment angles in [rad] and rotation periods in [s] and a
            dictionary containing the time evolution of these quantities for each neutron star (if the option to save
            the time evolution is enabled).
    """

    # Save the number of simulated neutron stars, which is flexible depending on whether
    # they are simulated all at once or one by one.
    n = len(t_age)

    # Initialization of a dictionary that will contain the evolution in time of B, chi and P.
    evolution_dictionary = {}

    # Initialization of the array for the three parameters.
    B_final = np.zeros(n)
    chi_final = np.zeros(n)
    P_final = np.zeros(n)

    # Initial conditions for the two parameters.
    y_initial = np.column_stack((chi_initial, P_initial))

    # Draw a random asymptotic value of the magnetic field at late times from a log-normal distribution.
    # This asymptotic value is based on the distribution of inferred magnetic fields for the old population
    # of millisecond pulsars.
    B_asymptotic = 10 ** np.random.normal(
        cfg["B_millisec_mean"], cfg["B_millisec_sigma"], n
    )

    for i in range(n):
        # Generating a time grid at which the solution is evaluated. We start to
        # evolve each star at its birth, corresponding to time 0, and do so for
        # its full age in steps of the time_step specified in the configuration
        # file. To obtain the magnetic field at the current time, we append the
        # current age value.
        time_grid = np.append(
            10
            ** np.arange(0, np.log10(t_age[i]), cfg["magrot_time_step_log10"]),
            t_age[i],
        )

        # To integrate the problem, we use scipy's odeint function.
        # We set tfirst=True to unify the structure of the input ODEs in order to be able
        # to compare different scipy functions to solve the ODEs.
        if cfg["magneto-thermal_model"] == "analytical":
            evol_output = np.array(
                odeint(
                    combined_derivatives_analytical,
                    y0=y_initial[i],
                    t=time_grid,
                    args=(
                        B_initial[i],
                        B_asymptotic[i],
                        cfg["L"],
                        cfg["sigma"],
                        cfg["n_e"],
                        cfg["NS_mass"],
                        cfg["NS_radius"],
                    ),
                    tfirst=True,
                )
            )
        else:
            evol_output = np.array(
                odeint(
                    combined_derivatives_fit,
                    y0=y_initial[i],
                    t=time_grid,
                    args=(
                        B_initial[i],
                        B_asymptotic[i],
                        cfg["a1"],
                        cfg["a2"],
                        cfg["A1"],
                        cfg["A2"],
                        cfg["b1"],
                        cfg["b2"],
                        cfg["tau_late"],
                        cfg["a_late"],
                        cfg["NS_mass"],
                        cfg["NS_radius"],
                    ),
                    tfirst=True,
                )
            )

        if cfg["save_magrot_evolution"]:
            # Evaluate the magnetic field evolution.
            if cfg["magneto-thermal_model"] == "analytical":
                B_t = mfev.magnetic_field_evolution_analytical_numpy(
                    B_initial[i],
                    time_grid,
                    B_asymptotic[i],
                    cfg["L"],
                    cfg["sigma"],
                    cfg["n_e"],
                )
            else:
                B_t = mfev.magnetic_field_evolution_fit_numpy(
                    B_initial[i],
                    time_grid,
                    B_asymptotic[i],
                    cfg["a1"],
                    cfg["a2"],
                    cfg["A1"],
                    cfg["A2"],
                    cfg["b1"],
                    cfg["b2"],
                    cfg["tau_late"],
                    cfg["a_late"],
                )

            # Save the evolution output of the i-th neutron star in a dictionary.
            evolution = {
                i: {
                    "t": time_grid.tolist(),
                    "B(t)": B_t.tolist(),
                    "chi(t)": evol_output[:, 0].tolist(),
                    "P(t)": evol_output[:, 1].tolist(),
                }
            }
            # Update the dictionary containing the evolution information.
            evolution_dictionary = {**evolution_dictionary, **evolution}

        # Save the final values of the magnetic field, inclination angle and spin period.
        if cfg["magneto-thermal_model"] == "analytical":
            B_final[i] = mfev.magnetic_field_evolution_analytical(
                B_initial[i],
                t_age[i],
                B_asymptotic[i],
                cfg["L"],
                cfg["sigma"],
                cfg["n_e"],
            )
        else:
            B_final[i] = mfev.magnetic_field_evolution_fit(
                B_initial[i],
                t_age[i],
                B_asymptotic[i],
                cfg["a1"],
                cfg["a2"],
                cfg["A1"],
                cfg["A2"],
                cfg["b1"],
                cfg["b2"],
                cfg["tau_late"],
                cfg["a_late"],
            )

        chi_final[i] = evol_output[-1, 0]
        P_final[i] = evol_output[-1, 1]

    return B_final, chi_final, P_final, evolution_dictionary

mlpoppyns.simulator.magneto_rotational_physics.misalignment_angle_derivative

Time derivative of the pulsar misalignment angle.

Authors:

Vanessa Graber (graber@ice.csic.es)
Michele Ronchi (ronchi@ice.csic.es)

misalignment_angle_derivative(B, chi, P, NS_mass, NS_radius)

This function determines the change in the misalignment angle, i.e., the angle between the magnetic dipolar moment and the rotation axis of a pulsar. It is taken from eq. (71) of Pons & Vigano (2019). For more details see, e.g., Spitkovsky (2006) or Philippov et al. (2014), who determine the coefficients k_0, k_1, k_2 (defined in the configuration file) for a pulsar embedded in a force-free and resistive magnetosphere from numerical simulations. Note that all three input parameters are time-dependent.

Parameters:

Name Type Description Default
B float

Values of the dipolar component of the magnetic field at the magnetic pole for the sample of simulated neutron stars, measured in [G].

 required
chi float

Angles between the magnetic dipolar moment, i.e., the magnetic field axis, and the rotation axis for all simulated pulsars, measured in [rad].

 required
P float

Spin periods of simulated pulsars, measured in [s].

 required
NS_mass float

Neutron star mass, measured in [g].

 required
NS_radius float

Neutron star radius, measured in [cm].

 required

Returns:

Type Description
float

Misalignment angle derivatives for all simulated pulsars in [rad/yr].
Source code in mlpoppyns/simulator/magneto_rotational_physics/misalignment_angle_derivative.py 
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
@jit(float64(float64, float64, float64, float64, float64))
def misalignment_angle_derivative(
    B: float, chi: float, P: float, NS_mass: float, NS_radius: float
) -> float:
    """
    This function determines the change in the misalignment angle, i.e., the angle between the
    magnetic dipolar moment and the rotation axis of a pulsar. It is taken from eq. (71) of
    Pons & Vigano (2019). For more details see, e.g., Spitkovsky (2006) or Philippov et al.
    (2014), who determine the coefficients k_0, k_1, k_2 (defined in the configuration file)
    for a pulsar embedded in a force-free and resistive magnetosphere from numerical simulations.
    Note that all three input parameters are time-dependent.

    Args:
        B (float): Values of the dipolar component of the magnetic field at the
            magnetic pole for the sample of simulated neutron stars, measured in [G].
        chi (float): Angles between the magnetic dipolar moment, i.e., the magnetic
            field axis, and the rotation axis for all simulated pulsars, measured in [rad].
        P (float): Spin periods of simulated pulsars, measured in [s].
        NS_mass (float): Neutron star mass, measured in [g].
        NS_radius (float): Neutron star radius, measured in [cm].

    Returns:
        (float): Misalignment angle derivatives for all simulated pulsars in [rad/yr].
    """

    # Canonical neutron star moment of inertia in [g cm^2] assuming a perfect solid sphere.
    NS_inertia = 2.0 / 5.0 * NS_mass * NS_radius**2

    # Auxiliary quantity beta as defined in eq. (72) of Pons & Vigano (2019).
    beta = np.pi**2 * NS_radius**6 / (NS_inertia * const.C**3)

    # Misalignment angle derivative.
    chi_deriv = (
        -k_coefficients_2
        * beta
        * B**2
        / (P**2)
        * np.sin(chi)
        * np.cos(chi)
    ) * const.YR_TO_S

    return chi_deriv

mlpoppyns.simulator.magneto_rotational_physics.period_derivative

Time derivative of the pulsar period.

Authors:

Vanessa Graber (graber@ice.csic.es)
Michele Ronchi (ronchi@ice.csic.es)

period_derivative(B, chi, P, NS_mass, NS_radius)

This function determines the change in the rotation period of a pulsar. It is taken from eq. (70) of Pons & Vigano (2019). For more details see, e.g., Spitkovsky (2006) or Philippov et al. (2014), who determine the coefficients k_0, k_1, k_2 (defined in the configuration file) for a pulsar embedded in a force-free and resistive magnetosphere from numerical simulations. Note that all three input parameters are time-dependent.

Parameters:

Name Type Description Default
B float

Value of the dipolar component of the magnetic field at the magnetic pole for a simulated neutron star, measured in [G].

 required
chi float

Angle between the magnetic dipolar moment, i.e., the magnetic field axis, and the rotation axis for a simulated pulsar, measured in [rad].

 required
P float

Spin period of a simulated pulsar, measured in [s].

 required
NS_mass float

Neutron star mass, measured in [g].

 required
NS_radius float

Neutron star radius, measured in [cm].

 required

Returns:

Type Description
float

Period derivative of a simulated pulsar in [s/yr].
Source code in mlpoppyns/simulator/magneto_rotational_physics/period_derivative.py 
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
@jit(float64(float64, float64, float64, float64, float64))
def period_derivative(
    B: float, chi: float, P: float, NS_mass: float, NS_radius: float
) -> float:
    """
    This function determines the change in the rotation period of a pulsar. It is taken
    from eq. (70) of Pons & Vigano (2019). For more details see, e.g., Spitkovsky (2006)
    or Philippov et al. (2014), who determine the coefficients k_0, k_1, k_2 (defined in
    the configuration file) for a pulsar embedded in a force-free and resistive magnetosphere
    from numerical simulations. Note that all three input parameters are time-dependent.

    Args:
        B (float): Value of the dipolar component of the magnetic field at the
            magnetic pole for a simulated neutron star, measured in [G].
        chi (float): Angle between the magnetic dipolar moment, i.e., the magnetic
            field axis, and the rotation axis for a simulated pulsar, measured in [rad].
        P (float): Spin period of a simulated pulsar, measured in [s].
        NS_mass (float): Neutron star mass, measured in [g].
        NS_radius (float): Neutron star radius, measured in [cm].

    Returns:
        (float): Period derivative of a simulated pulsar in [s/yr].
    """

    # Canonical neutron star moment of inertia in [g cm^2] assuming a perfect solid sphere.
    NS_inertia = 2.0 / 5.0 * NS_mass * NS_radius**2

    # Auxiliary quantity beta as defined in eq. (72) of Pons & Vigano (2019).
    beta = np.pi**2 * NS_radius**6 / (NS_inertia * const.C**3)

    # Period derivative.
    P_deriv = (
        beta
        * B**2
        / P
        * (k_coefficients_0 + k_coefficients_1 * np.sin(chi) ** 2)
    ) * const.YR_TO_S

    return P_deriv

period_derivative_numpy(B, chi, P, NS_mass, NS_radius)

This is the numpy array version of the function above. It determines the change in the rotation period of a pulsar. It is taken from eq. (70) of Pons & Vigano (2019). For more details see, e.g., Spitkovsky (2006) or Philippov et al. (2014), who determine the coefficients k_0, k_1, k_2 (defined in the configuration file) for a pulsar embedded in a force-free and resistive magnetosphere from numerical simulations. Note that all three input parameters are time-dependent.

Parameters:

Name Type Description Default
B ndarray

Array of values of the dipolar component of the magnetic field at the magnetic pole for the simulated neutron stars, measured in [G].

 required
chi ndarray

Array of angles between the magnetic dipolar moment, i.e., the magnetic field axis, and the rotation axis for the simulated pulsars, measured in [rad].

 required
P ndarray

Array of spin period of a simulated pulsars, measured in [s].

 required
NS_mass float

Neutron star mass, measured in [g].

 required
NS_radius float

Neutron star radius, measured in [cm].

 required

Returns:

Type Description
ndarray

Array of period derivatives of the simulated pulsars in [s/s].
Source code in mlpoppyns/simulator/magneto_rotational_physics/period_derivative.py 
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
def period_derivative_numpy(
    B: np.ndarray,
    chi: np.ndarray,
    P: np.ndarray,
    NS_mass: float,
    NS_radius: float,
) -> np.ndarray:
    """
    This is the numpy array version of the function above. It determines the change in the
    rotation period of a pulsar. It is taken from eq. (70) of Pons & Vigano (2019).
    For more details see, e.g., Spitkovsky (2006) or Philippov et al. (2014), who determine
    the coefficients k_0, k_1, k_2 (defined in the configuration file) for a pulsar embedded
    in a force-free and resistive magnetosphere from numerical simulations.
    Note that all three input parameters are time-dependent.

    Args:
        B (np.ndarray): Array of values of the dipolar component of the magnetic field at the
            magnetic pole for the simulated neutron stars, measured in [G].
        chi (np.ndarray): Array of angles between the magnetic dipolar moment, i.e., the magnetic
            field axis, and the rotation axis for the simulated pulsars, measured in [rad].
        P (np.ndarray): Array of spin period of a simulated pulsars, measured in [s].
        NS_mass (float): Neutron star mass, measured in [g].
        NS_radius (float): Neutron star radius, measured in [cm].

    Returns:
        (np.ndarray): Array of period derivatives of the simulated pulsars in [s/s].
    """

    # Canonical neutron star moment of inertia in [g cm^2] assuming a perfect solid sphere.
    NS_inertia = 2.0 / 5.0 * NS_mass * NS_radius**2

    # Auxiliary quantity beta as defined in eq. (72) of Pons & Vigano (2019).
    beta = np.pi**2 * NS_radius**6 / (NS_inertia * const.C**3)

    # Period derivative.
    P_deriv = (
        beta
        * B**2
        / P
        * (k_coefficients_0 + k_coefficients_1 * np.sin(chi) ** 2)
    )

    return P_deriv