In my experience there is another very simple offsetting method that is better for "simple" curves. I think it is (or was?) used in paper.js.

You take the intersection point of the curves's handles and offset curve's handles. Then you calculate the distances from these points to the curve's/offset curve's end points and use these distances to scale the handles of the offset curve.

Here is the code, reformatted so you can basically copy it into your existing code. It should require fewer subdivisions than TH, but more than shape control and your impressive solution.

    // This function needs to go into the Point class
    scale(s) {
        return new Point(this.x * s, this.y * s);
    }

    handle_scaling() {
        // p: points on original curve
        // q: points on offset curve
        const deriv = this.c.deriv();
        const p0 = this.c.p0();
        const p1 = this.c.p1();
        const p2 = this.c.p2();
        const p3 = this.c.p3();
        const q0 = p0.plus(this.eval_offset(0));
        const q3 = p3.plus(this.eval_offset(1));
        const tan0 = deriv.eval(0);
        const tan3 = deriv.eval(1);
        // s: intersection of handle vectors 
        const sp = ray_intersect(p0, tan0, p3, tan3);
        const sq = ray_intersect(q0, tan0, q3, tan3);
        // r: ratios of distances from handle intersection to end points
        const r0 = sq.dist(q0) / sp.dist(p0);
        const r3 = sq.dist(q3) / sp.dist(p3);
        // calculate control points of offset curve
        let q1;
        let q2;
        if (r0 > 0 && r3 > 0) {
            q1 = q0.plus(p1.minus(p0).scale(sq.dist(q0) / sp.dist(p0)));
            q2 = q3.plus(p2.minus(p3).scale(sq.dist(q3) / sp.dist(p3)));
        } else {
            q1 = p1.minus(p0).plus(q0);
            q2 = p2.minus(p3).plus(q3);
        }
        return CubicBez.from_pts(q0, q1, q2, q3);
    }