diff --git a/extras/cc_nnlo_ic/.gitignore b/extras/cc_nnlo_ic/.gitignore
new file mode 100644
index 00000000..1364d379
--- /dev/null
+++ b/extras/cc_nnlo_ic/.gitignore
@@ -0,0 +1,6 @@
+# Kirill
+2601.02916*.pdf
+CF_HQI_CC_DIS
+# GPL
+0410259*.pdf
+Multiple-Polylogarithm
diff --git a/extras/cc_nnlo_ic/README.md b/extras/cc_nnlo_ic/README.md
new file mode 100644
index 00000000..c445945f
--- /dev/null
+++ b/extras/cc_nnlo_ic/README.md
@@ -0,0 +1,4 @@
+# CC NNLO intrinsic
+Implementation of [2601.02916](https://arxiv.org/abs/2601.02916).
+
+Clone the results repo here: `git clone git@github.com:Qdashkin/CF_HQI_CC_DIS.git`.
diff --git a/extras/cc_nnlo_ic/mma.py b/extras/cc_nnlo_ic/mma.py
new file mode 100644
index 00000000..2e9a310b
--- /dev/null
+++ b/extras/cc_nnlo_ic/mma.py
@@ -0,0 +1,71 @@
+"""Python interface to Mathematica via the CLI."""
+
+import io
+import shutil
+import subprocess
+import threading
+from typing import Self
+
+
+def _thread_writer(p: subprocess.Popen, res: io.StringIO) -> None:
+    """
+    Thread worker.
+
+    Reads from the stdout of the subprocess until it is closed or an `@` is
+    encountered.
+
+    Parameters
+    ----------
+        p : subprocess.Popen
+            subprocess
+        res : stream
+            stream to which the output is copied
+    """
+    while True:
+        # print("read data: ")
+        data = p.stdout.read(1).decode("utf-8")
+        if not data or data == "@":
+            break
+        res.write(data)
+        res.flush()
+
+
+class MmaRunner:
+    """Call Mathematica interactively from Python via the CLI."""
+
+    p: subprocess.Popen
+    """CLI thread"""
+
+    def __enter__(self: Self) -> Self:
+        path = shutil.which("math")
+        assert path is not None
+        self.p = subprocess.Popen(
+            [path, "-noprompt"], stdin=subprocess.PIPE, stdout=subprocess.PIPE
+        )
+        return self
+
+    def send(self: Self, code: str) -> str:
+        """Send Mathematica code to the running CLI.
+
+        An explicit `@` is added as the EOF signal to the thread worker.
+        """
+        # reading needs to be on a seperate stream - idea from here:
+        # https://stackoverflow.com/questions/19880190/interactive-input-output-using-python/53312631#53312631
+        stream = io.StringIO()
+        writer = threading.Thread(target=_thread_writer, args=(self.p, stream))
+        writer.start()
+
+        self.p.stdin.write((code + 'Print["@"];\n').encode())
+        self.p.stdin.flush()
+        writer.join()
+        stream.seek(0)
+        s = stream.read()
+        return s[1:-2].strip()
+
+    def __exit__(self: Self, exc_type: type, _exc_value, _traceback):
+        """Close the interactive CLI, by sending the `Exit[]` command."""
+        if exc_type is not None:
+            return
+        if self.p is not None:
+            self.send("Exit[];")
+            self.p.terminate()
diff --git a/extras/cc_nnlo_ic/run.py b/extras/cc_nnlo_ic/run.py
new file mode 100644
index 00000000..ca8e2098
--- /dev/null
+++ b/extras/cc_nnlo_ic/run.py
@@ -0,0 +1,64 @@
+"""Check and import Kirill."""
+
+import numpy as np
+from mma import MmaRunner
+
+from yadism.coefficient_functions.intrinsic import f2_cc
+
+PATH = "CF_HQI_CC_DIS/strFun.mx"
+
+INIT = f'Get["{PATH}"];' + r"""
+colorRule = {cF -> 4/3, cA -> 3}
+ruleU = {u -> Sqrt[y1/(y1 - 4 y2)]};
+loResult = {Splus/2, (Splus xBj)/(1 - \[Beta]), 2 Rplus};
+simpleGRule = {G[0,x_] :> Log[x]};
+y2SingularRule = {y2 -> m2/(Q2+m2)};
+"""
+
+
+class MockESF:
+    Q2: float
+
+
+def check_lo(r: MmaRunner) -> None:
+    """Check LO is a delta function."""
+    for sf in [1, 2, 3]:
+        lo = r.send(rf"Print[L[{sf}, 0, 1] /. LCoefRules];")
+        print(f"F_{sf}|LO = {lo} which is {lo == 'PD[0]'}")
+
+
+def check_nlo(r: MmaRunner) -> None:
+    """Check NLO."""
+    Q2 = 1.0
+    m2 = 1.0
+    nf = 3
+    esf = MockESF()
+    esf.Q2 = Q2
+    yad = f2_cc.Splus(esf, nf, m1sq=m2)
+    # x = 0.1
+    # yad_lo = yad.LO()
+    # print(yad_lo.loc(x,yad_lo.args["loc"]))
+    yad_nlo = yad.NLO()
+    z = 0.1
+    # TODO: there is a magic 2 still
+    yad_nlo_sing = yad_nlo.sing(z, yad_nlo.args["sing"]) / yad.lo / 2.0
+    sf = 2
+    kk_nlo_sing = r.send(rf"""
+    Block[{{sf={sf},cf,sing}},
+        cf = cF L[sf, 1, 1] /. LCoefRules;
+        sing = Coefficient[cf,PD@2] * Log[1-z]/(1-z) + Coefficient[cf,PD@1] * 1/(1-z);
+        sing = sing /. colorRule /. simpleGRule /. y2SingularRule /. {{m2 -> {m2}, Q2 -> {Q2}, z -> {z}}};
+        Print@sing;
+    ]""")
+    kk_nlo_sing = float(kk_nlo_sing)
+    print(
+        f"F_{sf}|NLO|sing = {kk_nlo_sing} which is {np.isclose(kk_nlo_sing, yad_nlo_sing)}"
+    )
+
+
+if __name__ == "__main__":
+    # initialize
+    with MmaRunner() as runner:
+        runner.send(INIT)
+        # check_lo(runner)
+        check_nlo(runner)