8
$\begingroup$

This question can be viewed as a follow-up of

What is the definition of Curl in Mathematica?

First argument of Grad can be an array, but what definition does Mathematica use for gradient of an array? Document of Grad is rather brief on this topic. The following seems to be the only related paragraph in Details section:

In $\text{Grad}[f,{x_1,…,x_n},chart]$, if $f$ is an array, it must have dimensions ${n,…,n}$. The components of $f$ are interpreted as being in the orthonormal basis associated with $chart$.

In Properties & Relations section the following hint is included:

The gradient of an array equals the gradient of its components only in Cartesian coordinates.

And I can't find anything more in the document.

After reading some material about tensor (I found The Poor Man’s Introduction to Tensors a good one BTW), I guess it's something related to covariant derivative, but I fail to go further.

$\endgroup$
  • 2
    $\begingroup$ Oh, didn't notice this is 100th (undeleted) question of mine. $\endgroup$ – xzczd Jul 11 at 3:41
4
$\begingroup$

Natas' answer is almost correct, and I gave it an up vote. However, technically what Grad computes is the raised covariant derivative $\nabla^b T^{cd\ldots} = g^{ba}\partial_aT^{cd\ldots} + \Gamma^{bc}_aT^{ad\ldots} + \ldots $. The beauty of orthonormal bases, and the reason they're the only ones exposed in System` functionality, is that components are independent of raising and lowering (in Euclidean signature metrics). However, if we dig into the lower-level package and use the coordinate basis instead of the orthonormal one, you can see the difference:

Grad[
    SymbolicTensors`Tensor[
        {fr[r,θ],fθ[r,θ]},
        {SymbolicTensors`TangentBasis[{r,θ}]}
    ],
    {r,θ},
    "Polar"
]


(* SymbolicTensors`Tensor[
       {
           {Derivative[1, 0][fr][r, θ],  ((-r)*fθ[r, θ] + Derivative[0, 1][fr][r, θ])/r^2},
           {fθ[r, θ]/r + Derivative[1, 0][fθ][r, θ], fr[r, θ]/r + Derivative[0, 1][fθ][r, θ])/r^2}
       }, 
       {SymbolicTensors`TangentBasis[{r, θ}], SymbolicTensors`TangentBasis[{r, θ}]}
 ]*)

If the Grad were truly the covariant derviative, the new index would of type SymbolicTensors`CotangentBasis[{r, θ}] instead.

| improve this answer | |
$\endgroup$
  • $\begingroup$ I knew you'd come :) . I think this can be added to e.g. the Details section of document of Grad? As to the definition of raised covariant derivative, just to confirm, $\nabla^b T^{cd} = g^{ba} \nabla_a T^{cd}$, right? $\endgroup$ – xzczd Sep 12 at 2:56
  • $\begingroup$ I echo xzczd's plea that an explanation like this ought to be added somewhere in the docs for Grad[]. $\endgroup$ – J. M.'s discontentment Sep 12 at 6:01
  • 3
    $\begingroup$ @xzczd Yes, that is exactly the raised covariant derivative. We had some discussion about whether to add this, and there was some concern defining it exactly this way given how it relates to what's exposed in the system, but we batted around some ideas. I've put this on my list of pages to try to polish for 12.2. Stay tuned. $\endgroup$ – Itai Seggev Sep 18 at 21:19
7
+50
$\begingroup$

Indeed, Grad does compute the covariant derivative. This can be seen from the following example given in the documentation

In a curvilinear coordinate system, a vector with constant components may have a nonzero gradient:

Grad[{1, 1, 1}, {r, θ, ϕ}, "Spherical"]
(* {{0, -(1/r), -(1/r)}, {0, 1/r, -(Cot[θ]/r)}, {0, 0, (
  Csc[θ] (Cos[θ] + Sin[θ]))/r}} *)

Note that the result (and input) is always to be understood with respect to a particular basis. In the example above, it is my understanding that the input {1, 1, 1} is a vector

$$ \mathbf{v} = v^r \mathbf{e}_{r} + v^\theta \mathbf{e}_{\theta} + v^\phi \mathbf{e}_{\phi} $$

with $v^i = 1$ for $i = r, \theta, \phi$ and the $\mathbf{e}_i$ are orthonormal. (Note that it is commonplace in differential geometry to work with unnormalized bases.)

With respect to higher rank tensors I think that Mathematica treats all components as contravariant (upper indices), e.g. plugging in an array of rank two will be understood as

$$ \mathrm{Grad}_k A^{ij} = \partial_k A^{ij} + \Gamma^{i}_{kl} A^{lj} + \Gamma^{j}_{kl} A^{il}$$ where $\Gamma^{i}_{jk}$ are the coefficients of the Christoffel connection (with respect to the chosen basis).

Update: An example

Note furthermore that the point of the covariant derivative is to obtain an object which is a tensor (and therefore transforms as a tensor). As an example consider the following

vecGradBuiltin = Grad[{Subscript[f, 1][r, θ], Subscript[f, 2][r, θ]}, {r, θ}, "Polar"]

which is the covariant derivative of the vector with components $(f_1, f_2)$ with respect to the orthonormal basis. With respect to the coordinate basis ($\partial_r$, $\partial_\theta$) the components are $(f_1, \frac{1}{r} f_2)$ since the bases are related by

$$ \left( \begin{matrix} \mathbf{e}_r\\ \mathbf{e}_\theta \end{matrix} \right) = \left( \begin{matrix} 1 & 0 \\ 0 & \frac{1}{r} \end{matrix} \right) \left( \begin{matrix} \partial_r\\ \partial_\theta \end{matrix} \right) $$

and, of course, for two different coordinate systems

$$ \mathbf{v} = v^a \mathbf{e}_a = v'{}^a \mathbf{e}'_a $$

Using the code from (224280) to compute the Christoffel symbols, a possible implementation of a vector gradient is

ChristoffelSymbol[g_, xx_] := 
 Block[{n, ig, res}, n = Length[xx]; ig = Inverse[g];
  res = Table[(1/2)*
     Sum[ig[[i, s]]*(-D[g[[j, k]], xx[[s]]] + D[g[[j, s]], xx[[k]]] + 
         D[g[[s, k]], xx[[j]]]), {s, 1, n}], {i, 1, n}, {j, 1, n}, {k,
      1, n}];
  Simplify[res]]
vectorGrad[vec_, g_, coord_] := 
 With[{n = Length[coord], Γ = 
    ChristoffelSymbol[g, coord]}, 
  Table[D[vec[[b]], coord[[a]]] + 
    Sum[Γ[[b, a, c]] vec[[c]], {c, 1, n}], {b, 1, 
    n}, {a, 1, n}]]

with this you can compute vector gradient in the coordinate basis

g = DiagonalMatrix[{1, r^2}];
coord = {r, θ};
j = DiagonalMatrix[{1, 1/r}];
vec = Array[Subscript[f, #][r, θ]&, 2];
vecGradHomebrew = 
 vectorGrad[j.vec, g, coord]

Now you have the components of the vector gradient once with respect to the orthonormal basis (vecGradBuiltin) and the components with respect to the coordinate basis (vecGradHomebrew).

As mentioned before, the crucial point is that the vector gradient is a tensor. Therefore the components transform as a tensor. Since we know how the two bases are related to each other, one can verify (note that the vector gradient has one covariant and one contravariant index, the way it is defined here, the first index is contravariant and the second is covariant)

Inverse[j].vecGradHomebrew.j == vecGradBuiltin // Simplify
(* True *)

Update: Spherical example

It is straightforward to do this in three dimensions, e.g. with the spherical coordinate system

g = DiagonalMatrix[{1, r^2, r^2 Sin[θ]^2}];
coord = {r, θ, ϕ};
j = DiagonalMatrix[{1, 1/r, 1/(r Sin[θ])}];
vec = Array[Subscript[f, #][r, θ] &, 3];
vecGradHomebrew = vectorGrad[j.vec, g, coord]
Inverse[j].vecGradHomebrew.j == 
  Grad[vec, coord, "Spherical"] // Simplify
| improve this answer | |
$\endgroup$
  • 1
    $\begingroup$ I think $\Gamma^{j}_{kl} A^{jl}$ should be $\Gamma^{j}_{kl} A^{il}$? Also, how to extend the example to gradient of 2D array? I tried Inverse[j].Inverse[j].vecGradHomebrew.j but this doesn't seem to be correct. $\endgroup$ – xzczd Jul 10 at 5:13
  • $\begingroup$ Thanks, you are right about the indices. Note that the function vectorGrad I defined only works for vector (rank 1 contravariant tensors). (I certainly would not reimplement the general case but use one of the available packages like xAct.) I don't think that your proposed expression works. Dot will always operator on the "outermost" index. I think usually this type of index contraction is achieved by first Transposeing (with specification of the structure), then Dotting and then Transposeing back. $\endgroup$ – Natas Jul 10 at 6:48
  • 1
    $\begingroup$ Furthermore, if you want to define the components of the same rank 2 tensor with respect to the orthonormal basis and the coordinate basis you have to think about how the components should be "rescaled", i.e. the invariant object is given by $\mathbf{A} = A^{ij} \mathbf{e}_i \otimes \mathbf{e}_j$. $\endgroup$ – Natas Jul 10 at 7:06
  • $\begingroup$ I think I manage to. Hopefully it's not a coincidence. Key step: array = Table[Subscript[f, i, k][r, θ] j[[i, i]] j[[k, k]], {i, 2}, {k, 2}]; Inverse[j].Transpose[Inverse[j].Transpose@vecGradHomebrew].j == vecGradBuiltin // Simplify. (vectorGrad has been re-defined of course. ) Am I right? $\endgroup$ – xzczd Jul 10 at 7:17
  • 2
    $\begingroup$ Yes, that makes sense to me. Let me just add that for general coordinate changes (i.e. mixing of the coordinates) the matrix j is the Jacobian of the coordinate transformation (or inverse thereof, depending on your point of view) and a lot more care needs to be taken when contracting indices, i.e. you will also have to think about when to Transpose j. (And of course you should reconsider renaming the variables since now this is the "gradient" of a higher rank tensor. ;)) $\endgroup$ – Natas Jul 10 at 7:29
-2
$\begingroup$

In a Mathematica Array you need to define the function in the Array as a pure function with the symbols common to the coordinates system in use.

Grad[Array[f, 10], {x, y, z}]

{{0, 0, 0}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}, {0, 0, 0}}

f is a constant

Grad[Array[x y z &, 10], {x, y, z}]

{{y z, x z, x y}, {y z, x z, x y}, {y z, x z, x y}, {y z, x z, x y}, {y z, x z, x y}, {y z, x z, x y}, {y z, x z, x y}, {y z, x z, x y}, {y z, x z, x y}, {y z, x z, x y}}

This sticks to Function in the defition:

body& or Function[body]
is a pure (or "anonymous") function. The formal parameters are # (or #1), #2, etc. 

This works for all coordinate systems and even in relativistic problems.

| improve this answer | |
$\endgroup$
  • 3
    $\begingroup$ (-1) This is not even related to my question. Please read my question more carefully. $\endgroup$ – xzczd Jul 12 at 20:32

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.