EoM in Schwarzschild geometry: geodesic v Hamilton formalism

[Matter_Matters]

Hi there guys,

Currently writing and comparing two separate Mathematica scripts which can be found here and also here. The first one I've slightly modified to suit my needs and the second one is meant to reproduce the same results.

Both scripts are attempting to simulate the trajectory of a timelike geodesic using

1. the geodesic equations of motion.
2. The Hamiltonian formalism for equations of motion.

Now, the results should be equivelent. However, my problem is as follows: The Hamiltonian produces a timelike geodesic that falls radially inward directly to the event horizon, where the geodesic equation approach produces a beautiful precessing orbit. Now, the approach using the geodesic equation is correct, however, can anyone see why this is different in the second link above? I have no idea why they are not producing the same results!

 

[Paul Colby]

Do the two approaches produce the same equations of motion? This could be directly checked in Mathematica.

 

[Matter_Matters]

Do the two approaches produce the same equations of motion? This could be directly checked in Mathematica.

Yes the output should be complete identical however it's not! I have been messing around for a couple of hours with no success unfortunately!

 

[Paul Colby]

Yes the output should be complete identical however it's not!

Wait, same equations, same solver, same boundary conditions yield a different answer? Hum, I doubtful. Something must not be identical. My suggestion was to reduce the Hamiltonian equations of motion down to the geodesic ones and test for difference before solving them. I suspect that different equations result and the hamiltonian derived ones are not equivalent to the geodesic ones. Your problem is very likely a programming or problem setup one not one and not with the physics.

in the hamiltonian approach the momentum is given by $p = \dot q$ which is not $p = \frac{\partial L}{\partial\dot{q}}$

 

[pervect]

Wait, same equations, same solver, same boundary conditions yield a different answer? Hum, I doubtful. Something must not be identical. My suggestion was to reduce the Hamiltonian equations of motion down to the geodesic ones and test for difference before solving them. I suspect that different equations result and the hamiltonian derived ones are not equivalent to the geodesic ones. Your problem is very likely a programming or problem setup one not one and not with the physics.

in the hamiltonian approach the momentum is given by $p = \dot q$ which is not $p = \frac{\partial L}{\partial\dot{q}}$

As I recall, to derive the Hamiltonian approach, we can start out with some Lagrangian $L(q,\dot{q})$ and we DEFINE $p = \partial L / \partial \dot{q}$. The we compute the energy function $h(q,\dot{q})$, and use it to find the Hamiltonian H(p,q) with a change of variables. (Alternatively, we can define H as the Legendre transformation of L). So I'm not sure what point you're making. We write one of Hamilton's equations as $dq/dt = \partial H / \partial p$, but this isn't the same as what you wrote above. I'm not sure where you got $p = \dot{q}$, this doesn't look right. Wiki seems to bear me out, <<link>>.

 

[Paul Colby]

So I'm not sure what point you're making.

In the link provided in the OP (the second "here"), $p = \dot{q}$ which is very much not $\frac{\partial L}{\partial\dot{q}}$ at least by my reading of the Mathematica code provided. My point was perhaps this is an error.

 

[pervect]

$p = \frac{\partial S}{\partial q}$ could make sense in the Hamilton-Jacobi approach (henceforth the HJ approach), but that doesn't seem to be what the script is saying, if I'm interpreting it correctly it seems indeed saying that $p = dq/ds$. This equation isn't making any sense to me, sorry.

I'd suggest looking up the equations of motion in a text, rather than trying to reverse engineer them from a pair of Mathematica scripts, one of which isn't working correctly.

Attempting to follow my own advice, I do see a fair amount of discussion of the Hamilton-Jacobi equations of motion in my text (MTW), but nothing that mentions the "Hamiltionian" method.

I'm not sufficiently familiar with either the HJ method or Mathematica to be of any more help, unfortunately - if indeed the OP is interested in the HJ method ( something I'm not totally convinced of either).

I'm fairly familiar with the geodesic method, but that seems to be working correctly and there aren't any questions about it.

 

[Matter_Matters]

Hmm so I've got the following code. However, the output is still not corresponding to the geodesic equations approach:

[NOPARSE]

(*Symplectic schemes*)
Clear[m];
q = {t[s], r[s], θ[s], ϕ[s]};
p = {pt[s], pr[s], pθ[s], pϕ[s]};                
n = Length[q];
tt = 1 - 2 m/r[s];
rr = -1/tt;
θθ = -r[s]^2;
ϕϕ = -(r[s] Sin[θ[s]])^2;
metric = {{tt, 0, 0, 0}, {0, rr, 0, 0}, {0, 0, θθ, 0}, {0, 0, 0, ϕϕ}};
inversemetric = Simplify[Inverse[metric]];
ivs = {1.3, 0, 0, 0.088}; ics = {0, 6.5, π/2, 0};
m = 1;

(*hamiltonian=1/2 (D[q,s].inversemetric.D[q,s])-1;*)

hamiltonian = 1/2 (p.inversemetric.p) - 1;        (* define H in terms of p, q *)
pdot = Table[-D[hamiltonian, q[[i]]], {i, 1, n}];
qdot = Table[D[hamiltonian, p[[i]]], {i, 1, n}];

eqs1 = {{D[q, s] == qdot, D[p, s] == pdot},       (* adjust equations to new defs *)
  {(q /. s -> 0) == ics, (p /. s -> 0) == ivs}};

invariants = hamiltonian;
time = {s, 0, 100};                               (* change T to s to match q etc. *)
solee = NDSolve[eqs1, q, time, Method -> "ExplicitRungeKutta",
  StartingStepSize -> 0.1]
[/NOPARSE]

It's physical output is a particle that falls almost directly into the event horizon of the Schwarzschild black hole. I have no idea why it does this though.

 

[Paul Colby]

when $p$ is eliminated from the Hamiltonian system one must obtain the geodesic equation up to a possible parameter change. Once this is verified symbolically (not as a numerical solution but in terms of symbols) the answers will be the same. If they aren't symbolically equivalent then there is an error.

 

[Matter_Matters]

when $p$ is eliminated from the Hamiltonian system one must obtain the geodesic equation up to a possible parameter change. Once this is verified symbolically (not as a numerical solution but in terms of symbols) the answers will be the same. If they aren't symbolically equivalent then there is an error.

I believe the code that I posted above should be doing the job. The EoM are given by Hamilton's equations and I've edited the code in line with the above comments. Still though, I'm encountering the singularity almost at the beginning of the simulation which I can't seem to understand.

 

[Paul Colby]

Well, debugging others code is a good method to learn a language. Your inverse metric is poorly behaved at $\theta = 0$. Okay you're at $\theta = \pi/2$

 

[Paul Colby]

Okay, I get a stiff system warning at $s=10.8$ which isn't too surprising.

 

[Matter_Matters]

Okay, I get a stiff system warning at $s=10.8$ which isn't too surprising.

Paul, I don't think this physically makes sense which is how I know there is an issue somewhere. It should have a stable orbit with those IC's.

 

[Paul Colby]

Well, the $\phi(s)$ changes linearly with $s$ which is clearly wrong since angular momentum should lead to an ever increasing orbital rate. Also, there are constraints on the on the initial world velocity which seems broken. I did some rummaging around in the geodesic approach (your first "here" link) looking for the metric you are using. I couldn't find it. I agree with the form of the hamiltonian used but if the metric is non-physical all bets are off. Again, I suggest you forget solving and compare the equations after removing the $p$ dependence with the geodesic ones. My bet is this will fail. The reason I think the hamiltonian ones are broken because the $r(s)$ curve doesn't seem to be effected by changing $p\phi(0)$.

 

[pervect]

If you're just interested in the answer, see for instance see https://www.fourmilab.ch/gravitation/orbits/.- Basically, two of the geodesic equatoins for the Shwarzschild metric can be manipulated into the form.

$$E(\tau)= \left(1-\frac{2M}{r} \right) \frac{dt}{d\tau} \quad \frac{dE}{d\tau}=0 \quad L(\tau) = r^2 \, \frac{d\phi}{d\tau} \quad \frac{dL}{d\tau}=0 $$

If you expand dE and dL using the chain rule, you basically get two of the geodesic equations. You need one more equation to solve for the equations of motion, which can be taken as setting the magnitude of the 4 velocity to -1 (with a -+++ metric signature).

 

[Matter_Matters]

If you're just interested in the answer, see for instance see https://www.fourmilab.ch/gravitation/orbits/.- Basically, two of the geodesic equatoins for the Shwarzschild metric can be manipulated into the form.

$$E(\tau)= \left(1-\frac{2M}{r} \right) \frac{dt}{d\tau} \quad \frac{dE}{d\tau}=0 \quad L(\tau) = r^2 \, \frac{d\phi}{d\tau} \quad \frac{dL}{d\tau}=0 $$

If you expand dE and dL using the chain rule, you basically get two of the geodesic equations. You need one more equation to solve for the equations of motion, which can be taken as setting the magnitude of the 4 velocity to -1 (with a -+++ metric signature).

Thanks for the information. I'm trying to write it in terms of hamiltonian's specifically to later use symplectic integration schemes. However, this is proving to be quite the task! Well for my coding capabilities anyways.

 

[Paul Colby]

q = {t[s],r[s],θ[s],ϕ[s]};
p = {pt[s], pr[s], pθ[s], pϕ[s]};
tt = 1-2 M /r[s];
rr = -1/tt;
θθ = -r[s]^2;
ϕϕ = -(r[s] Sin[θ[s]])^2;
metric = {{tt,0,0,0},{0,rr,0,0},{0,0,θθ,0},{0,0,0,ϕϕ}};
inversemetric = Simplify[Inverse[metric]];
ivs = {1.3,0,0,0.088}; ics = {0,6.5, π/2, 0};
hamiltonian = 1/2 (p . inversemetric . p - 1);
M=1;
hamiltonian /. Join[Table[q[[n]]-> ics[[n]],{n,1,4}],Table[p[[n]]-> ivs[[n]],{n,1,4}]]
0.720464

The way the Hamiltonian is written its value at the initial conditions should be 0, not 0.720464. It's not a physical particle so...?

 

[Paul Colby]

For what it is worth, fixing the initial canonical momenta to be $p^\mu p_\mu = 1$ does little to change the straight into the hole issue. It is interesting that $p_\phi$ and $\theta$ are all constant which is the correct physics. I don't know what to make of $p_t$ being constant? This is clearly implied by the equations.

 

[Paul Colby]

Decreasing the mass is amusing. Stable orbits result but they are highly elliptical orbits. I get a stable orbit for $M=0.01$ but it really zings around the mass at closest approach. It's possible straight into the hole is well deserved for your initial conditions?

 

[Paul Colby]

For completeness here is a "working" code

Clear["Global`*"]
q = {t[s],r[s],θ[s],ϕ[s]};
p = {pt[s], pr[s], pθ[s], pϕ[s]};
tt = 1-2 M /r[s];
rr = -1/tt;
θθ = -r[s]^2;
ϕϕ = -(r[s] Sin[θ[s]])^2;
metric = {{tt,0,0,0},{0,rr,0,0},{0,0,θθ,0},{0,0,0,ϕϕ}};
inversemetric = Simplify[Inverse[metric]];
ivs = {1,0,0,0.088}; ics = {0,6.5, π/2, 0};
hamiltonian = 1/2 (p . inversemetric . p - 1);
M=0.01
0.01
norm = Sqrt[p . inversemetric . p /.
   Join[Table[q[[n]]-> ics[[n]],{n,1,4}],
        Table[p[[n]]-> ivs[[n]],{n,1,4}]]]
1.00145
ivs = ivs / norm;
ivs
{0.998552,0.,0.,0.0878725}
hamiltonian /.
   Join[Table[q[[n]]-> ics[[n]],{n,1,4}],
        Table[p[[n]]-> ivs[[n]],{n,1,4}]]
1.06604*10^-16
qdot = Table[D[hamiltonian,p[[n]]],{n,1,4}];
pdot = Table[-D[hamiltonian,q[[n]]],{n,1,4}];
qequs = Table[D[q[[n]],s]==qdot[[n]],{n,1,4}];
qbcs = Table[(q[[n]]/. s-> 0) == ics[[n]],{n,1,4}];
pbcs = Table[(p[[n]]/. s-> 0) == ivs[[n]],{n,1,4}];
pequs = Table[D[p[[n]],s]==pdot[[n]],{n,1,4}];
eqs1 = Join[qequs,pequs,qbcs, pbcs];
time = {s,0,1000};
pbcs
{pt(0)==0.998552,pr(0)==0.,pθ(0)==0.,pϕ(0)==0.0878725}
sol = NDSolve[eqs1,Join[p,q],time,StartingStepSize->0.01];

 

[vanhees71]

Let's start in clear formulae not some mediocre computer code. The most simple method is to use the quadratic form of the action principle, i.e.,
$$L=\frac{1}{2} g_{\mu \nu} \dot{x}^{\mu} \dot{x}^{\nu}.$$
It automatically leads to the parametrization with an affine parameter since
$$H=p_{\mu} \dot{x}^{\mu}-L=L$$
is conserved, i.e.,
$$g_{\mu \nu} \dot{x}^{\mu} \dot{x}^{\nu}=\text{const}.$$
A massive particle has a time-like trajectory, so that you can choose the overall scale of the parameter such that
$$\text{const}=1,$$
and then the affine parameter is the proper time $\tau$ of the particle.

Further you have
$$p_{\mu}=\frac{\partial L}{\partial \dot{x}^{\mu}}=g_{\mu \nu} \dot{x}^{\nu},$$
which implies that the Hamiltonian in terms of the canonical momenta (as is needed for the Hamilton canonical equations) is
$$H=\frac{1}{2} g^{\mu \nu} p_{\mu} p_{\nu}.$$
Now let's check the equations of motion

(a) Lagrange formalism

$$\dot{p}_{\mu}=\frac{\mathrm{d}}{\mathrm{d} \lambda} [g_{\mu \nu} \dot{x}^{\nu}]=g_{\mu \nu} \dot{x}^{\nu} + \partial_{\alpha} g_{\mu \nu} \dot{x}^{\alpha} \dot{x}^{\alpha} = \frac{\partial L}{\partial x^{\mu}}=\frac{1}{2} \partial_{\mu} g_{\alpha \nu} \dot{x}^{\alpha} \dot{x}^{\nu},$$
or in the more familiar form
$$g_{\mu \nu} \ddot{x}^{\nu}+\Gamma_{\mu \alpha \nu} \dot{x}^{\alpha} \dot{x}^{\nu},$$
where
$$\Gamma_{\mu \alpha \nu}=\frac{1}{2} (\partial_{\alpha} g_{\mu \nu} + \partial_{\nu} g_{\mu \alpha}-\partial_{\mu} g_{\alpha \nu})=g_{\mu \rho} {\Gamma^{\rho}}_{\alpha \nu},$$
with the usual Christoffel symbols for the Lorentzian spacetime manifold.

(b) Hamilton formalism

Here the equations of motion are given by the Hamilton canonical equations,
$$\dot{x}^{\mu}=\frac{\partial H}{\partial p_{\mu}}=g^{\mu \nu} p_{\nu}, \quad \dot{p}_{\mu}=-\frac{\partial H}{\partial x^{\mu}}=-\frac{1}{2} \partial_{\mu} g^{\alpha \nu} p_{\alpha} p_{\nu}.$$
From the first equation you get
$$p_{\mu} = g_{\mu \nu} \dot{x}^{\nu}.$$
Taking the derivative leads to
$$\dot{p}_{\mu}=g_{\mu \nu} \ddot{x}^{\nu}+\partial_{\alpha} g_{\mu \nu} \dot{x}^{\mu} \dot{x}^{\alpha}=-\frac{1}{2} (\partial_{\mu} g^{\alpha \nu}) g_{\alpha \rho} g_{\nu \sigma} \dot{x}^{\rho} \dot{x}^{\sigma}=+\frac{1}{2} g^{\alpha \nu} \partial_{\mu} g_{\alpha \rho} g_{\nu \sigma} \dot{x}^{\rho} \dot{x}^{\sigma}=\frac{1}{2}\partial_{\mu} g_{\alpha \rho} \dot{x}^{\alpha} \dot{x}^{\rho} .$$
This is obviously the same equation we derived above within the Lagrange formalism. That's no surprise, because the Lagrange and Hamilton formalisms derive the equations of motion from the same action principle and thus lead the same equations of motion.

 

[Paul Colby]

That's not surprise, because the Lagrange and Hamilton formalisms derive the equations of motion from the same action principle and thus lead the same equations of motion.

Totally true. I use this argument all the time in discussions with quality assurance about code bugs. Never seems to work.

 

[pervect]

Let's start in clear formulae not some mediocre computer code. The most simple method is to use the quadratic form of the action principle, i.e.,
$$L=\frac{1}{2} g_{\mu \nu} \dot{x}^{\mu} \dot{x}^{\nu}.$$

Possibly I'm confused, but something seems wrong here.

Let's take the flat space case when $g_{\mu\nu} = \eta_{\mu\nu}$ = diag(-1,1,1,1), so we have a -+++ metric signature.

We know that the relativistic Lagrangian for a free particle in the above flat space-time should be $L = -m\,c^2 \sqrt{1-\beta^2}$, where $\beta = v/c$. If we use geometric units so that c=1, we can write this as:

$L(t,\textbf{x},\frac{d \textbf{x}}{dt}) = -m \, \sqrt{ -g_{\mu\nu} \frac{dx^\mu}{dt} \frac{dx^\nu}{dt}}$

where $\textbf{x}$ is a spatial vector, and the implied sum is only over the spatial terms. I am using the non-covariant formalism of the Lagrangian, but when I look up the covariant formalism of the Lagrangian in Goldstein, I still see square roots.

I don't see how to reconcile this with what you wrote. I've seen in the past Einbein forrmalisms that eliminate the square root from my version of the Lagrangian above by introducing another variable, but if that's what's going on, a bit more explanation would be quite helpful.

Another way of putting this -a change in the action $\Delta S = L \Delta t$ should be proportional to change in proper time $\Delta \tau$, so that extremizing the action of a free particle is the same as extremizing proper time, and in flat space time this means that L is proportional to $\frac{d\tau}{dt} = 1/\gamma = \sqrt{1-\beta^2}$. And we know that even in curved space-time, test particles move along worldlines that extremize proper time. But I don't quite see how that applies to what you (vanhees71) wrote.

 

[PeterDonis]

where $\textbf{x}$ is a spatial vector, and the implied sum is only over the spatial terms

I think you left out a $1$ under the square root.

when I look up the covariant formalism of the Lagrangian in Goldstein, I still see square roots

Yes, you just do $g_{\mu \nu} U^\mu U^\nu$ under the square root where $U^\mu$ is now the 4-velocity vector, i.e., a 4-vector, not a spatial 3-vector. (And of course there will be no $1$ under the square root this time.)

 

[pervect]

I think you left out a $1$ under the square root.
Yes, you just do $g_{\mu \nu} U^\mu U^\nu$ under the square root where $U^\mu$ is now the 4-velocity vector, i.e., a 4-vector, not a spatial 3-vector. (And of course there will be no $1$ under the square root this time.)

Good points - I edited my post to just strike out that portion of my response - and fix up a few other things that bothered me on re-reading. I still think there's a square root missing from vanhees71 version of the Lagrangian.

 

[PeterDonis]

I still think there's a square root missing from vanhees71 version of the Lagrangian.

I don't think so; there is actually a quadratic form of the Lagrangian. I have seen it discussed in QFT texts; for example, Zee's Quantum Field Theory in a Nutshell does so.

 

[Paul Colby]

It seems clear (ha right) that if a curve is an extremum of $s$ that it is also an extremum of $s^2$

 

[Paul Colby]

I guess I take issue with what I said in #27. An extremum of $L(q,\dot{q})$ is found,

$\frac{d}{dt}(\frac{\partial L}{\partial \dot{q}}) - \frac{\partial L}{\partial q} = 0$
​

right? Now let $L = \frac{1}{2}A^2$. The equation becomes,

$\frac{d}{dt}(A\frac{\partial A}{\partial \dot{q}}) - A\frac{\partial A}{\partial q} = 0$
​

doesn't look to be the same equation. Still as shown #21 the equations are the same in this case.

 

[vanhees71]

Possibly I'm confused, but something seems wrong here.

Let's take the flat space case when $g_{\mu\nu} = \eta_{\mu\nu}$ = diag(-1,1,1,1), so we have a -+++ metric signature.

We know that the relativistic Lagrangian for a free particle in the above flat space-time should be $L = -m\,c^2 \sqrt{1-\beta^2}$, where $\beta = v/c$. If we use geometric units so that c=1, we can write this as:

$L(t,\textbf{x},\frac{d \textbf{x}}{dt}) = -m \, \sqrt{ -g_{\mu\nu} \frac{dx^\mu}{dt} \frac{dx^\nu}{dt}}$

where $\textbf{x}$ is a spatial vector, and the implied sum is only over the spatial terms. I am using the non-covariant formalism of the Lagrangian, but when I look up the covariant formalism of the Lagrangian in Goldstein, I still see square roots.

I don't see how to reconcile this with what you wrote. I've seen in the past Einbein forrmalisms that eliminate the square root from my version of the Lagrangian above by introducing another variable, but if that's what's going on, a bit more explanation would be quite helpful.

Another way of putting this -a change in the action $\Delta S = L \Delta t$ should be proportional to change in proper time $\Delta \tau$, so that extremizing the action of a free particle is the same as extremizing proper time, and in flat space time this means that L is proportional to $\frac{d\tau}{dt} = 1/\gamma = \sqrt{1-\beta^2}$. And we know that even in curved space-time, test particles move along worldlines that extremize proper time. But I don't quite see how that applies to what you (vanhees71) wrote.

If you want the 3+1-formulation, an analysis using Noether's theorem indeed leads to the square-root form (let me use my HEP convenction with signatur (+---) since otherwise I get confused myself). This is the parametrization-invariant formulation, i.e., you can use any parametrization of the particles' worldline you like. That's the case if the Lagrangian is homogeneous of 1st degree in the $\dot{x}^{\mu}$. Then, of course you get $H \equiv 0$, independent of whether you insert the solutions of the equations of motion or not, i.e., you cannot use the Hamilton formalism easily. The reason is indeed the reparametrization invariance.

Of course, you can also use affine parameters like, in the case of a massive particle, proper time. Then you envoke the constraint
$$g_{\mu \nu} \dot{x}^{\mu} \dot{x}^{\nu}=c^2$$
to choose the proper time as world-line parameter. Working this idea out, using the Lagrange-multiplier method leads to the square form of the Lagrangian, which I used, because it's particularly simple, leading to a nice Hamilton formalism without problems related to parametrization invariance. For details see

http://th.physik.uni-frankfurt.de/~hees/pf-faq/srt.pdf

p. 42 ff.

This can be extended to GR in a straight-forward way by substituting $g_{\mu \nu}$ for $\eta_{\mu \nu}$ in the action.

 

[vanhees71]

I guess I take issue with what I said in #27. An extremum of $L(q,\dot{q})$ is found,

$\frac{d}{dt}(\frac{\partial L}{\partial \dot{q}}) - \frac{\partial L}{\partial q} = 0$
​

right? Now let $L = \frac{1}{2}A^2$. The equation becomes,

$\frac{d}{dt}(A\frac{\partial A}{\partial \dot{q}}) - A\frac{\partial A}{\partial q} = 0$
​

doesn't look to be the same equation. Still as shown #21 the equations are the same in this case.

It becomes the same equation, given that for an affine parameter $A=\text{const}$. That you have an affine parameter is worked into the square-form Lagrangian, because it implies that $H=L=\text{const}$, because $L$ doesn't explicitly depend on the parameter (which you can choose, e.g., as the proper time in the case of a massive particle).

 

[Paul Colby]

It becomes the same equation, given that for an affine parameter A=constA=constA=\text{const}.

Ah, the missing fact, thanks.

 

[stevendaryl]

Possibly I'm confused, but something seems wrong here.

Let's take the flat space case when $g_{\mu\nu} = \eta_{\mu\nu}$ = diag(-1,1,1,1), so we have a -+++ metric signature.

We know that the relativistic Lagrangian for a free particle in the above flat space-time should be $L = -m\,c^2 \sqrt{1-\beta^2}$, where $\beta = v/c$. If we use geometric units so that c=1, we can write this as:

$L(t,\textbf{x},\frac{d \textbf{x}}{dt}) = -m \, \sqrt{ -g_{\mu\nu} \frac{dx^\mu}{dt} \frac{dx^\nu}{dt}}$

It's one of those lucky miracles that you get the same equations of motion using:

$L = -m \sqrt{g_{\mu \nu} \frac{dx^\mu}{ds} \frac{dx^\nu}{ds}}$

or

$L =\frac{m}{2} g_{\mu \nu} \frac{dx^\mu}{ds} \frac{dx^\nu}{ds}$

(Actually, the first gives more solutions than the second, but every solution to the first can be converted into a solution of the second by reparametrizing $s$)