geomstats.numerics package#

Submodules#

geomstats.numerics.bvp module#

Boundary value problem solvers implementation.

class geomstats.numerics.bvp.ScipySolveBVP(tol=0.001, max_nodes=1000, bc_tol=None, save_result=False)[source]#

Bases: object

Wrapper for scipy.integrate.solve_bvp.

integrate(fun, bc, x, y, fun_jac=None, bc_jac=None)[source]#

Solve a boundary value problem for a system of ODEs.

geomstats.numerics.geodesic module#

Geodesic solvers implementation.

class geomstats.numerics.geodesic.ExpODESolver(integrator=None)[source]#

Bases: ExpSolver

Geodesic initial value problem solver.

Parameters:

integrator (ODEIVPSolver) – Instance of ODEIVP integrator.

exp(space, tangent_vec, base_point)[source]#

Exponential map.

Parameters:
  • space (Manifold) – Equipped manifold.

  • tangent_vec (array-like, shape=[…, *space.shape]) – Tangent vector at the base point.

  • base_point (array-like, shape=[…, *space.shape]) – Point on the manifold.

Returns:

end_point (array-like, shape=[…, *space.shape]) – Point on the manifold.

geodesic_ivp(space, tangent_vec, base_point)[source]#

Geodesic curve for initial value problem.

Parameters:
  • space (Manifold) – Equipped manifold.

  • tangent_vec (array-like, shape=[…, *space.shape]) – Tangent vector at the base point.

  • base_point (array-like, shape=[…, *space.shape]) – Point on the manifold.

Returns:

path (callable) – Time parametrized geodesic curve. f(t).

class geomstats.numerics.geodesic.ExpSolver[source]#

Bases: ABC

Abstract class for geodesic initial value problem solvers.

abstract exp(space, tangent_vec, base_point)[source]#

Exponential map.

Parameters:
  • space (Manifold) – Equipped manifold.

  • tangent_vec (array-like, shape=[…, dim]) – Tangent vector at the base point.

  • base_point (array-like, shape=[…, dim]) – Point on the manifold.

Returns:

end_point (array-like, shape=[…, dim]) – Point on the manifold.

abstract geodesic_ivp(space, tangent_vec, base_point)[source]#

Geodesic curve for initial value problem.

Parameters:
  • space (Manifold) – Equipped manifold.

  • tangent_vec (array-like, shape=[…, dim]) – Tangent vector at the base point.

  • base_point (array-like, shape=[…, dim]) – Point on the manifold.

Returns:

path (callable) – Time parametrized geodesic curve. f(t).

class geomstats.numerics.geodesic.LogODESolver(n_nodes=10, integrator=None, initialization=None, use_jac=True)[source]#

Bases: _LogBatchMixins, LogSolver

Geodesic boundary value problem using an ODE solver.

Parameters:
  • n_nodes (Number of mesh nodes.)

  • integrator (ScipySolveBVP) – Instance of ScipySolveBVP.

  • initialization (callable) – Function to provide initial solution. f(space, point, base_point). Defaults to linear initialization.

geodesic_bvp(space, point, base_point)[source]#

Geodesic curve for boundary value problem.

Parameters:
  • space (Manifold) – Equipped manifold.

  • end_point (array-like, shape=[…, dim]) – Point on the manifold.

  • base_point (array-like, shape=[…, dim]) – Point on the manifold.

Returns:

path (callable) – Time parametrized geodesic curve. f(t). 0 <= t <= 1.

class geomstats.numerics.geodesic.LogShootingSolver(optimizer=None, initialization=None, flatten=True)[source]#

Bases: object

Geodesic boundary value problem solver using shooting.

Parameters:
  • optimizer (ScipyMinimize) – Instance of ScipyMinimize.

  • initialization (callable) – Function to provide initial solution. f(space, point, base_point). Defaults to linear initialization.

  • flatten (bool) – If True, the optimization problem is solved together for all the points.

class geomstats.numerics.geodesic.LogSolver[source]#

Bases: ABC

Abstract class for geodesic boundary value problem solvers.

abstract geodesic_bvp(space, point, base_point)[source]#

Geodesic curve for boundary value problem.

Parameters:
  • space (Manifold) – Equipped manifold.

  • end_point (array-like, shape=[…, dim]) – Point on the manifold.

  • base_point (array-like, shape=[…, dim]) – Point on the manifold.

Returns:

path (callable) – Time parametrized geodesic curve. f(t).

abstract log(space, point, base_point)[source]#

Logarithm map.

Parameters:
  • space (Manifold) – Equipped manifold.

  • end_point (array-like, shape=[…, dim]) – Point on the manifold.

  • base_point (array-like, shape=[…, dim]) – Point on the manifold.

Returns:

tangent_vec (array-like, shape=[…, dim]) – Tangent vector at the base point.

geomstats.numerics.ivp module#

Initial value problem solvers implementation.

class geomstats.numerics.ivp.GSIVPIntegrator(n_steps=10, step_type='euler', save_result=False)[source]#

Bases: ODEIVPSolver

In-house ODE integrator.

Parameters:
  • n_steps (int) – Number of steps to perform.

  • step_type (str) – Type of integration step. Possible values are euler, rk2, rk4.

  • save_result (bool) – If True, result is stored after calling integrate or integrate_t.

integrate(force, initial_state, end_time=1.0)[source]#

Integrate force.

Parameters:
  • force (callable) – Function to integrate.

  • initial_state (array-like) – Initial state.

  • end_time (float or None) – Integration end time.

Returns:

result (OdeResult)

integrate_t(force, initial_state, t_eval)[source]#

Integrate force at t_eval points.

Parameters:
  • force (callable) – Function to integrate.

  • initial_state (array-like) – Initial state.

  • t_eval (array-like) – Times at which to store the computed solution.

Returns:

result (OdeResult)

property step_type#

Integrator step type.

class geomstats.numerics.ivp.ODEIVPSolver(save_result=False, state_is_raveled=False, tfirst=False)[source]#

Bases: ABC

Abstract class for ode ivp solvers.

Parameters:
  • save_result (bool) – If True, result is stored after calling integrate or integrate_t.

  • state_is_raveled (bool) – If True, state is represented by a 1d array. Else, it is represented by (n_vars, dim).

  • tfirst (bool) – Declares function signature. If True f(t, y), else f(y, t).

abstract integrate(force, initial_state, end_time)[source]#

Integrate force.

Parameters:
  • force (callable) – Function to integrate.

  • initial_state (array-like) – Initial state.

  • end_time (float or None) – Integration end time.

Returns:

result (OdeResult)

abstract integrate_t(force, initial_state, t_eval)[source]#

Integrate force while choosing evaluating points.

Parameters:
  • force (callable) – Function to integrate.

  • initial_state (array-like) – Initial state.

  • t_eval (array-like) – Times at which to store the computed solution.

Returns:

result (OdeResult)

class geomstats.numerics.ivp.OdeResult[source]#

Bases: OptimizeResult

Bunch object (follows scipy).

Its purposes is to homogenize output of different integrators.

get_last_y()[source]#

Get value for last y.

Allows to have y represented as an gs.array or list[gs.array] (latter for cases where they have different shapes).

Assumes last t is the same.

class geomstats.numerics.ivp.ScipySolveIVP(method='RK45', save_result=False, **options)[source]#

Bases: ODEIVPSolver

Wrapper for scipy.integrate.solve_ivp.

Parameters:
  • method (str) – Integration method.

  • save_result (bool) – If True, result is stored after calling integrate or integrate_t.

integrate(force, initial_state, end_time=1.0)[source]#

Integrate force.

Parameters:
  • force (callable) – Function to integrate.

  • initial_state (array-like) – Initial state.

  • end_time (float or None) – Integration end time.

Returns:

result (OdeResult)

integrate_t(force, initial_state, t_eval)[source]#

Integrate force at t_eval points.

Parameters:
  • force (callable) – Function to integrate.

  • initial_state (array-like) – Initial state.

  • t_eval (array-like) – Times at which to store the computed solution.

Returns:

result (OdeResult)

geomstats.numerics.optimizers module#

Optimizers implementations.

class geomstats.numerics.optimizers.ScipyMinimize(method='L-BFGS-B', jac=None, hess=None, bounds=None, constraints=(), tol=None, callback=None, options=None, save_result=False)[source]#

Bases: object

Wrapper for scipy.optimize.minimize.

Only jac differs from scipy: if autodiff, then gs.autodiff.value_and_grad is used to compute the jacobian.

minimize(fun, x0, fun_jac=None, fun_hess=None, hessp=None)[source]#

Minimize objective function.

Parameters:
  • fun (callable) – The objective function to be minimized.

  • x0 (array-like) – Initial guess.

  • fun_jac (callable) – If not None, jac is ignored.

  • fun_hess (callable) – If not None, hess is ignored.

  • hessp (callable)

Module contents#