For future reference, here’s a tridiagonal solver in vex (very much not optimised, but it’s based on some Mathematica code) @jake rice 😮
“` #include
function complex cdiv(complex a, b) { complex conjB = conj(b); complex num = cmult(a, conjB); complex den = cmult(b, conjB); complex ans = cmult(num, CMPLX(1.0/den.real, 0)); return ans; }
function complex[] TDMAsolver(complex a[], b[], c[], d[]) { int nf = len(d); // number of equations complex p[]; complex q[]; for (int i=0;i