Add spins to BBH ID input file

support/Pipelines/Bbh/InitialData.py

@@ -3,9 +3,10 @@
 
 import logging
 from pathlib import Path
-from typing import Optional, Union
+from typing import Optional, Sequence, Union
 
 import click
+import numpy as np
 from rich.pretty import pretty_repr
 
 from spectre.support.Schedule import schedule, scheduler_options
@@ -17,6 +18,8 @@
 
 def id_parameters(
     mass_ratio: float,
+    dimensionless_spin_a: Sequence[float],
+    dimensionless_spin_b: Sequence[float],
     separation: float,
     orbital_angular_velocity: float,
     radial_expansion_velocity: float,
@@ -29,6 +32,8 @@ def id_parameters(
 
     Arguments:
       mass_ratio: Defined as q = M_A / M_B >= 1.
+      dimensionless_spin_a: Dimensionless spin of the larger black hole, chi_A.
+      dimensionless_spin_b: Dimensionless spin of the smaller black hole, chi_B.
       separation: Coordinate separation D of the black holes.
       orbital_angular_velocity: Omega_0.
       radial_expansion_velocity: adot_0.
@@ -44,6 +49,14 @@ def id_parameters(
     M_B = 1.0 / (1.0 + mass_ratio)
     x_A = separation / (1.0 + mass_ratio)
     x_B = x_A - separation
+    chi_A = np.asarray(dimensionless_spin_a)
+    r_plus_A = M_A * (1.0 + np.sqrt(1 - np.dot(chi_A, chi_A)))
+    Omega_A = -0.5 * chi_A / r_plus_A
+    Omega_A[2] += orbital_angular_velocity
+    chi_B = np.asarray(dimensionless_spin_b)
+    r_plus_B = M_B * (1.0 + np.sqrt(1 - np.dot(chi_B, chi_B)))
+    Omega_B = -0.5 * chi_B / r_plus_B
+    Omega_B[2] += orbital_angular_velocity
     return {
         "MassRight": M_A,
         "MassLeft": M_B,
@@ -53,6 +66,18 @@ def id_parameters(
         "ExcisionRadiusLeft": 0.89 * 2.0 * M_B,
         "OrbitalAngularVelocity": orbital_angular_velocity,
         "RadialExpansionVelocity": radial_expansion_velocity,
+        "DimensionlessSpinRight_x": chi_A[0],
+        "DimensionlessSpinRight_y": chi_A[1],
+        "DimensionlessSpinRight_z": chi_A[2],
+        "DimensionlessSpinLeft_x": chi_B[0],
+        "DimensionlessSpinLeft_y": chi_B[1],
+        "DimensionlessSpinLeft_z": chi_B[2],
+        "HorizonRotationRight_x": Omega_A[0],
+        "HorizonRotationRight_y": Omega_A[1],
+        "HorizonRotationRight_z": Omega_A[2],
+        "HorizonRotationLeft_x": Omega_B[0],
+        "HorizonRotationLeft_y": Omega_B[1],
+        "HorizonRotationLeft_z": Omega_B[2],
         # Resolution
         "L": refinement_level,
         "P": polynomial_order,
@@ -61,6 +86,8 @@ def id_parameters(
 
 def generate_id(
     mass_ratio: float,
+    dimensionless_spin_a: Sequence[float],
+    dimensionless_spin_b: Sequence[float],
     separation: float,
     orbital_angular_velocity: float,
     radial_expansion_velocity: float,
@@ -109,6 +136,8 @@ def generate_id(
     # Determine initial data parameters from options
     id_params = id_parameters(
         mass_ratio=mass_ratio,
+        dimensionless_spin_a=dimensionless_spin_a,
+        dimensionless_spin_b=dimensionless_spin_b,
         separation=separation,
         orbital_angular_velocity=orbital_angular_velocity,
         radial_expansion_velocity=radial_expansion_velocity,
@@ -137,6 +166,22 @@ def generate_id(
     help="Mass ratio of the binary, defined as q = M_A / M_B >= 1.",
     required=True,
 )
+@click.option(
+    "--dimensionless-spin-A",
+    "--chi-A",
+    type=float,
+    nargs=3,
+    help="Dimensionless spin of the larger black hole, chi_A.",
+    required=True,
+)
+@click.option(
+    "--dimensionless-spin-B",
+    "--chi-B",
+    type=float,
+    nargs=3,
+    help="Dimensionless spin of the smaller black hole, chi_B.",
+    required=True,
+)
 @click.option(
     "--separation",
     "-D",

support/Pipelines/Bbh/InitialData.yaml

@@ -29,12 +29,18 @@ Background: &background
     ObjectLeft: &kerr_left
       KerrSchild:
         Mass: {{ MassLeft }}
-        Spin: [0., 0., 0.]
+        Spin:
+          - {{ DimensionlessSpinLeft_x }}
+          - {{ DimensionlessSpinLeft_y }}
+          - {{ DimensionlessSpinLeft_z }}
         Center: [0., 0., 0.]
     ObjectRight: &kerr_right
       KerrSchild:
         Mass: {{ MassRight }}
-        Spin: [0., 0., 0.]
+        Spin:
+          - {{ DimensionlessSpinRight_x }}
+          - {{ DimensionlessSpinRight_y }}
+          - {{ DimensionlessSpinRight_z }}
         Center: [0., 0., 0.]
     AngularVelocity: {{ OrbitalAngularVelocity }}
     Expansion: {{ RadialExpansionVelocity }}
@@ -53,7 +59,10 @@ DomainCreator:
         ExciseWithBoundaryCondition:
           ApparentHorizon:
             Center: [*x_right, 0., 0.]
-            Rotation: [0., 0., 0.]
+            Rotation:
+              - {{ HorizonRotationRight_x }}
+              - {{ HorizonRotationRight_y }}
+              - {{ HorizonRotationRight_z }}
             Lapse: *kerr_right
             NegativeExpansion: *kerr_right
       UseLogarithmicMap: True
@@ -65,7 +74,10 @@ DomainCreator:
         ExciseWithBoundaryCondition:
           ApparentHorizon:
             Center: [*x_left, 0., 0.]
-            Rotation: [0., 0., 0.]
+            Rotation:
+              - {{ HorizonRotationLeft_x }}
+              - {{ HorizonRotationLeft_y }}
+              - {{ HorizonRotationLeft_z }}
             Lapse: *kerr_left
             NegativeExpansion: *kerr_left
       UseLogarithmicMap: True

tests/support/Pipelines/Bbh/Test_InitialData.py

@@ -6,6 +6,7 @@
 import unittest
 from pathlib import Path
 
+import numpy.testing as npt
 import yaml
 from click.testing import CliRunner
 
@@ -28,6 +29,8 @@ def tearDown(self):
     def test_generate_id(self):
         params = id_parameters(
             mass_ratio=1.5,
+            dimensionless_spin_a=[0.1, 0.2, 0.3],
+            dimensionless_spin_b=[0.4, 0.5, 0.6],
             separation=20.0,
             orbital_angular_velocity=0.01,
             radial_expansion_velocity=-1.0e-5,
@@ -42,6 +45,22 @@ def test_generate_id(self):
         self.assertAlmostEqual(params["ExcisionRadiusLeft"], 0.712)
         self.assertEqual(params["OrbitalAngularVelocity"], 0.01)
         self.assertEqual(params["RadialExpansionVelocity"], -1.0e-5)
+        self.assertEqual(
+            [params[f"DimensionlessSpinRight_{xyz}"] for xyz in "xyz"],
+            [0.1, 0.2, 0.3],
+        )
+        self.assertEqual(
+            [params[f"DimensionlessSpinLeft_{xyz}"] for xyz in "xyz"],
+            [0.4, 0.5, 0.6],
+        )
+        npt.assert_allclose(
+            [params[f"HorizonRotationRight_{xyz}"] for xyz in "xyz"],
+            [-0.043236994315732, -0.086473988631464, -0.119710982947196],
+        )
+        npt.assert_allclose(
+            [params[f"HorizonRotationLeft_{xyz}"] for xyz in "xyz"],
+            [-0.337933017966707, -0.422416272458383, -0.49689952695006],
+        )
         self.assertEqual(params["L"], 1)
         self.assertEqual(params["P"], 5)
         # COM is zero
@@ -55,6 +74,14 @@ def test_cli(self):
         common_args = [
             "--mass-ratio",
             "1.5",
+            "--chi-A",
+            "0.1",
+            "0.2",
+            "0.3",
+            "--chi-B",
+            "0.4",
+            "0.5",
+            "0.6",
             "--separation",
             "20",
             "--orbital-angular-velocity",

tests/support/Pipelines/Bbh/Test_Inspiral.py

@@ -27,6 +27,8 @@ def setUp(self):
         self.bin_dir = Path(unit_test_build_path(), "../../bin").resolve()
         generate_id(
             mass_ratio=1.5,
+            dimensionless_spin_a=[0.0, 0.0, 0.0],
+            dimensionless_spin_b=[0.0, 0.0, 0.0],
             separation=20.0,
             orbital_angular_velocity=0.01,
             radial_expansion_velocity=-1.0e-5,

tests/support/Pipelines/Bbh/Test_Ringdown.py

@@ -28,6 +28,8 @@ def setUp(self):
         self.bin_dir = Path(unit_test_build_path(), "../../bin").resolve()
         generate_id(
             mass_ratio=1.5,
+            dimensionless_spin_a=[0.0, 0.0, 0.0],
+            dimensionless_spin_b=[0.0, 0.0, 0.0],
             separation=20.0,
             orbital_angular_velocity=0.01,
             radial_expansion_velocity=-1.0e-5,
@@ -53,7 +55,7 @@ def setUp(self):
     def tearDown(self):
         shutil.rmtree(self.test_dir, ignore_errors=True)
 
-    def test_id_parameters(self):
+    def test_ringdown_parameters(self):
         with open(self.inspiral_dir / "Inspiral.yaml") as open_input_file:
             _, inspiral_input_file = yaml.safe_load_all(open_input_file)
         params = ringdown_parameters(