Function cart_to_oblate_velocity

examples/pde/bhrt.cpp:366–397

* @brief Transform the velocity vectors from cartesian to oblate coordinates * * @param vel * @param pos * @return af::array */

 * @return af::array
 */
af::array cart_to_oblate_velocity(const af::array& vel, const af::array& pos) {
    if (vel.dims() != pos.dims())
        throw make_error("Arrays must have the same dimensions");
    else if (pos.dims()[0] != 3)
        throw make_error("Arrays must have 3 principal coordintes");

    af::array x = pos(0, af::span);
    af::array y = pos(1, af::span);
    af::array z = pos(2, af::span);

    auto a    = J / M;
    auto diff = x * x + y * y + z * z - a * a;

    af::array r =
        af::sqrt((diff + af::sqrt(diff * diff + z * z * a * a * 4.0)) / 2.0);
    af::array o = af::acos(z / r);
    af::array p = af::atan2(y, x);

    af::array ux = vel(0, af::span);
    af::array uy = vel(1, af::span);
    af::array uz = vel(2, af::span);

    af::array ra = r * r + a * a;
    af::array ur = ((ux * x + uy * y) * r + uz * ra * z / r) /
                   (r * r + af::pow(a * af::cos(o), 2.0));
    af::array up = (uy * x - ux * y) / (x * x + y * y);
    af::array uo = ((ux * x + uy * y) / af::tan(o) - uz * z * af::tan(o)) /
                   (r * r + af::pow(a * af::cos(o), 2.0));
    af::array transformed_vel = af::join(0, ur, uo, up);

    return transformed_vel;
}

/**
 * @brief Transform the position vectors from spherical to cartesian coordinates