Как вы используете Apache Commons Math 3.0+ для решения ОДУ с помощью якобиана?

Что касается конкретной библиотеки Java для решения ОДУ с приближением Якоби, предыдущая библиотека Apache Commons Math 2.2 интуитивно понятна, и также есть несколько наглядных примеров:

Интеграция обыкновенных дифференциальных уравнений

Но 2.2 сейчас устарела, а недавняя библиотека 3.0+

Apache Commons Math 3.0

особенно пакеты:

org.apache.commons.math3.ode

org.apache.commons.math3.ode.events

org.apache.commons.math3.ode.nonstiff

являются заменами, но кажутся гораздо более сложными. Я думаю, что новые версии ODE были мотивированы обобщением структур данных и методов, но трудно понять, как некоторые из классов предназначены для совместного использования (например, MainStateJacobianProvider, ExpandableStatefulODE, ParameterizedODE), и нигде нет примеров этого, пока Apache не обновил свои информационные страницы. Документация по API не очень информативна и тоже кажется неоднозначной.

Поэтому вместо того, чтобы звучать так, будто я жалуюсь на все подряд, я задам конкретный вопрос: Как бы вы перекодировали ту же проблему из одного из примеров, продемонстрированных выше, с использованием 2.2, используя библиотеку 3.0?

Вот код:

public class BasicCircleODE implements ParameterizedODE {

    private double[] c;
    private double omega;

    public BasicCircleODE(double[] c, double omega) {
        this.c = c;
        this.omega = omega;
    }

    public int getDimension() {
        return 2;
    }

    public void computeDerivatives(double t, double[] y, double[] yDot) {
        yDot[0] = omega * (c[1] - y[1]);
        yDot[1] = omega * (y[0] - c[0]);
    }   

    public int getParametersDimension() {
        // we are only interested in the omega parameter
        return 1;
    }

    public void setParameter(int i, double value) {
        omega = value;
    }

}

double[] hY = new double[] { 0.001, 0.001 };
double[] hP = new double[] { 1.0e-6 };
FirstOrderIntegratorWithJacobians integrator = new FirstOrderIntegratorWithJacobians(dp853, ode, hY, hP);
integrator.integrate(t0, y0, dy0dp, t, y, dydy0, dydp);