Bubbz's profile - activity

For example in the code that is in my post, whatever the impact parameter al maybe, however large, if b=0 meaning that w

Thanks for your response. My issue is not with the horizon but with the axis of the object. For example in the edge-on c

Switching coordinate charts mid-integration with solve_across_charts I am trying to utilize solve_across_charts to avoid

As it turns out, even though the error is raised, the restriction does pass on the coordinate_chart which is even weirde

Coordinate chart restriction raising a floating-point-overflow error I am using the 9.3 version and what i am trying to

Coordinate chart restriction raising a floating-point-overflow error I am using the 9.3 version and what i am trying to

I am trying to figure out how to load a previously computed and saved integrated geodesic and find the tangent vector using geod.tangent_vector_eval_at at a later time. The problem is that i can not seem to find a way to refer a solution_key= to the solution after i've loaded it in order to interpolate it. Is there a way to do that at a time different than the one that the solution is being computed? Thanks in advance.

Evaluating the tangent vector of a previously saved and loaded integrated geodesic I am trying to figure out how to load

Workaround works just fine ! thank you very much!

Hi there. I am trying to evaluate the tangent vector at a specific value of the affine parameter. I have tried everything from https://doc.sagemath.org/html/en/refe... but it raises an error i can't find a solution to.

edit: Just to be precise and avoid taking up anyone's time. I am actually trying to create a poincare section of my orbit so i think i could use the numerical solution by calculating the velocity where the sol meets the criteria. Maybe there's an easier way (?). Couldnt think of anything else that can help me find the poincare section after and not during the integration.

M = Manifold(4, 'M', latex_name=r'\mathcal{M}', structure='Lorentzian')
X.<t,r,th,ph> = M.chart(r't r:(2.0001,+inf) th:(0,pi):\theta ph:(0,2*pi):periodic\phi') 
X.coord_range()

var('m,b,a', domain='real')
g = M.metric()
m=1

g[0,0] = -(1-2*m*r/(r^2 + (a*cos(th))^2))
g[0,3] = -2*a*m*r*sin(th)^2/(r^2 + (a*cos(th))^2)
g[1,1] = (r^2 + (a*cos(th))^2)/(r^2 -2*m*r + a^2)
g[2,2] = (r^2 + (a*cos(th))^2)
g[3,3] = (r^2+a^2+2*m*r*(a*sin(th))^2/(r^2 + (a*cos(th))^2))*sin(th)^2
g.display()

def g00(m,r,th,a):return g[0,0](m,r,th,a)
def g03(m,r,th,a):return g[0,3](m,r,th,a)
def g11(m,r,th,a):return g[1,1](m,r,th,a)
def g22(m,r,th,a):return g[2,2](m,r,th,a)
def g33(m,r,th,a):return g[3,3](m,r,th,a)
def D(m,r,th,a): return (g03(m,r,th,a))^2-g00(m,r,th,a)*g33(m,r,th,a)

E.<x,y,z> = EuclideanSpace()
phi = M.diff_map(E, [r*sin(th)*cos(ph), r*sin(th)*sin(ph), r*cos(th)])
phi.display()

def initial_vector(r0, b,al, ph0=0, E=1, inward=False):
    t0,th0=0,pi/2
    L = -b*E
    vt0=(E*g33(m,r0,th0,a)+L*g03(m,r0,th0,a))/(D(m,r0,th0,a))
    vr0=sqrt(((1/D(m,r0,th0,a))*((L^2)*g00(m,r0,th0,a)+2*E*L*g03(m,r0,th0,a)+(E^2)*g33(m,r0,th0,a)))/g11(m,r0,th0,a))
    if inward:
        vr0 = - vr0
    vth0 = al/(100*r0)
    vph0 = -(1/D(m,r0,th0,a))*(E*g03(m,r0,th0,a)+L*(g00(m,r0,th0,a)))
    p0 = M((t0, r0, th0, ph0), name='p_0')
    return M.tangent_space(p0)((vt0, vr0, vth0, vph0), name='v_0')

r0 = 100 
th0=pi/2
s = var('s')

v0 = initial_vector(r0, al=1,b=1, ph0=0, inward=True)
geod = M.integrated_geodesic(g, (s, 0, 500), v0, across_charts=True)
sol = geod.solve_across_charts(step=0.01,method='odeint',parameters_values={a:1},
                               solution_key='sol',verbose=True) 

interp = geod.interpolate(solution_key='sol',
                          interpolation_key='interp 1', verbose=True)

where using:

v = geod.tangent_vector_eval_at(200, verbose=True)

raises :

Evaluating tangent vector components from the interpolation associated with the key 'interp 1' by default...

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-24-6fbce2379e9a> in <module>
----> 1 v = geod.tangent_vector_eval_at(Integer(200), verbose=True)

/opt/sagemath-9.2/local/lib/python3.7/site-packages/sage/manifolds/differentiable/integrated_curve.py in tangent_vector_eval_at(self, t, interpolation_key, verbose)
   2275             # partial test, in case future interpolation objects do not
   2276             # contain lists of instances of the Spline class
-> 2277             raise TypeError("unexpected type of interpolation object")
   2278 
   2279         interpolated_coordinates=[coordinate_curve_spline(t)

TypeError: unexpected type of interpolation object

'solve_across_charts()' works way better for me when dealing with geodesics that reach the event horizon. A photon that

'solve_across_charts()' works way better for me when dealing with geodesics that reach the event horizon. A photon that

i've included everything on the original post! thanks for your time! really appreciate it.

Geodesics/tangent vector evaluation interpolation error Hi there. I am trying to evaluate the tangent vector at a specif

Geodesics/tangent vector evaluation interpolation error Hi there. I am trying to evaluate the tangent vector at a specif

Geodesics/tangent vector evaluation interpolation error Hi there. I am trying to evaluate the tangent vector at a specif

I am actually trying to create a poincare section of my orbit so i think i could use the numerical solution by calculati

Of course. I've edited both into the original post. Haven't really used anything that different than the SageManifolds'

Geodesics/tangent vector evaluation interpolation error Hi there. I am trying to evaluate the tangent vector at a specif

Tangent vector evaluation/Interpolating Geodesic Error Hi there. I am trying to evaluate the tangent vector at a specifi

Interpolating Geodesic Error Hi there. I am trying to evaluate the tangent vector at a specific value of the affine para

Interpolating Geodesic Error Hi there. I am trying to evaluate the tangent vector at a specific value of the affine para

It could work if in X.add_restrictions(r>R(th)) , R(th) can somehow be a callable function(?). That way i should be a

Coordinate charts, functions of one coordinate as bounds of another(?) Hi there. I have a Manifold (4,M) and its coordin

Coordinate charts, functions of one coordinate as bounds of another(?) Hi there. I have a Manifold (4,M) and its coordin

Coordinate charts, functions of one coordinate as bounds of another(?) Hi there. I have a Manifold (4,M) and its coordin

I've added on the original post, one of the surfaces (r,th) that i need to pass somehow as a restriction for r. I've tri

Coordinate charts, functions of one coordinate as bounds of another(?) Hi there. I have a Manifold (4,M) and its coordin

Thanks for your answer! I'm sorry if I wasn't clear, r=th^2+0.52 was just an example and it does actually work with your

Coordinate charts, functions of one variable as bounds of another(?) Hi there. I have a Manifold (4,M) and its coordinat

Coordinate charts, functions of one variable as bounds of another(?) Hi there. I have a Manifold (4,M) and its coordinat

Is it possible to integrate a geodesic equation backwards in time in SageMath? I'm using the following code in order to ray trace null geodesics:

v0 = initial_vector(r0, b, al, ph0=0, inward=False)
geod = M.integrated_geodesic(g, (s, 50, 0), v0, across_charts=True)
sol = geod.solve_across_charts(step=0.1, parameters_values={a: 0, q: 0},
                               verbose=True)

where 's' is my affine parameter.

From what i've understood, solve_across_charts does use scipy.odeint which supports backwards integration but when i try to run it i get the following error:

---------------------------------------------------------------------------
ValueError                                Traceback (most recent call last)
<ipython-input-58-13e5246d5cee> in <module>
     14         geod = M.integrated_geodesic(g, (s, Integer(50), Integer(0)), v0, across_charts=True)
     15         sol = geod.solve_across_charts(step=RealNumber('0.1'),parameters_values={a:Integer(0),q:Integer(0)},
---> 16                                        verbose=True) 
     17 
     18 

/opt/sagemath-9.2/local/lib/python3.7/site-packages/sage/manifolds/differentiable/integrated_curve.py in solve_across_charts(self, charts, step, solution_key, parameters_values, verbose, **control_param)
   1719         ics = initial_pt_coords + initial_tgt_vec_comps
   1720         times = np.linspace(t_min, t_max, int((t_max - t_min) / step) + 1,
-> 1721                             endpoint=True)
   1722         nt = len(times)
   1723 

<__array_function__ internals> in linspace(*args, **kwargs)

/opt/sagemath-9.2/local/lib/python3.7/site-packages/numpy/core/function_base.py in linspace(start, stop, num, endpoint, retstep, dtype, axis)
    113     num = operator.index(num)
    114     if num < 0:
--> 115         raise ValueError("Number of samples, %s, must be non-negative." % num)
    116     div = (num - 1) if endpoint else num
    117 

ValueError: Number of samples, -499, must be non-negative.

Which does make sense. Any workaround? I've been also trying to set my affine parameter as a negative but this also raises an error.

Thanks in advance.

It does seem to do the trick for the time being. Simulations still running but so far so good. Thanks a lot!

Integrating geodesics backwards in time Is it possible to integrate a geodesic equation backwards in time in sagemath ?