I have written code that does carries out the following task:
Construct a fully symmetric rank-$P$ "tensor" out of a collection of user-input rank-1 vectors, and a rank-2 symmetric tensor $h$. The number of times each vector is to be repeated, and the indices of the final rank-$P$ tensor are user-input.
Example #1: Three vectors, three indices-- Given input vectors $u$ repeated once, $v$ repeated twice, and indices {$a,b,c$}, the output should be the symmetric tensor:
$$u_a v_b v_c + u_b v_a v_c + u_c v_a v_b$$
Example #2: Two vectors, four indices-- This time, take input vector $u$ repeated twice, but with four indices {$a,b,c,d$}. Since there are two excess indices, they must be assigned to $h$. The output:
$$u_a u_b h_{c d}+u_a u_c h_{b d}+u_a u_d h_{b c}+u_b u_c h_{a d}+u_b u_d h_{a c}+u_c u_d h_{a b}$$
Note: User must input such that number of excess indices must be even, and can never be less than the number of rank-1 tensors.
My code here as follows, starting with definitions and input:
(*Tells engine that \[DoubleStruckH] is symmetric *)
Clear[\[DoubleStruckH]];
Subscript[\[DoubleStruckH], i_, j_] ^:= Subscript[\[DoubleStruckH], j, i] /; ! OrderedQ[{i, j}];
Here are the input variables:
(*Example 1*)
sym = {u, v}; (*symbols of rank-1 vectors*)
n = {1, 2}; (*repetition: u once, v twice*)
indices = {a, b, c}; (*Indices of rank-3 tensor.*)
(*Example 2*)
sym = {u}; (*symbols of rank-1 vectors*)
n = {2}; (*repetition: u twice*)
indices = {a, b, c, d}; (*Indices of rank-4 tensor.*)
Useful definitions:
P = Length[indices]; (*rank of final tensor*)
r = (P - Total[n])/2; (*number of \[DoubleStruckH] in tensor*)
So now my clumsy code: strategy is to begin with one-representative term of the symmetric tensor (which needs some help).
representative = Apply[Function[#,
Product[1/
n[[j]]! Product[Subscript[
sym[[j]], #[[i]]], {i, 1 + Sum[n[[k]], {k, 1, j - 1}],
Sum[n[[k]], {k, 1, j}]}], {j, 1, Length[sym]}]*1/(r! 2^r)
Product[
Subscript[\[DoubleStruckH], #[[i]], #[[i + 1]]], {i,
Total[n] + 1, P, 2}]] &, {Table[Unique[], {P}]}]
In the case of (*Example 1*)
, the output of representative
is:
I know I need help here, since representative
is supposed to be a pure function, with the arguments being the slots where the indices go. But the output of this reveals that it doesn't fully evaluate the products.
Anyway, carrying on. Now all I need is to add up all possible permutations:
answer = Plus @@ With[{func = representative}, Map[Apply[func, #] &, Permutations[indices]]]
It works, but it slow for large ($P\gtrsim 6$) tensors.
The two places where I need help is:
I need to make this work somewhat quickly for large tensors (with 6 or 8 indices). Clearly the code does more work than it needs to because of the factorials that I use in
representative
to divide out extra terms.My
representative
function needs help, too (see above).