TODO rename to macroscopic electromagnetism. Is this only macroscopic? Does it exist as a property of spacetime via Kaluza-Klein (esp if we relax the cylindrical constraint)?
(from https://en.wikipedia.org/wiki/Covariant_formulation_of_classical_electromagnetism, but adjusting sign convention)
Let $D^i = \epsilon_0 E^i + P^i$
Units of $D_i$ and $P_i$ are $\frac{C^2 \cdot s^2}{kg \cdot m^3} \frac{kg \cdot m}{C \cdot s^2} = \frac{C}{m^2}$
Let $H_i = \frac{1}{\mu_0} B_i - M_i$
Units of $H_i$ and $M_i$ are in $\frac{A}{m} = \frac{C}{m \cdot s}$
$P^i = \epsilon_0 \chi_e E^i$
For $\chi_e = $ material's electric volume susceptibillity.
$M_i = \chi_m H_i$
For $\chi_m =$ material's magnetic volume susceptibillity.
$D^i = \epsilon_0 E^i + P^i = \epsilon_0 (1 + \chi_e) E^i = \epsilon E^i$
For $\epsilon = $ electric permittivity
So $\frac{\epsilon}{\epsilon_0} = \epsilon_r = 1 + \chi_e$
Where $\epsilon_r$ is the relative permittivity, which is unitless
$H_i = \frac{1}{\mu_0} B_i - \chi_m H_i$
$\mu_0 (1 + \chi_m) H_i = \mu H_i = B_i$
For $\mu = $ magnetic permeability
So $\frac{\mu}{\mu_0} = \mu_r = 1 + \chi_m$
Where $\mu_r$ is the relative permeability, which is unitless
Notice that $\chi_e, \chi_m, \epsilon, \mu, \epsilon_r, \mu_r$ are operators, not necessarily scalars, not necessarily real
Also:
Let $v_p$ be the phase velocity of light such that $v_p \le c, \frac{1}{(v_p)^2} = \mu \epsilon$, in units $[\frac{m}{s}]$
$\frac{1}{\mu_0} F_{ab} = \downarrow a(i) \overset{\rightarrow b(j)}{\pmatrix{
0 & -\epsilon_0 c E_j \\
\epsilon_0 c E_i & \frac{1}{\mu_0} {\bar\epsilon_{ij}}^k B_k
}}$
Let $M_{uv} = \pmatrix{ 0 & c P_j \\ -c P_i & {\bar\epsilon_{ij}}^k M_k} = \pmatrix{
0 & c \epsilon_0 \chi_e E_j \\
-c \epsilon_0 \chi_e E_i & {\bar\epsilon_{ij}}^k \chi_m \mu^{-1} B_k
}$
$M_{uv}$ has units $\frac{C}{m \cdot s}$
Let $D_{uv} = \pmatrix{ 0 & -c D_j \\ c D_i & {\bar\epsilon_{ij}}^k H_k} = \pmatrix{
0 & -c \epsilon E_j \\
c \epsilon E_i & {\bar\epsilon_{ij}}^k \mu^{-1} B_k
}$
$D_{uv}$ has units $\frac{C}{m \cdot s}$
So $D_{uv} = \frac{1}{\mu_0} F_{uv} - M_{uv}$
$(J_{free})^u = {D^{uv}}_{,v}$
$(J_{free})^0 = {D^{0u}}_{,u} = c D_{j,j} = c \rho_{free}$
$(J_{free})^i = -D_{i,0} + \bar\epsilon^{ijk} H_{k,j}$
$(J_{free})^u$ is in units of $\frac{C}{m^2 \cdot s}$
${(J_{free})^u}_{,u} = D_{j,j0} - D_{i,0i} + \bar\epsilon^{ijk} H_{k,ji} = 0$
${(J_{free})^u}_{,u}$ is in units of $\frac{C}{m^3 \cdot s}$
$(J_{bound})^u = {M^{uv}}_{,v}$
$(J_{bound})^u$ is in units of $\frac{C}{m^2 \cdot s}$
$(J_{bound})^0 = {M^{0v}}_{,v} = -c P_{j,j} = c \rho_{bound}$
$(J_{bound})^i = P_{i,0} + \bar\epsilon^{ijk} M_{k,j}$
$(J_{bound})^u$ is in units of $\frac{C}{m^2 \cdot s}$
${(J_{bound})^u}_{,u} = -P_{j,j0} + P_{i,i0} + \bar\epsilon^{ijk} M_{k,ji} = 0$
${(J_{bound})^u}_{,u}$ is in units of $\frac{C}{m^3 \cdot s}$
$J^a = (J_{free})^a + (J_{bound})^a$ is in units of $\frac{C}{m^3 \cdot s}$
$(J_{free})^a = \pmatrix{ c \rho_{free} \\ {J^i}_{free} } = \pmatrix{ c D_{j,j} \\ -D_{i,0} + {\bar\epsilon_i}^{jk} H_{k,j} }$
$(J_{bound})^a = \pmatrix{ c \rho_{bound} \\ {J^i}_{bound} } = \pmatrix{ -c P_{j,j} \\ P_{i,0} + {\bar\epsilon_i}^{jk} M_{k,j} }$
For linear susceptibility:
${J^i}_{free} = \sigma E^i$
As an exterior derivative in special relativity:
$\frac{1}{2} \delta D = \frac{1}{2} \star d \star D = J$
$J_u dx^u = \frac{1}{2} \partial_a (\star D)_{bc} \star (dx^a \wedge dx^b \wedge dx^c)$
becomes:
$J_0 dx^0 = \frac{1}{2} ((\star D)_{23,1} + (\star D)_{31,2} + (\star D)_{12,3} - (\star D)_{32,1} - (\star D)_{13,2} - (\star D)_{21,3}) \star (dx^1 \wedge dx^2 \wedge dx^3)$
$J_1 dx^1 = \frac{1}{2} ((\star D)_{23,0} + (\star D)_{30,2} + (\star D)_{02,3} - (\star D)_{32,0} - (\star D)_{03,2} - (\star D)_{20,3}) \star (dx^0 \wedge dx^2 \wedge dx^3)$
$J_2 dx^2 = \frac{1}{2} ((\star D)_{13,0} + (\star D)_{30,1} + (\star D)_{01,3} - (\star D)_{31,0} - (\star D)_{03,1} - (\star D)_{10,3}) \star (dx^0 \wedge dx^1 \wedge dx^3)$
$J_3 dx^3 = \frac{1}{2} ((\star D)_{12,0} + (\star D)_{20,1} + (\star D)_{01,2} - (\star D)_{21,0} - (\star D)_{02,1} - (\star D)_{10,2}) \star (dx^0 \wedge dx^1 \wedge dx^2)$
becomes:
$J_0 dx^0 = \frac{1}{2} ((\star D)_{23,1} + (\star D)_{31,2} + (\star D)_{12,3} - (\star D)_{32,1} - (\star D)_{13,2} - (\star D)_{21,3}) (-dx^0)$
$J_1 dx^1 = \frac{1}{2} ((\star D)_{23,0} + (\star D)_{30,2} + (\star D)_{02,3} - (\star D)_{32,0} - (\star D)_{03,2} - (\star D)_{20,3}) dx^1$
$J_2 dx^2 = \frac{1}{2} ((\star D)_{13,0} + (\star D)_{30,1} + (\star D)_{01,3} - (\star D)_{31,0} - (\star D)_{03,1} - (\star D)_{10,3}) (-dx^2)$
$J_3 dx^3 = \frac{1}{2} ((\star D)_{12,0} + (\star D)_{20,1} + (\star D)_{01,2} - (\star D)_{21,0} - (\star D)_{02,1} - (\star D)_{10,2}) dx^3$
becomes:
$c \rho_{free} dx^0 = c (D_{1,1} + D_{2,2} + D_{3,3}) dx^0$
$J_1 dx^1 = ( -c D_{1,0} + H_{3,2} - H_{2,3}) dx^1$
$J_2 dx^2 = ( -c D_{2,0} + H_{1,3} - H_{3,1}) dx^2$
$J_3 dx^3 = ( -c D_{3,0} + H_{2,1} - H_{1,2}) dx^3$
becomes:
$D_{1,1} + D_{2,2} + D_{3,3} = \rho_{free}$
$D_{1,t} - H_{3,2} + H_{2,3} = -J_1$
$D_{2,t} - H_{1,3} + H_{3,1} = -J_2$
$D_{3,t} - H_{2,1} + H_{1,2} = -J_3$
special relativity electromagnetic stress-energy in medium:
$T_{uv} = F_{ua} {D_v}^a - \frac{1}{4} g_{uv} F_{ab} D^{ab}$
...has units $\frac{kg}{C \cdot s} \cdot \frac{C}{m \cdot s} = \frac{kg}{m \cdot s^2}$
${T^u}_v = -({F^u}_a {D^a}_v)^{TF}$
${F^u}_v = \downarrow u(i) \overset{\rightarrow v(j)}{\pmatrix{
0 & \frac{1}{c} E_j \\
\frac{1}{c} E_i & {\bar\epsilon_{ij}}^k B_k
}}$
${D^u}_v = \downarrow u(i) \overset{\rightarrow v(j)}{\pmatrix{
0 & c D_j \\
c D_i & {\bar\epsilon_{ij}}^k H_k
}}$
${F^u}_a {D^a}_v = \downarrow u(i) \overset{\rightarrow v(j)}{\pmatrix{
E_k D^k &
\frac{1}{c} E_k {\bar\epsilon_{kj}}^l H_l \\
c {\bar\epsilon_{ik}}^l D^k B_l &
E_k D^k + {\bar\epsilon_{ik}}^l B_l {\bar\epsilon_{kj}}^m H_m
}}$
${F^u}_a {D^a}_v = \downarrow u(i) \overset{\rightarrow v(j)}{\pmatrix{
E_k \epsilon E^k &
-\frac{1}{c} {\bar\epsilon_{jk}}^l E_k \mu^{-1} B_l \\
c {\bar\epsilon_{ik}}^l \epsilon E^k B_l &
E_j \epsilon E^i + B_j \mu^{-1} B_i - \delta^i_j B_m \mu^{-1} B_m
}}$
$tr (F \cdot D) = 2 (E_k \epsilon E^k - B_k \mu^{-1} B^k)$
$({F^u}_a {D^a}_v)^{TF} = {F^u}_a {D^a}_v - \frac{1}{4} \delta^u_v {F^b}_a {D^a}_b$
$= \downarrow u(i) \overset{\rightarrow v(j)}{\pmatrix{
\frac{1}{2} (E_k \epsilon E^k + B_k \mu^{-1} B^k) &
-\frac{1}{c} {\bar\epsilon_{jk}}^l E_k \mu^{-1} B_l \\
c {\bar\epsilon_{ik}}^l \epsilon E^k B_l &
E_j \epsilon E^i + B_j \mu^{-1} B_i - \frac{1}{2} \delta^i_j (E_k \epsilon E^k + B_k \mu^{-1} B^k)
}}$
${T^u}_v = -{((F \cdot D)^{TF})^u}_v = \downarrow u(i) \overset{\rightarrow v(j)}{\pmatrix{
-\frac{1}{2} (E_k \epsilon E^k + B_k \mu^{-1} B^k) &
\frac{1}{c} {\bar\epsilon_{jk}}^l E_k \mu^{-1} B_l \\
-c {\bar\epsilon_{ik}}^l \epsilon E^k B_l &
\frac{1}{2} \delta^i_j (E_k \epsilon E^k + B_k \mu^{-1} B^k) - E_j \epsilon E^i - B_j \mu^{-1} B_i
}}$
$T^{uv} = -((F \cdot D)^{TF})^{uv} = \downarrow u(i) \overset{\rightarrow v(j)}{\pmatrix{
\frac{1}{2} (E_k \epsilon E^k + B_k \mu^{-1} B^k) &
\frac{1}{c} {\bar\epsilon_{jk}}^l E_k \mu^{-1} B_l \\
c {\bar\epsilon_{ik}}^l \epsilon E^k B_l &
\frac{1}{2} \eta^{ij} (E_k \epsilon E^k + B_k \mu^{-1} B^k) - E_j \epsilon E^i - B_j \mu^{-1} B_i
}}$
In the case that $M_{uv} = 0$ we find the vacuum stress-energy tensor, which is symmetric.
In terms of $B$ and $D = \epsilon E$:
$T^{uv} = \downarrow u(i) \overset{\rightarrow v(j)}{\pmatrix{
\frac{1}{2} (D_k \epsilon^{-1} D^k + B_k \mu^{-1} B^k) &
\frac{1}{c} {\bar\epsilon_{jk}}^l \epsilon^{-1} D_k \mu^{-1} B_l \\
c {\bar\epsilon_{ik}}^l D^k B_l &
\frac{1}{2} \eta^{ij} (D_k \epsilon^{-1} D^k + B_k \mu^{-1} B^k) - D_j \epsilon^{-1} D^i - B_j \mu^{-1} B_i
}}$
There, looks better. $D$ and $B$ match up with $\epsilon$ and $\mu$
What is $T^{uv}$, assuming that $\epsilon$ and $\mu$ are scalars?
$T^{uv} = \downarrow u(i) \overset{\rightarrow v(j)}{\pmatrix{
\frac{1}{2} (\epsilon E^2 + \frac{1}{\mu} B^2) &
\frac{1}{c \mu} \bar\epsilon_{jkl} E_k B_l \\
c \epsilon \bar\epsilon_{ikl} E_k B_l &
\frac{1}{2} \eta^{ij} (\epsilon E^2 + \frac{1}{\mu} B^2) - \epsilon E^i E^j - \frac{1}{\mu} B^i B^j
}}$
Scalar susceptibility in materials, with $D$ instead of $E$:
$T^{uv} = \downarrow u(i) \overset{\rightarrow v(j)}{\pmatrix{
\frac{1}{2} (\frac{1}{\epsilon} D^2 + \frac{1}{\mu} B^2) &
c (\frac{v_p}{c})^2 \bar\epsilon_{jkl} D_k B_l \\
c \bar\epsilon_{ikl} D_k B_l &
\frac{1}{2} \eta^{ij} (\frac{1}{\epsilon} D^2 + \frac{1}{\mu} B^2) - \frac{1}{\epsilon} D^i D^j - \frac{1}{\mu} B^i B^j
}}$
$T^{[uv]} = \downarrow u(i) \overset{\rightarrow v(j)}{\pmatrix{
0 &
-\frac{1}{2} (1 - (\frac{v_p}{c})^2) c \bar\epsilon_{jkl} D_k B_l
\\
\frac{1}{2} (1 - (\frac{v_p}{c})^2) c \bar\epsilon_{ikl} D_k B_l
&
0
}}$
Now to simplify it a bit...
$\frac{1}{2} (1 - (\frac{v_p}{c})^2) c \bar\epsilon_{ikl} D_k B_l$
Let $n = \frac{c}{v_p}$ be the index of refraction.
$= \frac{1}{2} (1 - \frac{1}{n^2}) c \bar\epsilon_{ikl} D_k B_l$
For most materials (except metamaterials), we have the index of refraction $n \ge 1$, so $\frac{1}{n} \in (0,1]$, so $1 - \frac{1}{n^2} \in [0,1)$.
Also, $\vec{D} \times \vec{B} = \epsilon \vec{E} \times \mu \vec{H} = \frac{1}{(v_p)^2} \vec{S} = \frac{1}{c^2} n^2 \vec{S}$, for $\vec{S}$ the Poynting vector.
$= \frac{1}{2} \frac{1}{c} (n^2 - 1) S_i$
$T^{[uv]} = \downarrow u(i) \overset{\rightarrow v(j)}{\pmatrix{
0 & -\frac{1}{2} \frac{1}{c} (n^2 - 1) S_j \\
\frac{1}{2} \frac{1}{c} (n^2 - 1) S_i & 0
}}$
$\vec{D} \times \vec{B}$ has units $\frac{kg}{m^2 \cdot s}$
$\vec{E} \times \vec{H} = \vec{S}$ has units $\frac{kg}{s^3}$
$T^{[i0]} = \frac{1}{2} \frac{1}{c} (n^2 - 1) S_i$ has units $\frac{kg}{m \cdot s^2}$
For curiousity's sake, what is the symmetric component?
$T^{(uv)} = \downarrow u(i) \overset{\rightarrow v(j)}{\pmatrix{
\frac{1}{2} (\frac{1}{\epsilon} D^2 + \frac{1}{\mu} B^2) &
\frac{1}{2} (c (\frac{v_p}{c})^2 \bar\epsilon_{jkl} D_k B_l + c \bar\epsilon_{jkl} D_k B_l) \\
\frac{1}{2} (c (\frac{v_p}{c})^2 \bar\epsilon_{ikl} D_k B_l + c \bar\epsilon_{ikl} D_k B_l) &
\frac{1}{2} \eta^{ij} (\frac{1}{\epsilon} D^2 + \frac{1}{\mu} B^2) - \frac{1}{\epsilon} D^i D^j - \frac{1}{\mu} B^i B^j
}}$
Simplifying...
$T^{(i0)} = \frac{1}{2} (c (\frac{v_p}{c})^2 \bar\epsilon_{ikl} D_k B_l + c \bar\epsilon_{ikl} D_k B_l)$
$= \frac{1}{2} c (1 + (\frac{v_p}{c})^2) \bar\epsilon_{ikl} D_k B_l$
$= \frac{1}{2} \frac{1}{c} (n^2 + 1) S_i$
$T^{(uv)} = \downarrow u(i) \overset{\rightarrow v(j)}{\pmatrix{
\frac{1}{2} (\frac{1}{\epsilon} D^2 + \frac{1}{\mu} B^2) &
\frac{1}{2} \frac{1}{c} (n^2 + 1) S_j \\
\frac{1}{2} \frac{1}{c} (n^2 + 1) S_i &
\frac{1}{2} \eta^{ij} (\frac{1}{\epsilon} D^2 + \frac{1}{\mu} B^2) - \frac{1}{\epsilon} D^i D^j - \frac{1}{\mu} B^i B^j
}}$
What is the difference between $T^{uv}$ and $T^{uv}_{vacuum}$, assuming $\epsilon$ and $\mu$ are scalars?
(I'm thinking of redoing this or getting rid of it...)
$T^{uv} - T^{uv}_{vacuum} = \pmatrix{
\frac{1}{2} (\epsilon_0 (1 + \chi_e) E^2 + \frac{1}{\mu_0 (1 + \chi_m)} B^2) &
\frac{1}{c \mu_0 (1 + \chi_m)} \bar\epsilon_{jkl} E_k B_l \\
c \epsilon_0 (1 + \chi_e) \bar\epsilon_{ikl} E_k B_l &
\frac{1}{2} \eta^{ij} (\epsilon_0 (1 + \chi_e) E^2 + \frac{1}{\mu_0 (1 + \chi_m)} B^2) - \epsilon_0 (1 + \chi_e) E^i E^j - \frac{1}{\mu_0 (1 + \chi_m)} B^i B^j
}$ $- \pmatrix{
\frac{1}{2} (\epsilon_0 E^2 + \frac{1}{\mu_0} B^2) &
\frac{1}{c \mu_0} \bar\epsilon_{jkl} E_k B_l \\
c \epsilon_0 \bar\epsilon_{ikl} E_k B_l &
\frac{1}{2} \eta^{ij} (\epsilon_0 E^2 + \frac{1}{\mu_0} B^2) - \epsilon_0 E^i E^j - \frac{1}{\mu_0} B^i B^j
}$
$= \pmatrix{
\frac{1}{2} (\epsilon_0 \chi_e E^2 - \frac{\chi_m}{\mu} B^2) &
-\frac{1}{c} \frac{\chi_m}{\mu} \bar\epsilon_{jkl} E_k B_l \\
-\frac{1}{c} \frac{\chi_m}{\mu} \bar\epsilon_{ikl} E_k B_l &
\frac{1}{2} (\epsilon_0 \chi_e E^2 - \frac{\chi_m}{\mu} B^2) - \epsilon_0 \chi_e E^i E^j + \frac{\chi_m}{\mu} B^i B^j
}$
Notice the difference is still symmetric. So for scalar $\epsilon, \mu$ we still have a symmetric stress-energy tensor, which means no spacetime torsion.
stress-energy divergence of electromagnetism in medium with scalar susceptibility in special relativity:
$T^{uv} = \pmatrix{
\frac{1}{2} (\frac{1}{\epsilon} D^2 + \frac{1}{\mu} B^2) &
c (\frac{v_p}{c})^2 \bar\epsilon_{jkl} D_k B_l \\
c \bar\epsilon_{ikl} D_k B_l &
\frac{1}{2} \eta^{ij} (\frac{1}{\epsilon} D^2 + \frac{1}{\mu} B^2) - \frac{1}{\epsilon} D^i D^j - \frac{1}{\mu} B^i B^j
}$
And now for ${T^{uv}}_{,u}$...
${T^{uv}}_{,v} = $...
${T^{0v}}_{,v} = \frac{1}{2} \frac{1}{c} (\frac{1}{\epsilon} D^2 + \frac{1}{\mu} B^2)_{,0}
+ \frac{1}{c} \bar\epsilon_{jkl} (\frac{1}{\mu\epsilon} D_k B_l)_{,j}$
${T^{iv}}_{,v} = \bar\epsilon_{ikl} (D_k B_l)_{,0}
+ (\frac{1}{2} \eta^{ij} (\frac{1}{\epsilon} D^2 + \frac{1}{\mu} B^2) - \frac{1}{\epsilon} D^i D^j - \frac{1}{\mu} B^i B^j)_{,j}$
$ = \bar\epsilon_{ikl} (D_k B_l)_{,0}
+ \frac{1}{2} (\frac{1}{\epsilon})_{,i} D^2
- (\frac{1}{\epsilon})_{,j} D_i D_j
+ \frac{1}{\epsilon} D_j D_{j,i}
- \frac{1}{\epsilon} D_j D_{i,j}
- \frac{1}{\epsilon} D_i D_{j,j}
+ \frac{1}{2} (\frac{1}{\mu})_{,i} B^2
- (\frac{1}{\mu})_{,j} B_i B_j
+ \frac{1}{\mu} B_j B_{j,i}
- \frac{1}{\mu} B_j B_{i,j}
- \frac{1}{\mu} B^i B_{j,j}
$
electromagnetic stress-energy in medium in general relativity
$F_{uv} = \pmatrix{
0 & -\frac{1}{c} \alpha \mathcal{E}_j \\
\frac{1}{c} \alpha \mathcal{E}_i & \epsilon_{ijk} \mathcal{B}^k
}$
$D_{uv} = \pmatrix{
0 & -c \alpha \mathcal{D}_j \\
c \alpha \mathcal{D}_i & \epsilon_{ijk} \mathcal{H}^k
}$
$(F \cdot D)_{uv} = \pmatrix{
0 & -\frac{1}{c} \alpha \mathcal{E}_k \\
\frac{1}{c} \alpha \mathcal{E}_i & \epsilon_{ikl} \mathcal{B}^l
}
\pmatrix{
-1/\alpha^2 & \beta^m / \alpha^2 \\
\beta^k / \alpha^2 & \gamma^{km} - \beta^k \beta^m / \alpha^2
}
\pmatrix{
0 & -c \alpha \mathcal{D}_j \\
c \alpha \mathcal{D}_m & \epsilon_{mjn} \mathcal{H}^n
}$
$= \pmatrix{
-\frac{1}{\alpha} \frac{1}{c} \mathcal{E}_k \beta^k
&
- \alpha \frac{1}{c} \mathcal{E}^m
+ \frac{1}{\alpha} \frac{1}{c} \mathcal{E}_k \beta^k \beta^m
\\
-\frac{1}{\alpha} \frac{1}{c} \mathcal{E}_i
+ \frac{1}{\alpha}^2 \epsilon_{ikl} \mathcal{B}^l \beta^k
&
\frac{1}{\alpha} \frac{1}{c} \mathcal{E}_i \beta^m
+ {\epsilon_i}^{ml} \mathcal{B}_l
- \frac{1}{\alpha^2} \epsilon_{ikl} \mathcal{B}^l \beta^k \beta^m
}
\pmatrix{
0 & -c \alpha \mathcal{D}_j \\
c \alpha \mathcal{D}_m & \epsilon_{mjn} \mathcal{H}^n
}$
$= \pmatrix{
- \alpha^2 \mathcal{E}^k \mathcal{D}_k
+ \mathcal{E}_k \beta^k \mathcal{D}_l \beta^l
&
\mathcal{D}_j \mathcal{E}_k \beta^k
+ \frac{1}{c} \alpha \epsilon_{jkl} \mathcal{E}^k \mathcal{H}^l
+ \frac{1}{c} \frac{1}{\alpha} \epsilon_{jkl} \mathcal{H}^k \beta^l \mathcal{E}_m \beta^m
\\
\mathcal{E}_i \mathcal{D}_k \beta^k
+ c \alpha \epsilon_{ikl} \mathcal{D}^k \mathcal{B}^l
+ c \frac{1}{\alpha} \epsilon_{ikl} \mathcal{B}^k \beta^l \mathcal{D}_m \beta^m
&
\mathcal{E}_i \mathcal{D}_j
+ \mathcal{H}_i \mathcal{B}_j
- \gamma_{ij} \mathcal{B}_k \mathcal{H}^k
+ c \frac{1}{\alpha} \mathcal{D}_j \epsilon_{ikl} \mathcal{B}^k \beta^l
+ \frac{1}{c} \frac{1}{\alpha} \mathcal{E}_i \epsilon_{jkl} \mathcal{H}^k \beta^l
+ \frac{1}{\alpha^2} \epsilon_{ikl} \mathcal{B}^k \beta^l \epsilon_{jmn} \mathcal{H}^m \beta^n
}$
$tr(F \cdot D) = (F \cdot D)_{uv} g^{uv}$
$=
-1/\alpha^2 (
- \alpha^2 \mathcal{E}^k \mathcal{D}_k
+ \mathcal{E}_k \beta^k \mathcal{D}_l \beta^l
)
+ \frac{1}{\alpha^2} \beta^j (
\mathcal{D}_j \mathcal{E}_k \beta^k
+ \frac{1}{c} \alpha \epsilon_{jkl} \mathcal{E}^k \mathcal{H}^l
+ \frac{1}{c} \frac{1}{\alpha} \epsilon_{jkl} \mathcal{H}^k \beta^l \mathcal{E}_m \beta^m
)
+ \frac{1}{\alpha^2} \beta^i (
\mathcal{E}_i \mathcal{D}_k \beta^k
+ c \alpha \epsilon_{ikl} \mathcal{D}^k \mathcal{B}^l
+ c \frac{1}{\alpha} \epsilon_{ikl} \mathcal{B}^k \beta^l \mathcal{D}_m \beta^m
)
+ (\gamma^{ij} - \frac{1}{\alpha^2} \beta^i \beta^j) (
\mathcal{E}_i \mathcal{D}_j
+ \mathcal{H}_i \mathcal{B}_j
- \gamma_{ij} \mathcal{B}_k \mathcal{H}^k
+ c \frac{1}{\alpha} \mathcal{D}_j \epsilon_{ikl} \mathcal{B}^k \beta^l
+ \frac{1}{c} \frac{1}{\alpha} \mathcal{E}_i \epsilon_{jkl} \mathcal{H}^k \beta^l
+ \frac{1}{\alpha^2} \epsilon_{ikl} \mathcal{B}^k \beta^l \epsilon_{jmn} \mathcal{H}^m \beta^n
)
$
$=
2 \mathcal{E}^k \mathcal{D}_k
- 2 \mathcal{B}_k \mathcal{H}^k
+ 2 \frac{1}{\alpha^2} \beta^i \beta_i \mathcal{B}_k \mathcal{H}^k
- 2 \frac{1}{\alpha^2} \mathcal{H}_i \beta^i \mathcal{B}_j \beta^j
+ 2 \frac{1}{c} \frac{1}{\alpha} \epsilon_{jkl} \mathcal{E}^j \mathcal{H}^k \beta^l
+ 2 c \frac{1}{\alpha} \epsilon_{ikl} \mathcal{D}^i \mathcal{B}^k \beta^l
$
For $\mathcal{D} = \epsilon \mathcal{E}$ and $\mathcal{H} = \frac{1}{\mu} \mathcal{B}$:
$tr(F \cdot D) =
2 \epsilon \mathcal{E}^2
- 2 \frac{1}{\mu} \mathcal{B}^2
+ 2 \frac{1}{\mu} \frac{1}{\alpha^2} \beta^2 \mathcal{B}^2
- 2 \frac{1}{\mu} \frac{1}{\alpha^2} (\mathcal{B}_k \beta^k)^2
+ 2 \frac{1}{\mu} \frac{1}{c} \frac{1}{\alpha} \epsilon_{jkl} \mathcal{E}^j \mathcal{B}^k \beta^l
+ 2 \epsilon c \frac{1}{\alpha} \epsilon_{ikl} \mathcal{E}^i \mathcal{B}^k \beta^l
$
Continuing with the scalar assumption of $\epsilon$ and $\mu$...
$T_{uv} = -(F \cdot D)^{TF}_{uv} = -(F \cdot D)_{uv} - \frac{1}{4} g_{uv} tr (F \cdot D)$
$= \pmatrix{
\epsilon \alpha^2 \mathcal{E}^2
- \epsilon (\mathcal{E}_k \beta^k)^2
- \frac{1}{2} \epsilon \alpha^2 \mathcal{E}^2
+ \frac{1}{2} \frac{1}{\mu} \alpha^2 \mathcal{B}^2
- \frac{1}{2} \frac{1}{\mu} \beta^2 \mathcal{B}^2
+ \frac{1}{2} \frac{1}{\mu} (\mathcal{B}_k \beta^k)^2
- \frac{1}{\mu} \frac{1}{c} \alpha \frac{1}{2} \epsilon_{klm} \mathcal{E}^k \mathcal{B}^l \beta^m
- \epsilon c \alpha \frac{1}{2} \epsilon_{klm} \mathcal{E}^k \mathcal{B}^l \beta^m
+ \frac{1}{2} \epsilon \beta^2 \mathcal{E}^2
- \frac{1}{2} \frac{1}{\mu} \beta^2 \mathcal{B}^2
+ \frac{1}{2} \frac{1}{\mu} \frac{1}{\alpha^2} \beta^4 \mathcal{B}^2
- \frac{1}{2} \frac{1}{\mu} \frac{1}{\alpha^2} \beta^2 (\mathcal{B}_k \beta^k)^2
+ \frac{1}{2} \frac{1}{\mu} \frac{1}{c} \frac{1}{\alpha} \beta^2 \epsilon_{klm} \mathcal{E}^k \mathcal{B}^l \beta^m
+ \frac{1}{2} \epsilon c \frac{1}{\alpha} \beta^2 \epsilon_{klm} \mathcal{E}^k \mathcal{B}^l \beta^m
&
- \epsilon \mathcal{E}_j \mathcal{E}_k \beta^k
- \frac{1}{\mu} \frac{1}{c} \alpha \epsilon_{jkl} \mathcal{E}^k \mathcal{B}^l
- \frac{1}{\mu} \frac{1}{c} \frac{1}{\alpha} \epsilon_{jkl} \mathcal{B}^k \beta^l \mathcal{E}_m \beta^m
+ \frac{1}{2} \epsilon \mathcal{E}^2 \beta_j
- \frac{1}{2} \frac{1}{\mu} \mathcal{B}^2 \beta_j
+ \frac{1}{2} \frac{1}{\mu} \frac{1}{\alpha^2} \beta^2 \mathcal{B}^2 \beta_j
- \frac{1}{2} \frac{1}{\mu} \frac{1}{\alpha^2} (\mathcal{B}_k \beta^k)^2 \beta_j
+ \frac{1}{2} \frac{1}{\mu} \frac{1}{c} \frac{1}{\alpha} \epsilon_{klm} \mathcal{E}^k \mathcal{B}^l \beta^m \beta_j
+ \frac{1}{2} \epsilon c \frac{1}{\alpha} \epsilon_{klm} \mathcal{E}^k \mathcal{B}^l \beta^m \beta_j
\\
- \epsilon \mathcal{E}_i \mathcal{E}_k \beta^k
- \epsilon c \alpha \epsilon_{ikl} \mathcal{E}^k \mathcal{B}^l
- \epsilon c \frac{1}{\alpha} \epsilon_{ikl} \mathcal{B}^k \beta^l \mathcal{E}_m \beta^m
+ \frac{1}{2} \epsilon \mathcal{E}^2 \beta_i
- \frac{1}{2} \frac{1}{\mu} \mathcal{B}^2 \beta_i
+ \frac{1}{2} \frac{1}{\mu} \frac{1}{\alpha^2} \beta^2 \mathcal{B}^2 \beta_i
- \frac{1}{2} \frac{1}{\mu} \frac{1}{\alpha^2} (\mathcal{B}_k \beta^k)^2 \beta_i
+ \frac{1}{2} \frac{1}{\mu} \frac{1}{c} \frac{1}{\alpha} \epsilon_{klm} \mathcal{E}^k \mathcal{B}^l \beta^m \beta_i
+ \frac{1}{2} \epsilon c \frac{1}{\alpha} \epsilon_{klm} \mathcal{E}^k \mathcal{B}^l \beta^m \beta_i
&
- \epsilon \mathcal{E}_i \mathcal{E}_j
- \frac{1}{\mu} \mathcal{B}_i \mathcal{B}_j
+ \frac{1}{\mu} \gamma_{ij} \mathcal{B}_k \mathcal{B}^k
- \epsilon c \frac{1}{\alpha} \mathcal{E}_j \epsilon_{ikl} \mathcal{B}^k \beta^l
- \frac{1}{\mu} \frac{1}{c} \frac{1}{\alpha} \mathcal{E}_i \epsilon_{jkl} \mathcal{B}^k \beta^l
- \frac{1}{\mu} \frac{1}{\alpha^2} \epsilon_{ikl} \mathcal{B}^k \beta^l \epsilon_{jmn} \mathcal{B}^m \beta^n
+ \frac{1}{2} \epsilon \mathcal{E}^2 \gamma_{ij}
- \frac{1}{2} \frac{1}{\mu} \mathcal{B}^2 \gamma_{ij}
+ \frac{1}{2} \frac{1}{\mu} \frac{1}{\alpha^2} \beta^2 \mathcal{B}^2 \gamma_{ij}
- \frac{1}{2} \frac{1}{\mu} \frac{1}{\alpha^2} (\mathcal{B}_k \beta^k)^2 \gamma_{ij}
+ \frac{1}{2} \frac{1}{\mu} \frac{1}{c} \frac{1}{\alpha} \epsilon_{klm} \mathcal{E}^k \mathcal{B}^l \beta^m \gamma_{ij}
+ \frac{1}{2} \epsilon c \frac{1}{\alpha} \epsilon_{klm} \mathcal{E}^k \mathcal{B}^l \beta^m \gamma_{ij}
}$
TODO do this, the matter SR, and the vaccuum GR / SR all in the same co/contra-variance, which ever one you choose.
Antisymmetric component of general-relativistic stress-energy:
$T_{[uv]} = \pmatrix{
0 &
\frac{1}{2} \frac{1}{c} \frac{1}{\mu} (n^2 - 1) (
\alpha \epsilon_{jkl} \mathcal{E}^k \mathcal{B}^l
+ \frac{1}{\alpha} \epsilon_{jkl} \mathcal{B}^k \beta^l \mathcal{E}_m \beta^m
)
\\
-\frac{1}{2} \frac{1}{c} \frac{1}{\mu} (n^2 - 1) (
\alpha \epsilon_{ikl} \mathcal{E}^k \mathcal{B}^l
+ \frac{1}{\alpha} \epsilon_{ikl} \mathcal{B}^k \beta^l \mathcal{E}_m \beta^m
)
&
\frac{1}{c} \frac{1}{\mu} (n^2 - 1) \frac{1}{\alpha} \mathcal{E}_{[i} \epsilon_{j]kl} \mathcal{B}^k \beta^l
}$
Let $\mathcal{S}_i = {\epsilon_i}^{jk} \mathcal{E}_j \mathcal{H}_k$
$T_{[uv]} = \pmatrix{
0 &
\frac{1}{2} \frac{1}{c} (n^2 - 1) (
\alpha \mathcal{S}_j
+ \frac{1}{\alpha} \epsilon_{jkl} \mathcal{H}^k \beta^l \mathcal{E}_m \beta^m
)
\\
-\frac{1}{2} \frac{1}{c} (n^2 - 1) (
\alpha \mathcal{S}_i
+ \frac{1}{\alpha} \epsilon_{ikl} \mathcal{H}^k \beta^l \mathcal{E}_m \beta^m
)
&
\frac{1}{c} (n^2 - 1) \frac{1}{\alpha} \mathcal{E}_{[i} \epsilon_{j]kl} \mathcal{H}^k \beta^l
}$
So here is the antisymmetric component of our Ricci tensor.
Einstein Field equations, asymmetric terms:
$G_{[ab]} = 8 \pi \frac{G}{c^4} T_{[ab]}$
Notice since $G_{ab} = R_{ab} - \frac{n}{2} R g_{ab}$, so $G_{[ab]} = R_{[ab]}$ even if $T_{ab}$ wasn't trace-free (since $T_{[ab]}$ is always trace-free).
It just so happens that the electromagnetism stress-energy is trace-free anyways, but even if it wasn't this would hold true:
$R_{[ab]} = 8 \pi \frac{G}{c^4} T_{[ab]}$
Time to kick it up a notch from Einstein to Einstein-Cartan.
See my worksheet here for where I got this from.
Riemann curvature with torsion (but no commutation):
${R^d}_{cab} =
\partial_a ({\Gamma^d}_{bc})
- \partial_b ({\Gamma^d}_{ac})
+ {\Gamma^d}_{au} {\Gamma^u}_{bc}
- {\Gamma^d}_{bu} {\Gamma^u}_{ac}
$
Ricci tensor with torsion:
$R_{ab} = \partial_u ({\Gamma^u}_{ba})
- \partial_b ({\Gamma^u}_{ua})
+ {\Gamma^u}_{uv} {\Gamma^v}_{ba}
- {\Gamma^u}_{bv} {\Gamma^v}_{ua}
$
Antisymmetry of Ricci tensor with torsion:
$R_{[ab]} = \frac{3}{2} (
\nabla_{[u} {T^u}_{ba]} + {T^u}_{v[u} {T^v}_{ba]}
)$
Assuming zero shift ($\beta^i = 0$)...
I'm keeping $\alpha$ still, because most of the Earth's gravity comes from $\alpha$ and not $\beta$, but who knows, maybe $\beta$ becomes important in asymmetric stress-energies.
$T_{[uv]} = \pmatrix{
0 &
\frac{1}{2} \frac{1}{c} (n^2 - 1) \alpha \mathcal{S}_j \\
-\frac{1}{2} \frac{1}{c} (n^2 - 1) \alpha \mathcal{S}_i &
0
}$
$T_{[0i]} = \frac{1}{2} \frac{1}{c} (n^2 - 1) \alpha \mathcal{S}_i$
$R_{[0i]} = \frac{3}{2} ( \nabla_{[u} {T^u}_{i0]} + {T^u}_{v[u} {T^v}_{i0]})$
$S_i$ has units $\frac{kg}{s^3}$
$T_{[0i]}$ has units $\frac{kg}{m \cdot s^2}$
If $R_{ab}$ has units $\frac{1}{m^2}$ then that means $\nabla_a {T^b}_{cd}$ has units $\frac{1}{m^2}$ and ${T^a}_{bc}$ has units $\frac{1}{m}$
$4 \pi \frac{G}{c^5} (n^2 - 1) \alpha \mathcal{S}_i = \frac{3}{2} (\nabla_{[u} {T^u}_{i0]} + {T^u}_{v[u} {T^v}_{i0]}) $
$8 \pi \frac{G}{c^5} (n^2 - 1) \alpha \mathcal{S}_i =
\nabla_u {T^u}_{i0}
+ \nabla_i {T^u}_{0u}
+ \nabla_0 {T^u}_{ui}
+ {T^u}_{vu} {T^v}_{i0}
+ {T^u}_{vi} {T^v}_{0u}
+ {T^u}_{v0} {T^v}_{ui}
$
$8 \pi \frac{G}{c^5} (n^2 - 1) \alpha \mathcal{S}_i =
\nabla_j {T^j}_{i0}
+ \nabla_i {T^j}_{0j}
+ \nabla_0 {T^j}_{ji}
+ {T^0}_{i0} {T^j}_{0j}
+ {T^0}_{k0} {T^k}_{i0}
+ {T^j}_{kj} {T^k}_{i0}
$
... assuming steady-state, and that $({T^a}_{bc})^2 \approx 0$...
$8 \pi \frac{G}{c^5} (n^2 - 1) \alpha \mathcal{S}_i = \nabla_j {T^j}_{i0} - \nabla_i {T^j}_{j0}$
$8 \pi \frac{G}{c^5} (n^2 - 1) \alpha \mathcal{S}_i = \epsilon_{ijp} \epsilon^{mnp} \nabla_m {T^j}_{0n}$
... so the Poynting vector is equal to the curl of the torsion two-form component perpendicular to time (second index spatial, second index time) ... crossed with the torsion first index ...
...assume zero trace of the torsion...
$8 \pi \frac{G}{c^5} (n^2 - 1) \alpha \mathcal{S}_i = \nabla_j {T^j}_{i0}$
... one term falls off, and now the Poynting vector equals the divergence of the first index of the torsion...
The other term doesn't show up explicitly in the torsion influence of the geodesic anyways...
Now to assume zero-curl of the torsion first index and solve with Helmholtz decomposition:
${T^k}_{i0}(\vec{x}) = -\nabla^k \frac{1}{4 \pi} \int \frac{\nabla_j {T^j}_{i0}(\vec{y})}{|\vec{x} - \vec{y}|^2} d\vec{y}$
${T^k}_{i0}(\vec{x}) = -2 \frac{G}{c^5} (n^2 - 1) \nabla^k \int \frac{\alpha \mathcal{S}_i(\vec{x})}{|\vec{x} - \vec{y}|^2} d\vec{y}$
Now compare this to the geodesic equation influence from torsion.
See my worksheet here for where I got this from.
$\dot{u}^a = -{\Gamma^a}_{bc} u^b u^c + u_b {T^b}_{dc} u^c g^{da}$
...just looking at the torsion influence of gravitation on the spatial components...
$(\dot{u}^i)_{\Theta} = u_b {T^b}_{dc} u^c g^{di}$
...maintain that no-shift contraint: $\beta^i = 0$, and $u_0 \approx 1$...
$(\dot{u}^i)_{\Theta} = u_b {T^b}_{jc} u^c \gamma^{ji}$
$(\dot{u}^i)_{\Theta} =
-{T^0}_{j0} \gamma^{ji}
+ {T^0}_{jl} u^l \gamma^{ji}
- u_k {T^k}_{j0} \gamma^{ji}
+ u_k {T^k}_{jl} u^l \gamma^{ji}
$
... isolate the component related to the Poynting vector (with steady-state and zero trace torsion)...
$(\dot{u}^i)_{\hat\Theta} = -u_k {T^{ki}}_0$
...substitute ...
$(\dot{u}^i)_{\hat\Theta} = 2 \frac{G}{c^5} (n^2 - 1) u^k \nabla_k \int \frac{\alpha \mathcal{S}^i(\vec{x})}{|\vec{x} - \vec{y}|^2} d\vec{y} $
And there you have it:
The influence of gravitation (geodesic ficudial force) due to torsion, specifically from the Poynting field.
Example:
$u^a$ is unitless, $u^i = \frac{v^i}{c}$
$S$ is in units of $\frac{W}{m^2} = \frac{kg}{s^3}$
$\int S d\vec{x}$ is in units of $W m = \frac{kg \cdot m^3}{s^3}$
$u^k \nabla_k \int S d\vec{x}$ is in units of $W = \frac{kg \cdot m^2}{s^3}$
$\frac{G}{c^5} u^k \nabla_k \int S d\vec{x}$ is unitless.
So within the context of units:
$(\dot{v}^i)_{\hat\Theta} = 2 \frac{G}{c^5} (n^2 - 1) \gamma^{ji} v^k \nabla_k \int \frac{\alpha \mathcal{S}_j(\vec{x})}{|\vec{x} - \vec{y}|^2} d\vec{y} $
...for velocity $v^i = c u^i$
gravitational constant: $G = 6.67384 \cdot 10^{-11} \cdot \frac{m^3}{kg \cdot s^2}$
speed of light: $c = 299792458 \cdot \frac{m}{s}$
$|\dot{v}|_{\hat\Theta} = 2 \cdot \frac{6.67384 \cdot 10^{-11}}{(299792458)^5} (n^2 - 1) \cdot |v| $ ... times Poynting magnitude?
$= 5.5119119535891 \cdot 10^{-53} \cdot (n^2 - 1) \cdot |v| $ ... times Poynting magnitude?
All in all, this is a really longwinded way to come up with the relation between the geodesic acceleration, the torsion, and (specifically) the Poynting vector (times the spatial velocity, again).
Why bother, especially when the geodesic acceleration itself is directly related to the torsion components ${T^0}_{j0}$ without any scalar of the velocity.
${T^0}_{j0}$ is related to $S_j$ by $S_i \propto \nabla_0 {T^0}_{i0}$
So how does $S_j$ influence ${T^0}_{i0}$?
This is still scaled by $\frac{G}{c^5}$...
Tensor susceptibility:
(from https://en.wikipedia.org/wiki/Magnetic_susceptibility)
$P^i = \epsilon_0 {(\chi_e)^i}_j E^j$
$M_i = {(\chi_m)_i}^j H_j$ ... non-scalar in crystals, nonlinear in ferromagnetic crystals, complex in terms of AC fields
$D^i = \epsilon_0 E^i + P^i = \epsilon_0 (\delta^i_j + {(\chi_e)^i}_j) E^j$
$\mu_0 (\delta^i_j + {(\chi_m)^i}_j ) H_i = B_j$
${\mu^i}_j H_i = B_j$
Let $D^u = {\epsilon^u}_v E^v$
Let $H^u = {(\mu^{-1})^u}_v B^v$
$F_{uv} = \downarrow u(i) \overset{\rightarrow v(i)}{\pmatrix{
0 & -\frac{1}{c} E_j \\
\frac{1}{c} E_i & {\bar\epsilon_{ij}}^k B_k
}}$
$D_{uv} = \downarrow u(i) \overset{\rightarrow v(i)}{\pmatrix{
0 & -c D_j \\
c D_i & {\bar\epsilon_{ij}}^k H_k
}}$
In special relativity:
${\epsilon_a}^u F_{uv} {\epsilon_b}^v$
$= \pmatrix{
c^2 & 0 \\
0 & {\epsilon_i}^m
} \pmatrix{
0 & -\frac{1}{c} E_n \\
\frac{1}{c} E_m & {\bar\epsilon_{mn}}^k B_k
} \pmatrix{
c^2 & 0 \\
0 & {\epsilon_j}^n
}$
$ = \pmatrix{
0 & -c {\epsilon_j}^n E_n \\
c {\epsilon_i}^m E_m & {\epsilon_i}^m {\epsilon_j}^n {\bar\epsilon_{mn}}^k B_k
}$
$D_{uv} = \pmatrix{
0 & -c {\epsilon_j}^k E_k \\
c {\epsilon_i}^k E_k & {\bar\epsilon_{ij}}^k {(\mu^{-1})_k}^m B_m
}$
For ${\epsilon^i}_m {\epsilon^j}_n \bar\epsilon^{mnk} = \bar\epsilon^{ijm} {(\mu^{-1})^k}_m$
this would simplify further to $D_{ab} = {\epsilon_a}^u F_{uv} {\epsilon_b}^v$.
But that isn't often the case, because it is neglecting the use of the $M_{ab}$ tensor.
Maxwell's equations in a medium, in vector 4-form:
(from https://en.wikipedia.org/wiki/Covariant_formulation_of_classical_electromagnetism)
When $\chi$'s are just scalars, we get this:
$D^{uv} u_v = c^2 \epsilon F^{uv} u_v$
$\star D^{uv} u_v = \frac{1}{\mu} \star F^{uv} u_v$
expanded, in special relativity:
$D^{uv} u_v = \pmatrix{
0 & c \epsilon E_j \\
-c \epsilon E_i & {\bar\epsilon_{ij}}^k \mu^{-1} B_k
} \cdot \pmatrix{
u^t \\
u^j
} = \pmatrix{
c \epsilon E_j u^j \\
-c \epsilon E_i u^t + \mu^{-1} {\bar\epsilon_{ij}}^k B_k u^j
}$
$c^2 \epsilon F^{uv} u_v = c^2 \epsilon \pmatrix{
0 & \frac{1}{c} E_j \\
-\frac{1}{c} E_i & {\bar\epsilon_{ij}}^k B_k
} \cdot \pmatrix{
u^t \\
u^j
} = \pmatrix{
c \epsilon E_j u^j \\
-c \epsilon E_i u^t + c^2 \epsilon {\bar\epsilon_{ij}}^k B_k u^j
}$
...but this is only true when $\frac{1}{c^2} = \mu\epsilon$, which isn't always necessarily true, right?
combining non-dispersive materials with general relativity formulation of electromagnetism:
Notice that $E^u$'s and $B^u$'s are suppose to be spatial.
I'm suspicious that the $\epsilon_{ab}$'s and $(\mu^{-1})_{ab}$'s are also spatial as well
$F^{ab} = \frac{1}{c} n^a E^b - \frac{1}{c} n^b E^a + n_d \epsilon^{dabc} B_c$
$D^{ab} = c n^a D^b - c n^b D^a + n_d \epsilon^{dabc} H_c$
If susceptibility is a scalar:
$D^{ab} = c \epsilon (n^a E^b - n^b E^a) + n_d \epsilon^{dabc} \mu^{-1} B_c$
$\mu D^{ab} = c^2 \mu \epsilon \cdot \frac{1}{c} (n^a E^b - n^b E^a) + n_d \epsilon^{dabc} B_c$
$\mu D^{ab} = (\frac{c}{v_p})^2 \frac{1}{c} (n^a E^b - n^b E^a) + n_d \epsilon^{dabc} B_c$
$\mu D^{ab} = ((\frac{c}{v_p})^2 - 1) \frac{1}{c} (n^a E^b - n^b E^a) + F^{ab}$
or...
$\mu_0 D^{ab} = \frac{\epsilon}{\epsilon_0} \frac{1}{c} (n^a E^b - n^b E^a) + n_d \epsilon^{dabc} \frac{\mu_0}{\mu} B_c$
$\mu_0 D^{ab} = (1 + \chi_e) \frac{1}{c} (n^a E^b - n^b E^a) + (\frac{1}{1 + \chi_m}) n_d \epsilon^{dabc} B_c$
$D^{ab} = \frac{1}{\mu_0} \chi_e \frac{1}{c} (n^a E^b - n^b E^a) + \frac{1}{\mu_0} (\frac{1}{1 + \chi_m} - 1) n_d \epsilon^{dabc} B_c + \frac{1}{\mu_0} F^{ab}$
$D^{ab} = \frac{1}{\mu_0} F^{ab} - M^{ab}$
Therefore we get
$M^{ab} = -\frac{1}{\mu_0} (\chi_e \frac{1}{c} (n^a E^b - n^b E^a) + (\frac{1}{1 + \chi_m} - 1) n_d \epsilon^{dabc} B_c)$
$M^{ab} = -( c \chi_e \epsilon_0 (n^a E^b - n^b E^a) + \frac{1}{\mu_0} (\frac{\chi_m}{1 + \chi_m}) n_d \epsilon^{dabc} B_c)$
$M^{ab} = -( c (n^a P^b - n^b P^a) + n_d \epsilon^{dabc} M_c)$
...for $\epsilon_0 \chi_e E^a = P^a, \chi_m H_a = M_a$
If susceptibility is a tensor:
$D^{ab} = c (n^a {\epsilon^b}_u - n^b {\epsilon^a}_u) E^u + n_d \epsilon^{dabc} {(\mu^{-1})_c}^u B_u$
${F^a}_c {D^c}_b =
(\frac{1}{c} (n^a E_c - n_c E^a) + n_d {\epsilon^{da}}_{ce} B^e)
\cdot (c (n^c D_b - n_b D^c) + n_f {\epsilon^{fc}}_{bg} H^g)
$
${F^a}_c {D^c}_b =
(\frac{1}{c} (n^a E_c - n_c E^a) + n_d {\epsilon^{da}}_{ce} B^e)
\cdot
(c (n^c \epsilon_{bu} - n_b {\epsilon^c}_u) E^u + n_f {\epsilon^{fc}}_{bh} (\mu^{-1})^{hu} B_u)
$
${F^a}_c {D^c}_b =
- n^a E_c n_b {\epsilon^c}_u E^u
+ \frac{1}{c} n^a E_c n_f {\epsilon^{fc}}_{bh} (\mu^{-1})^{hu} B_u
- n_c E^a n^c \epsilon_{bu} E^u $ $
+ n_c E^a n_b {\epsilon^c}_u E^u
- \frac{1}{c} n_c E^a n_f {\epsilon^{fc}}_{bh} (\mu^{-1})^{hu} B_u
+ c n_d {\epsilon^{da}}_{ce} B^e n^c \epsilon_{bu} E^u $ $
- c n_d {\epsilon^{da}}_{ce} B^e n_b {\epsilon^c}_u E^u
+ n_d {\epsilon^{da}}_{ce} B^e n_f {\epsilon^{fc}}_{bh} (\mu^{-1})^{hu} B_u
$
TODO use projection tensors to apply $\epsilon$ and $\mu$ to only the spatial components of $E$ and $B$.
In fact,
If $F = E \wedge dt + \star (B \wedge dt)$
so $F_{[ab]} = 2 E_{[a} n_{b]} + {\epsilon_{ab}}^{cd} n_{[c} B_{d]}$
Then I'm suspicious that tensors $\star \epsilon = (\mu)^{-1}$ or something like that.
Maxwell's laws in general relativity using the constitutive tensor as a combination of permeabillity and permittivity tensors:
$D_{uv} = {Z_{uv}}^{ab} F_{ab}$
$D_{0i} = c^2 \epsilon F_{0i}$
$D_{i0} = c^2 \epsilon F_{i0}$
$D_{ij} = \frac{1}{\mu} F_{ij}$
Therefore
$c^2 \epsilon = \frac{\epsilon}{\epsilon_0} \frac{1}{\mu_0} = {Z_{0i}}^{0i} = {Z_{i0}}^{i0} = -{Z_{0i}}^{i0} = -{Z_{i0}}^{0i}$
$\frac{1}{\mu} = \frac{\mu_0}{\mu} \frac{1}{\mu_0} = {Z_{jk}}^{jk} = {Z_{kj}}^{kj} = -{Z_{jk}}^{kj} = -{Z_{kj}}^{jk}$
${Z_{aa}}^{bc} = {Z_{bc}}^{aa} = 0$
How suspicious, when the susceptibilities are scalars the tensor that relates the Faraday tensor with the displacement tensor fulfills all the properties of a Riemann curvature tensor.
${Z_{ab}}^{cd} = \frac{1}{\mu_0} \left[\begin{matrix}
& & 01 & 02 & 03 & 23 & 31 & 12 \\
\\
01 & & \epsilon_r & \cdot & \cdot & \cdot & \cdot & \cdot \\
02 & & \cdot & \epsilon_r & \cdot & \cdot & \cdot & \cdot \\
03 & & \cdot & \cdot & \epsilon_r & \cdot & \cdot & \cdot \\
23 & & \cdot & \cdot & \cdot & (\mu_r)^{-1} & \cdot & \cdot \\
31 & & \cdot & \cdot & \cdot & \cdot & (\mu_r)^{-1} & \cdot \\
12 & & \cdot & \cdot & \cdot & \cdot & \cdot & (\mu_r)^{-1}
\end{matrix}\right]$
for $\epsilon_r = \frac{\epsilon}{\epsilon_0}, \mu_r = \frac{\mu}{\mu_0}$
It also fulfills all the properties of a Levi-Civita permutation tensor, and of course when your coordinate system is rotation angles then the Riemann curvature tensor equals the Levi-Civita permutation tensor.
Maybe that's because the basis of the electromagnetic fields are rotations around the basis of the position vector?
Another thing to think about: in Kaluza-Klein the Gauss-Ampere constraint ${D_{au}}^{;u} = ({Z_{au}}^{bc} F_{bc})^{;u} \approx {Z_{au}}^{bc} {F_{bc}}^{;u} = 4 \pi J_a$ shows up as a 5th dimension right beside the Riemann curvature constraint in the Einstein Field Equations: $G_{ab} = {(\star R \star)^u}_{aub} = 8 \pi T_{ab}$.
It looks even more uncannily like the Riemann curvature tensor if we treat the approximation to be an equality: ${Z_{ab}}^{cd;u} = 0$, just like the Bianchi identity ${R_{ab}}^{[cd;e]} = 0$, though the Gauss-Ampere law itself appears in vacuum form in Kaluza-Klein with cylindrical constrints. Maybe we can relax that constraint to add in some extra terms?
Notice that potentials $\phi = A_0 = g_{50}$ is on the same footing as $\frac{1}{2} \phi_{grav} \approx 1 - g_{ii} \approx -1 - g_{00}$.
Notice as well the Lorentz force contributions $F = dA$ is on the same footing as the gravitational force contribution from the spacetime connection 1-form $de_b = {\omega^a}_b \otimes e_a = {\Gamma^a}_{cb} dx^c \otimes e_a$.