From 9a7e0e15f69ce7e2aa2b7ff9df998b2b5b937edd Mon Sep 17 00:00:00 2001
From: IndigoCarmine <tamanisikaminai@gmail.com>
Date: Tue, 17 Mar 2026 04:49:42 +0900
Subject: [PATCH 1/9] Add tetrapod model implementation in C and Python

---
 sasmodels/models/tetrapod.c  | 60 +++++++++++++++++++++++++++++
 sasmodels/models/tetrapod.py | 74 ++++++++++++++++++++++++++++++++++++
 2 files changed, 134 insertions(+)
 create mode 100644 sasmodels/models/tetrapod.c
 create mode 100644 sasmodels/models/tetrapod.py

diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c
new file mode 100644
index 00000000..622e4b83
--- /dev/null
+++ b/sasmodels/models/tetrapod.c
@@ -0,0 +1,60 @@
+const double A =
+    109.5 / 2.0 * M_PI / 180.0;  // half of the angle between arms in radians
+
+static double u_n(int n, double theta, double alpha) {
+  const double phi[4] = {0.0, M_PI_2, M_PI, 3.0 * M_PI_2};
+  return cos(A) * cos(theta + phi[n] * 2) +
+         sin(A) * sin(theta + phi[n] * 2) * cos(alpha - phi[n]);
+}
+
+static double A_n(double q, double u, double L, double R) {
+  double quL2 = q * u * L * 0.5;
+
+  double term1 = (fabs(quL2) > 1e-12) ? sin(quL2) / quL2 : 1.0;
+
+  double mu = sqrt(fmax(0.0, 1.0 - u * u));
+  double qmuR = q * mu * R;
+
+  double term2 = (fabs(qmuR) > 1e-12) ? 2.0 * sas_J1(qmuR) / qmuR : 1.0;
+
+  return term1 * term2;
+}
+
+static double form_volume(double L, double R) {
+  // V = 4 * \pi * R^2 * L
+  return 4.0 * M_PI * R * R * L;
+}
+
+static double Iq(double q, double L, double R, double sld_particle,
+                 double sld_solvent) {
+  double contrast = sld_particle - sld_solvent;
+  double total = 0.0;
+
+  for (int dtheta = 0; dtheta < GAUSS_N; dtheta++) {
+    double theta =
+        M_PI_2 * (GAUSS_Z[dtheta] + 1.0);  // map from [-1, 1] to [0, pi]
+    double w_theta =
+        GAUSS_W[dtheta] * M_PI_2;  // adjust weight for the new range
+
+    double integral_alpha = 0.0;
+    for (int dalpha = 0; dalpha < GAUSS_N; dalpha++) {
+      double alpha =
+          M_PI * (GAUSS_Z[dalpha] + 1.0);  // map from [-1, 1] to [0, pi]
+      double w_alpha =
+          GAUSS_W[dalpha] * M_PI;  // adjust weight for the new range
+
+      double sum_arms = 0.0;
+      for (int n = 0; n < 4; n++) {
+        for (int m = 0; m < 4; m++) {
+          double u = u_n(n, theta, alpha);
+          double v = u_n(m, theta, alpha);
+          sum_arms += A_n(q, u, L, R) * A_n(q, v, L, R) *
+                      cos(q * (u - v) * L / 2.0) * sin(theta);
+        }
+      }
+      integral_alpha += sum_arms * w_alpha;
+    }
+    total += integral_alpha * w_theta;
+  }
+  return contrast * contrast * total * 2 * M_PI;
+}
\ No newline at end of file
diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py
new file mode 100644
index 00000000..d0c51b94
--- /dev/null
+++ b/sasmodels/models/tetrapod.py
@@ -0,0 +1,74 @@
+r"""
+Definition
+----------
+
+Calculates the scattering from a tetrapod-shaped structure. A tetrapod consists
+of four cylindrical arms radiating from a central point, oriented along the
+(1,1,1), (-1,-1,1), (-1,1,-1), and (1,-1,-1) directions.
+
+The scattering intensity is calculated as:
+
+.. math::
+
+    I(q) = (\Delta \rho)^2 \left[ 4 I_{\text{arm}}(q) + 6 I_{\text{corr}}(q) \right]
+
+where $\Delta\rho$ is the scattering length density contrast between the
+tetrapod and the solvent.
+
+The arm scattering $I_{\text{arm}}(q)$ is obtained by integrating over all
+orientations of each cylindrical arm:
+
+.. math::
+
+    I_{\text{arm}}(q) = \int_0^{\pi/4} F_{\text{cylinder}}(q, \theta, L, R)^2 \sin(\theta) d\theta
+
+where $F_{\text{cylinder}}$ is the form factor of a cylinder.
+
+The correlation term $I_{\text{corr}}(q)$ accounts for coherent scattering
+between different arms.
+
+Geometry
+--------
+
+The four arms are positioned along the following directions:
+- Arm 1: (1, 1, 1)
+- Arm 2: (-1, -1, 1)
+- Arm 3: (-1, 1, -1)
+- Arm 4: (1, -1, -1)
+
+Each arm has length $L$ and radius $R$.
+
+References
+----------
+
+#. J. S. Pedersen, *J. Appl. Cryst.*, 33 (2000) 488-494
+
+Authorship and Verification
+----------------------------
+
+* **Author:** NIST IGOR/DANSE **Date:** pre 2010
+* **Last Modified by:** Paul Butler **Date:** March 20, 2016
+"""
+
+import numpy as np
+from numpy import cos, inf, pi, sin
+
+name = "tetrapod"
+title = "Tetrapod with four cylindrical arms"
+description = """
+    Calculates the scattering from a tetrapod structure with four cylindrical
+    arms radiating from a central point.
+"""
+category = "shape:cylinder"
+
+#             [ "name", "units", default, [lower, upper], "type", "description"],
+parameters = [
+    ["length", "Ang", 400, [0, inf], "volume", "Cylindrical arm length"],
+    ["radius", "Ang", 30, [0, inf], "volume", "Cylindrical arm radius"],
+    ["sld", "1e-6/Ang^2", 4, [-inf, inf], "sld", "Tetrapod scattering length density"],
+    ["sld_solvent", "1e-6/Ang^2", 1, [-inf, inf], "sld", "Solvent scattering length density"],
+]
+
+source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "tetrapod.c"]
+valid = "radius >= 0.0 && length >= 0.0"
+have_Fq = False

From 40769a5b2de74053d110c1571b9de735827262c4 Mon Sep 17 00:00:00 2001
From: IndigoCarmine <tamanisikaminai@gmail.com>
Date: Tue, 17 Mar 2026 05:37:05 +0900
Subject: [PATCH 2/9] Fix incorrect references to cited papers.

---
 sasmodels/models/tetrapod.py | 7 +++----
 1 file changed, 3 insertions(+), 4 deletions(-)

diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py
index d0c51b94..f2658268 100644
--- a/sasmodels/models/tetrapod.py
+++ b/sasmodels/models/tetrapod.py
@@ -41,13 +41,13 @@
 References
 ----------
 
-#. J. S. Pedersen, *J. Appl. Cryst.*, 33 (2000) 488-494
+#. Seoki Kyoo Seo *Korean J. Chem. Eng.* 34(2017) 1192-1198
 
 Authorship and Verification
 ----------------------------
 
-* **Author:** NIST IGOR/DANSE **Date:** pre 2010
-* **Last Modified by:** Paul Butler **Date:** March 20, 2016
+* **Author:** Yuhei Yamada
+* **Last Modified by:**
 """
 
 import numpy as np
@@ -70,5 +70,4 @@
 ]
 
 source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "tetrapod.c"]
-valid = "radius >= 0.0 && length >= 0.0"
 have_Fq = False

From c995274dbea32d1c2f195b7eb7e5ce886190b73f Mon Sep 17 00:00:00 2001
From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com>
Date: Mon, 16 Mar 2026 20:43:08 +0000
Subject: [PATCH 3/9] Initial plan


From 4399102559eccf0b2e925dd66ca8423682106ece Mon Sep 17 00:00:00 2001
From: "copilot-swe-agent[bot]" <198982749+Copilot@users.noreply.github.com>
Date: Mon, 16 Mar 2026 21:03:40 +0000
Subject: [PATCH 4/9] Fix u_n formula in tetrapod.c and update docstring

Co-authored-by: IndigoCarmine <79552713+IndigoCarmine@users.noreply.github.com>
---
 sasmodels/models/tetrapod.c  |  7 +++--
 sasmodels/models/tetrapod.py | 54 +++++++++++++++++++++++++-----------
 2 files changed, 42 insertions(+), 19 deletions(-)

diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c
index 622e4b83..f2bf7c22 100644
--- a/sasmodels/models/tetrapod.c
+++ b/sasmodels/models/tetrapod.c
@@ -3,8 +3,9 @@ const double A =
 
 static double u_n(int n, double theta, double alpha) {
   const double phi[4] = {0.0, M_PI_2, M_PI, 3.0 * M_PI_2};
-  return cos(A) * cos(theta + phi[n] * 2) +
-         sin(A) * sin(theta + phi[n] * 2) * cos(alpha - phi[n]);
+  const double sign[4] = {1.0, -1.0, 1.0, -1.0};
+  return sign[n] * cos(A) * cos(theta) +
+         sin(A) * sin(theta) * cos(alpha - phi[n]);
 }
 
 static double A_n(double q, double u, double L, double R) {
@@ -39,7 +40,7 @@ static double Iq(double q, double L, double R, double sld_particle,
     double integral_alpha = 0.0;
     for (int dalpha = 0; dalpha < GAUSS_N; dalpha++) {
       double alpha =
-          M_PI * (GAUSS_Z[dalpha] + 1.0);  // map from [-1, 1] to [0, pi]
+          M_PI * (GAUSS_Z[dalpha] + 1.0);  // map from [-1, 1] to [0, 2*pi]
       double w_alpha =
           GAUSS_W[dalpha] * M_PI;  // adjust weight for the new range
 
diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py
index f2658268..54d0cc10 100644
--- a/sasmodels/models/tetrapod.py
+++ b/sasmodels/models/tetrapod.py
@@ -6,35 +6,52 @@
 of four cylindrical arms radiating from a central point, oriented along the
 (1,1,1), (-1,-1,1), (-1,1,-1), and (1,-1,-1) directions.
 
-The scattering intensity is calculated as:
+The scattering intensity is calculated as a powder average over all
+orientations:
 
 .. math::
 
-    I(q) = (\Delta \rho)^2 \left[ 4 I_{\text{arm}}(q) + 6 I_{\text{corr}}(q) \right]
+    I(q) = \frac{(\Delta \rho)^2}{4\pi}
+        \int_0^{\pi} \int_0^{2\pi}
+        \left|\sum_{n=1}^{4} F_n(q, \theta, \varphi)\right|^2
+        \sin\theta \, d\theta \, d\varphi
 
-where $\Delta\rho$ is the scattering length density contrast between the
-tetrapod and the solvent.
-
-The arm scattering $I_{\text{arm}}(q)$ is obtained by integrating over all
-orientations of each cylindrical arm:
+where $\Delta\rho$ is the SLD contrast and $F_n$ is the form factor amplitude
+of the $n$-th cylindrical arm:
 
 .. math::
 
-    I_{\text{arm}}(q) = \int_0^{\pi/4} F_{\text{cylinder}}(q, \theta, L, R)^2 \sin(\theta) d\theta
+    F_n(q, \theta, \varphi) = \text{sinc}\!\left(\frac{q u_n L}{2}\right)
+        \cdot \frac{2 J_1(q \mu_n R)}{q \mu_n R}
+
+with $u_n = \hat{q} \cdot \hat{a}_n$ the projection of $\hat{q}$ onto the
+arm axis, $\mu_n = \sqrt{1 - u_n^2}$, $L$ the arm length, and $R$ the arm
+radius.  Expanding the squared modulus into a double sum and exploiting the
+reality of $F_n$ gives:
 
-where $F_{\text{cylinder}}$ is the form factor of a cylinder.
+.. math::
 
-The correlation term $I_{\text{corr}}(q)$ accounts for coherent scattering
-between different arms.
+    I(q) = \frac{(\Delta \rho)^2}{4\pi}
+        \int \sum_{n=1}^{4}\sum_{m=1}^{4}
+        F_n F_m \cos\!\left(\frac{q(u_n-u_m)L}{2}\right)
+        \sin\theta \, d\theta \, d\varphi
+
+The cosine factor is the interference term between the centres of arms $n$
+and $m$, which are displaced by $\tfrac{L}{2}\hat{a}_n$ from the junction.
 
 Geometry
 --------
 
-The four arms are positioned along the following directions:
-- Arm 1: (1, 1, 1)
-- Arm 2: (-1, -1, 1)
-- Arm 3: (-1, 1, -1)
-- Arm 4: (1, -1, -1)
+The four arms are oriented along tetrahedral directions.  With
+$A = 109.5°/2$ (the half-angle between arms), the arm unit vectors and the
+corresponding projections $u_n$ are
+
+.. math::
+
+    u_n = s_n \cos A \cos\theta + \sin A \sin\theta \cos(\varphi - \varphi_n)
+
+where $(s_n, \varphi_n) = (+1,\ 0),\ (-1,\ \pi/2),\ (+1,\ \pi),\ (-1,\ 3\pi/2)$
+for $n = 1, 2, 3, 4$ respectively.
 
 Each arm has length $L$ and radius $R$.
 
@@ -71,3 +88,8 @@
 
 source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "tetrapod.c"]
 have_Fq = False
+
+tests = [
+    # Default parameters: length=400, radius=30, sld=4, sld_solvent=1
+    [{}, 0.1, 1.0026479478514564e-3],
+]

From 03c63fbbc53501098c52d5d2353f6c3c754bffe2 Mon Sep 17 00:00:00 2001
From: IndigoCarmine <tamanisikaminai@gmail.com>
Date: Wed, 18 Mar 2026 12:24:34 +0900
Subject: [PATCH 5/9] Refactor Fq_n function in tetrapod.c and remove text
 generated through hallucination

---
 sasmodels/models/tetrapod.c  | 12 ++++--------
 sasmodels/models/tetrapod.py |  7 +------
 2 files changed, 5 insertions(+), 14 deletions(-)

diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c
index f2bf7c22..4ea80584 100644
--- a/sasmodels/models/tetrapod.c
+++ b/sasmodels/models/tetrapod.c
@@ -8,17 +8,13 @@ static double u_n(int n, double theta, double alpha) {
          sin(A) * sin(theta) * cos(alpha - phi[n]);
 }
 
-static double A_n(double q, double u, double L, double R) {
+static double Fq_n(double q, double u, double L, double R) {
   double quL2 = q * u * L * 0.5;
 
-  double term1 = (fabs(quL2) > 1e-12) ? sin(quL2) / quL2 : 1.0;
-
   double mu = sqrt(fmax(0.0, 1.0 - u * u));
   double qmuR = q * mu * R;
 
-  double term2 = (fabs(qmuR) > 1e-12) ? 2.0 * sas_J1(qmuR) / qmuR : 1.0;
-
-  return term1 * term2;
+  return sas_sinx_x(quL2) * sas_2J1x_x(qmuR);
 }
 
 static double form_volume(double L, double R) {
@@ -49,7 +45,7 @@ static double Iq(double q, double L, double R, double sld_particle,
         for (int m = 0; m < 4; m++) {
           double u = u_n(n, theta, alpha);
           double v = u_n(m, theta, alpha);
-          sum_arms += A_n(q, u, L, R) * A_n(q, v, L, R) *
+          sum_arms += Fq_n(q, u, L, R) * Fq_n(q, v, L, R) *
                       cos(q * (u - v) * L / 2.0) * sin(theta);
         }
       }
@@ -57,5 +53,5 @@ static double Iq(double q, double L, double R, double sld_particle,
     }
     total += integral_alpha * w_theta;
   }
-  return contrast * contrast * total * 2 * M_PI;
+  return contrast * contrast * total * 2 / M_PI;
 }
\ No newline at end of file
diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py
index 54d0cc10..0dff0b0b 100644
--- a/sasmodels/models/tetrapod.py
+++ b/sasmodels/models/tetrapod.py
@@ -6,7 +6,7 @@
 of four cylindrical arms radiating from a central point, oriented along the
 (1,1,1), (-1,-1,1), (-1,1,-1), and (1,-1,-1) directions.
 
-The scattering intensity is calculated as a powder average over all
+The scattering intensity is calculated as an average over all
 orientations:
 
 .. math::
@@ -88,8 +88,3 @@
 
 source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "tetrapod.c"]
 have_Fq = False
-
-tests = [
-    # Default parameters: length=400, radius=30, sld=4, sld_solvent=1
-    [{}, 0.1, 1.0026479478514564e-3],
-]

From 7037bea3789b90ca531953fa71cb9cac18acbfd2 Mon Sep 17 00:00:00 2001
From: IndigoCarmine <tamanisikaminai@gmail.com>
Date: Wed, 18 Mar 2026 12:26:16 +0900
Subject: [PATCH 6/9] Remove unused numpy imports in tetrapod.py

---
 sasmodels/models/tetrapod.py | 3 +--
 1 file changed, 1 insertion(+), 2 deletions(-)

diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py
index 0dff0b0b..fa9e3a80 100644
--- a/sasmodels/models/tetrapod.py
+++ b/sasmodels/models/tetrapod.py
@@ -67,8 +67,7 @@
 * **Last Modified by:**
 """
 
-import numpy as np
-from numpy import cos, inf, pi, sin
+from numpy import inf
 
 name = "tetrapod"
 title = "Tetrapod with four cylindrical arms"

From 6a2524746f6d4d6ce82bbd9eacc8f2d38d4fd91e Mon Sep 17 00:00:00 2001
From: IndigoCarmine <tamanisikaminai@gmail.com>
Date: Sat, 28 Mar 2026 02:16:06 +0900
Subject: [PATCH 7/9] Update Iq function to include form_volume in return value
 and enhance author attribution in tetrapod.py

---
 sasmodels/models/tetrapod.c  | 2 +-
 sasmodels/models/tetrapod.py | 3 ++-
 2 files changed, 3 insertions(+), 2 deletions(-)

diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c
index 4ea80584..4bcc14fb 100644
--- a/sasmodels/models/tetrapod.c
+++ b/sasmodels/models/tetrapod.c
@@ -53,5 +53,5 @@ static double Iq(double q, double L, double R, double sld_particle,
     }
     total += integral_alpha * w_theta;
   }
-  return contrast * contrast * total * 2 / M_PI;
+  return contrast * contrast * total * 2 / M_PI * form_volume(L, R);
 }
\ No newline at end of file
diff --git a/sasmodels/models/tetrapod.py b/sasmodels/models/tetrapod.py
index fa9e3a80..eae3521f 100644
--- a/sasmodels/models/tetrapod.py
+++ b/sasmodels/models/tetrapod.py
@@ -63,7 +63,7 @@
 Authorship and Verification
 ----------------------------
 
-* **Author:** Yuhei Yamada
+* **Author:** Yuhei Yamada (Github user name: Indigo Carmine, https://orcid.org/0009-0003-9780-4135)
 * **Last Modified by:**
 """
 
@@ -87,3 +87,4 @@
 
 source = ["lib/polevl.c", "lib/sas_J1.c", "lib/gauss76.c", "tetrapod.c"]
 have_Fq = False
+opencl = True

From f7de9f6e56d82ef82227a87dc522fe66ef64132e Mon Sep 17 00:00:00 2001
From: IndigoCarmine <tamanisikaminai@gmail.com>
Date: Sat, 28 Mar 2026 02:35:24 +0900
Subject: [PATCH 8/9] Fix Iq function to remove form_volume calculation and
 adjust return value

---
 sasmodels/models/tetrapod.c | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c
index 4bcc14fb..aeeb8c07 100644
--- a/sasmodels/models/tetrapod.c
+++ b/sasmodels/models/tetrapod.c
@@ -53,5 +53,5 @@ static double Iq(double q, double L, double R, double sld_particle,
     }
     total += integral_alpha * w_theta;
   }
-  return contrast * contrast * total * 2 / M_PI * form_volume(L, R);
+  return 1e-4 * contrast * contrast * total * 2 / M_PI;
 }
\ No newline at end of file

From 26b9812b552c099aaaab25091d7ac79d3c171d2c Mon Sep 17 00:00:00 2001
From: IndigoCarmine <tamanisikaminai@gmail.com>
Date: Sat, 28 Mar 2026 02:48:58 +0900
Subject: [PATCH 9/9] Refactor angle calculation for clarity in tetrapod.c

---
 sasmodels/models/tetrapod.c | 4 ++--
 1 file changed, 2 insertions(+), 2 deletions(-)

diff --git a/sasmodels/models/tetrapod.c b/sasmodels/models/tetrapod.c
index aeeb8c07..b0189fd9 100644
--- a/sasmodels/models/tetrapod.c
+++ b/sasmodels/models/tetrapod.c
@@ -1,5 +1,5 @@
-const double A =
-    109.5 / 2.0 * M_PI / 180.0;  // half of the angle between arms in radians
+// half of the angle between arms in radians
+const double A = acos(-1 / 3) / 2.0 * M_PI / 180.0;
 
 static double u_n(int n, double theta, double alpha) {
   const double phi[4] = {0.0, M_PI_2, M_PI, 3.0 * M_PI_2};