Local Magneto-shear instability in Newtonian gravity

Author: Thomas Celora

Clear["Global`*"]
systemInfo[] := Column[{Row[{Style["Last Modified: ", Bold], DateString[FileDate[NotebookFileName[]]]}]}]
clearMemory[] := {MemoryInUse[], Unprotect[In, Out], Clear[In, Out], Protect[In, Out], MemoryInUse[]};

systemInfo[]
clearMemory[];

Sheared background

Dispersion relation

Calculating the determinant

mat = {{-om, ro k1, ro k2, ro k3,0,0,0},{cs^2 / ro k1 ,-om -I s1,0,0,I1,b2 k1, b3 k1},{cs^2 / ro k2, 0,-om - I s2, 0,b1 k2,I2, b3 k2},{cs^2 / ro k3 ,0,0, -om +I(s1+s2),b1 k3, b2 k3,I3},{0,I1, b1 k2,b1 k3,-om +I s1,0,0},{0,b2 k1,I2,b2 k3,0, -om +I s2,0},{0,b3 k1 , b3 k2,I3,0,0,-om -I(s1 + s2)}};
mat // MatrixForm
shearmat = ReplaceAll[mat, {I1 -> -b2 k2 - b3 k3, I2 -> -b1 k1 - b3 k3, I3-> -b1 k1 - b2 k2}];
sheardet = Det[shearmat] ;

Computing the coefficients and writing them in a separate file

Collect[sheardet,om,Simplify];
a0 = ReplaceAll[Expand[sheardet], {om^7->0,om^6 ->0,om^5->0,om^4->0,om^3->0,om^2->0,om->0}] ;
a1 = ReplaceAll[Expand[sheardet] - a0, {om^7->0,om^6 ->0,om^5->0,om^4->0,om^3->0,om^2->0}] /. om->1 ;
a2 = ReplaceAll[Expand[sheardet]- a0-a1 om, {om^7->0,om^6 ->0,om^5->0,om^4->0,om^3->0}] /. om->1 ;
a3 = ReplaceAll[Expand[sheardet] - a0-a1 om-a2 om^2, {om^7->0,om^6 ->0,om^5->0,om^4->0}] /. om->1 ;
a4 = ReplaceAll[Expand[sheardet] - a0-a1 om-a2 om^2-a3 om^3, {om^7->0,om^6 ->0,om^5->0}] /. om->1 ;
a5 = ReplaceAll[Expand[sheardet] - a0-a1 om-a2 om^2-a3 om^3-a4 om^4, {om^7->0,om^6 ->0}] /. om->1 ;
a6 = ReplaceAll[Expand[sheardet]- a0-a1 om-a2 om^2-a3 om^3-a4 om^4-a5 om^5, {om^7->0}] /. om->1 ;
a7 = (Expand[sheardet] - a0-a1 om-a2 om^2-a3 om^3-a4 om^4-a5 om^5-a6 om^6) /. om ->1;

coefflist = {a0,a1,a2,a3,a4,a5,a6,a7};
For[i=1,i < Length[coefflist]+1,i++,If[coefflist[[i]]==0, Print["0: ",i-1]] ]
For[i=1,i < Length[coefflist]+1,i++,If[coefflist[[i]]==-1, Print["-1: ", i-1]] ]
For[i=1,i < Length[coefflist]+1,i++,If[coefflist[[i]]==1, Print["1: ",i-1]] ]

SetDirectory[NotebookDirectory[]];
DumpSave["shearbg_coeff.mx", {a0,a1,a2,a3,a4,a5}];

Coefficients decomposition

SetDirectory[NotebookDirectory[]];
<< "shearbg_coeff.mx";
B= {b1,b2,b3};
K = {k1,k2,k3};
U = Cross[B,K]; 
shear = {{s1,0,0},{0,s2,0},{0,0, -(s1+s2)}} ;
shear2 = shear . shear ;
shear3 = shear.shear.shear;

W =TensorContract[  Outer[Times, LeviCivitaTensor[3] , shear] ,{{3,4}}] ;
bkW = TensorContract[ Outer[Times, W, B, K] , {{1,4},{2,5}}] ;
Wbk = TensorContract [ Outer[Times, W, B, K] , {{2,4},{3,5}}];

Zero-th coefficient

Simplify[a0 - (-I cs^2)((B.K)^2(K.K * Det[shear] - 1/2 K.shear.K * Tr[shear2] + (B.K)^2 K.shear.K)+ Det[shear](K.shear2.K -1/2K.K * Tr[shear2]))]

there are no terms proportional to the speed of sound

a0nocs - a0 /. cs->0 ;
a0nocs

a0D = Simplify[a0 / (-I cs^2)];

partnoB = a0D /. Table[B[[i]]->0, {i,1,Length[B]}];

decomposition of the terms without B: quadratic in K and shear +

rest = Simplify[partnoB / Det[shear]];

shearvec = {s1,s2};
sterms = Outer[Times,shearvec, shearvec];
l = Length[sterms];
Sterms = sterms[[1]];
For[i=1,i<l,i++ ,
Sterms = Join[Sterms, sterms[[i+1]]]
]
DeleteDuplicates[Sterms];

kterms = Outer[Times,K, K];
l = Length[kterms];
Kterms = kterms[[1]];
For[i=1,i<l,i++ ,
Kterms = Join[Kterms, kterms[[i+1]]]
]
DeleteDuplicates[Kterms];

skterms = Outer[Times,Kterms,Sterms];
l = Length[skterms];
SKterms = skterms[[1]];
For[i=1,i<l,i++,
SKterms = Join[SKterms, skterms[[i+1]]];
]
TOTterms = DeleteDuplicates[SKterms];

sc1 = Expand[K.shear2.K];
sc2 = Expand[K.K * Tr[shear2]];

scalars = {sc1,sc2};
alphas = {c1,c2};

expr = Expand[ rest - scalars.alphas  ];


l = Length[TOTterms];
subs=  Table[TOTterms[[i]]->0 , {i,1,l}] ;
Subs = Table[Delete[subs,i],{i,1,l}];
Expr = Table[Collect[ ReplaceAll[expr,Subs[[i]]],TOTterms,Simplify],{i,1,l}];

Expr = ReplaceAll[Expr, {k1->1,k2->1,k3->1,s1->1,s2->1}];
Expr = DeleteCases[Expr,0];
Length[Expr] ;
If[ DuplicateFreeQ[Expr] , 
Expr = Expr , 
Expr = DeleteDuplicates[Expand[Expr]] ;
]
Length[Expr];
Solve[Table[Expr[[i]] == 0 , {i,1,Length[Expr]}], alphas]

Simplify[partnoB - Det[shear](K.shear2.K -1/2K.K * Tr[shear2])];

partB = a0csD  - (a0csD /. Table[B[[i]]->0, {i,1,Length[B]}]);

rest = Simplify[partB / (B.K)^2];
restlinS =  Expand[(Expand[ Coefficient[rest,s1]] /. {s1->0, s2->0 })*s1+ (Expand[ Coefficient[rest,s2]] /.{s1->0, s2->0 } )*s2] ;
restcubS = Expand[(Expand[ Coefficient[rest,s1^3]] /. s2->0 )*s1^3 +  (Expand[ Coefficient[rest,s1^2 s2]] /. s2->0 )*s1^2 s2+ (Expand[ Coefficient[rest,s1 s2^2]] /. s1->0 )*s1 s2^2  + (Expand[ Coefficient[rest,s2^3]] /. s1->0 )*s2^3] ;
Expand[rest - restlinS - restcubS];

decomposition of the part cubic in shear: quadratic in K

shearvec = {s1,s2};
sterms = Outer[Times,shearvec, shearvec];
l = Length[sterms];
Sterms = sterms[[1]];
For[i=1,i<l,i++ ,
Sterms = Join[Sterms, sterms[[i+1]]]
]
DeleteDuplicates[Sterms];

sterms = Outer[Times,shearvec,Sterms];
l = Length[sterms];
Sterms = sterms[[1]];
For[i=1,i<l,i++ ,
Sterms = Join[Sterms, sterms[[i+1]]]
]
DeleteDuplicates[Sterms];

kterms = Outer[Times,K, K];
l = Length[kterms];
Kterms = kterms[[1]];
For[i=1,i<l,i++ ,
Kterms = Join[Kterms, kterms[[i+1]]]
]
DeleteDuplicates[Kterms];

skterms = Outer[Times,Kterms,Sterms];
l = Length[skterms];
SKterms = skterms[[1]];
For[i=1,i<l,i++,
SKterms = Join[SKterms, skterms[[i+1]]];
]
TOTterms = DeleteDuplicates[SKterms];

sc1 = Expand[K.shear.K * Tr[shear2]];
sc2 = Expand[K.K * Det[shear]];

scalars = {sc1,sc2};
alphas = {c1,c2};

expr = Expand[restcubS- scalars.alphas  ];


l = Length[TOTterms];
subs=  Table[TOTterms[[i]]->0 , {i,1,l}] ;
Subs = Table[Delete[subs,i],{i,1,l}];
Expr = Table[Collect[ ReplaceAll[expr,Subs[[i]]],TOTterms,Simplify],{i,1,l}];

Expr = ReplaceAll[Expr, {k1->1,k2->1,k3->1,s1->1,s2->1}];
Expr = DeleteCases[Expr,0];
Length[Expr] ;
If[ DuplicateFreeQ[Expr] , 
Expr = Expr , 
Expr = DeleteDuplicates[Expand[Expr]] ;
]
Length[Expr];
Solve[Table[Expr[[i]] == 0 , {i,1,Length[Expr]}], alphas]

Simplify[restcubS - (K.K * Det[shear] - 1/2 K.shear.K * Tr[shear2])];

Simplify[restlinS - (B.K)^2 K.shear.K];

Simplify[partB -  (B.K)^2(K.K * Det[shear] - 1/2 K.shear.K * Tr[shear2] + (B.K)^2 K.shear.K)];

First coefficient

Simplify[a1 - cs^2((B.K)^4 K.K + (B.K)^2(K.shear2.K - Tr[shear2] * K.K)+ Det[shear]* (K.shear.K) + 1/4 Tr[shear2]^2 K.K - 1/2 (K.shear2.K) * Tr[shear2]) -(B.shear2.B* K.shear2.K + 1/2Tr[shear2] * (U.shear2.U) - 1/4(Tr[shear2]^2)* U.U -Det[shear] (U.shear.U + K.shear.K * B.B+ B.shear.B * K.K)- (B.K)^2 B.shear.B K.shear.K - Det[shear]^2)]

part with the speed of sound

a1cs = a1 -(a1/. cs->0);
a1csD = Simplify[a1cs / cs^2];

partnoS = a1csD /. {s1->0,s2->0};
partlinS =  Expand[(Expand[ Coefficient[a1csD,s1]] /. {s1->0, s2->0 })*s1+ (Expand[ Coefficient[a1csD,s2]] /.{s1->0, s2->0 } )*s2] ;
partquadS = (Expand[ Coefficient[a1csD,s1^2]] /. {s1->0, s2->0 } )* s1^2 + (Expand[ Coefficient[a1csD,s2^2]] /. {s1->0, s2->0 })* s2^2 + (Expand[ Coefficient[a1csD,s1 s2]] /. {s1^2 ->0, s1 s2->0, s2^2->0 })* s1 s2 ;
partcubS = Expand[(Expand[ Coefficient[a1csD,s1^3]] /. s2->0 )*s1^3 +  (Expand[ Coefficient[a1csD,s1^2 s2]] /. s2->0 )*s1^2 s2+ (Expand[ Coefficient[a1csD,s1 s2^2]] /. s1->0 )*s1 s2^2  + (Expand[ Coefficient[a1csD,s2^3]] /. s1->0 )*s2^3] ;
partquartS =Expand[ (Expand[ Coefficient[a1csD,s1^4]] )*s1^4  + (Expand[ Coefficient[a1csD,s1^3 s2]] )*s1^3 s2 + (Expand[ Coefficient[a1csD,s1^2 s2^2]] )*s1^2 s2^2 + (Expand[ Coefficient[a1csD,s1 s2^3]] )*s1 s2^3 + (Expand[ Coefficient[a1csD,s2^4]] )*s2^4  ];

Simplify[partnoS - (B.K)^4 K.K];

decomposition of the terms quadratic in shear

rest = Simplify[partquadS/ (B.K)^2];

kterms = Outer[Times, K, K];
l = Length[kterms] ;
Kterms = kterms[[1]];
For[i=1,i<l, i++, 
Kterms = Join[Kterms, kterms[[i+1]]]
]
Kterms = DeleteDuplicates[Kterms];

shearvec = {s1,s2};
sterms = Outer[Times,shearvec,shearvec];
l = Length[sterms];
Sterms = sterms[[1]];
For[i=1,i<l, i++, 
Sterms = Join[Sterms, sterms[[i+1]]]
]
Sterms = DeleteDuplicates[Sterms] ;

totterms = Outer[Times, Kterms, Sterms];
l = Length[totterms] ;
TOTterms = totterms[[1]];
For[i=1,i<l, i++, 
TOTterms = Join[TOTterms, totterms[[i+1]]]
]
TOTterms = DeleteDuplicates[TOTterms];

sc1 = Expand[K.shear2.K];
sc2 = Expand[(K.shear.K) *( K.shear.K)];
sc3 = Expand[Tr[shear2] * K.K] ;

alphas = {a,b,c} ;
scalars = {sc1,sc2,sc3};
expr =Expand[ rest- alphas.scalars] ;

l = Length[TOTterms];
subs=  Table[TOTterms[[i]]->0 , {i,1,l}] ;
Subs = Table[Delete[subs,i],{i,1,l}];
Expr = Table[Collect[ ReplaceAll[expr,Subs[[i]]],TOTterms,Simplify],{i,1,l}];

Expr = ReplaceAll[Expr, {k1->1,k2->1,k3->1,s1->1,s2->1}];
Expr = DeleteCases[Expr,0];
Length[Expr] ;
If[ DuplicateFreeQ[Expr] , 
 , 
Expr = DeleteDuplicates[Expand[Expr]] ;
]
Length[Expr]; Solve[Expr==0, alphas]

Simplify[partquadS- (B.K)^2(K.shear2.K - Tr[shear2] * K.K)];

decomposition of the terms quartic in shear: quadratic in K

kterms = Outer[Times, K, K];
l = Length[kterms] ;
Kterms = kterms[[1]];
For[i=1,i<l, i++, 
Kterms = Join[Kterms, kterms[[i+1]]]
]
Kterms = DeleteDuplicates[Kterms];

shearvec = {s1,s2};
sterms = Outer[Times,shearvec,shearvec];
l = Length[sterms];
Sterms = sterms[[1]];
For[i=1,i<l, i++, 
Sterms = Join[Sterms, sterms[[i+1]]]
]
Sterms = DeleteDuplicates[Sterms] ;

sterms = Outer[Times,Sterms,Sterms];
l = Length[sterms];
Sterms = sterms[[1]];
For[i=1,i<l, i++, 
Sterms = Join[Sterms, sterms[[i+1]]]
]
Sterms = DeleteDuplicates[Sterms] ;

totterms = Outer[Times, Kterms, Sterms];
l = Length[totterms] ;
TOTterms = totterms[[1]];
For[i=1,i<l, i++, 
TOTterms = Join[TOTterms, totterms[[i+1]]]
]
TOTterms = DeleteDuplicates[TOTterms];

sc1 = Expand[(K.shear2.K) * Tr[shear2]]  ;
sc2 = Expand[ Det[shear]* (K.shear.K)] ;
sc3 = Expand[Tr[shear2]^2 K.K];

scalars = {sc1, sc2,sc3};
alphas= {a,b,c};
expr = Expand[partquartS  - scalars.alphas  ];


l = Length[TOTterms];
subs=  Table[TOTterms[[i]]->0 , {i,1,l}] ;
Subs = Table[Delete[subs,i],{i,1,l}];
Expr = Table[Collect[ ReplaceAll[expr,Subs[[i]]],TOTterms,Simplify],{i,1,l}];

Expr = ReplaceAll[Expr, {k1->1,k2->1,k3->1,s1->1,s2->1}];
Expr = DeleteCases[Expr,0];
Length[Expr] ;
If[ DuplicateFreeQ[Expr] , 
Expr = Expr , 
Expr = DeleteDuplicates[Expand[Expr]] ;
]
Length[Expr];
Solve[Table[Expr[[i]] == 0 , {i,1,Length[Expr]}], alphas]

Simplify[partquartS - (Det[shear]* (K.shear.K) + 1/4 Tr[shear2]^2 K.K - 1/2 (K.shear2.K) * Tr[shear2])];

Simplify[a1cs - cs^2((B.K)^4 K.K + (B.K)^2(K.shear2.K - Tr[shear2] * K.K)+ Det[shear]* (K.shear.K) + 1/4 Tr[shear2]^2 K.K - 1/2 (K.shear2.K) * Tr[shear2])]

part without the speed of sound

a1nocs = a1 /. cs->0 ;

partnoB = a1nocs /. Table[B[[i]]->0, {i,1,Length[B]}] ;
Simplify[partnoB - (-Det[shear]^2)]

partB = a1nocs - (a1nocs /. Table[B[[i]]->0, {i,1,Length[B]}]);

partquadS = (Expand[ Coefficient[partB,s1^2]] /. {s1->0, s2->0 } )* s1^2 + (Expand[ Coefficient[partB,s2^2]] /. {s1->0, s2->0 })* s2^2 + (Expand[ Coefficient[partB,s1 s2]] /. {s1^2 ->0, s1 s2->0, s2^2->0 })* s1 s2 ;
partquartS =Expand[ (Expand[ Coefficient[partB,s1^4]] )*s1^4  + (Expand[ Coefficient[partB,s1^3 s2]] )*s1^3 s2 + (Expand[ Coefficient[partB,s1^2 s2^2]] )*s1^2 s2^2 + (Expand[ Coefficient[partB,s1 s2^3]] )*s1 s2^3 + (Expand[ Coefficient[partB,s2^4]] )*s2^4  ];
Expand[partB - partquadS - partquartS];

Simplify[partquadS - (- (B.K)^2 B.shear.B K.shear.K)] ;

decomposition of the terms quartic in shear

kterms = Outer[Times, K, K];
l = Length[kterms] ;
Kterms = kterms[[1]];
For[i=1,i<l, i++, 
Kterms = Join[Kterms, kterms[[i+1]]]
]
Kterms = DeleteDuplicates[Kterms];

bterms = Outer[Times, B, B];
l = Length[kterms] ;
Bterms = bterms[[1]];
For[i=1,i<l, i++, 
Bterms = Join[Bterms, bterms[[i+1]]]
]
Bterms = DeleteDuplicates[Bterms];

shearvec = {s1,s2};
sterms = Outer[Times,shearvec,shearvec];
l = Length[sterms];
Sterms = sterms[[1]];
For[i=1,i<l, i++, 
Sterms = Join[Sterms, sterms[[i+1]]]
]
Sterms = DeleteDuplicates[Sterms] ;

sterms = Outer[Times,Sterms,Sterms];
l = Length[sterms];
Sterms = sterms[[1]];
For[i=1,i<l, i++, 
Sterms = Join[Sterms, sterms[[i+1]]]
]
Sterms = DeleteDuplicates[Sterms] ;

totterms = Outer[Times, Kterms, Bterms];
l = Length[totterms] ;
TOTterms = totterms[[1]];
For[i=1,i<l, i++, 
TOTterms = Join[TOTterms, totterms[[i+1]]]
]
TOTterms = DeleteDuplicates[TOTterms];

totterms = Outer[Times, TOTterms, Sterms];
l = Length[totterms] ;
TOTterms = totterms[[1]];
For[i=1,i<l, i++, 
TOTterms = Join[TOTterms, totterms[[i+1]]]
]
TOTterms = DeleteDuplicates[TOTterms];

sc1 = Expand[B.shear2.B* K.shear2.K];
sc2 = Expand[B.shear2.B *K.K *Tr[shear2]];
sc3 = Expand[K.shear2.K *B.B *Tr[shear2]];
sc4 = Expand[Tr[shear2] * (K.shear.K )* (B.shear.B) ];
sc5 = Expand[Tr[shear2] * (K.shear.B)^2];
sc6 = Expand[Tr[shear2] * (U.shear2.U)];
sc7 = Expand[(Tr[shear2]^2)* B.B *K.K];
sc8 = Expand[(Tr[shear2]^2)* (B.K)^2];
sc9 = Expand[(Tr[shear2]^2)* U.U];
sc10 = Expand[Det[shear]*(U.shear.U)];
sc11 = Expand[Det[shear]*(K.shear.K) * B.B];
sc12 = Expand[Det[shear]*(B.shear.B) * K.K];

scalars = {sc1, sc2,sc3,0,sc5,sc6,0,sc8,sc9,sc10,sc11,sc12};
alphas= {c1,c2,c3,c4,c5,c6,c7,c8,c9,c10,c11,c12};
scalars.alphas;

expr = Expand[partquartS - scalars.alphas  ];

l = Length[TOTterms];
subs=  Table[TOTterms[[i]]->0 , {i,1,l}] ;
Subs = Table[Delete[subs,i],{i,1,l}];
Expr = Table[Collect[ ReplaceAll[expr,Subs[[i]]],TOTterms,Simplify],{i,1,l}];

Expr = ReplaceAll[Expr, {b1->1,b2->1,b3->1,k1->1,k2->1,k3->1,s1->1,s2->1}];
Expr = DeleteCases[Expr,0];
Length[Expr] ;
If[ DuplicateFreeQ[Expr] , 
Expr = Expr , 
Expr = DeleteDuplicates[Expand[Expr]] ;
]
Length[Expr];
Solve[Table[Expr[[i]] == 0 , {i,1,Length[Expr]}], alphas]

Simplify[partquartS - (B.shear2.B* K.shear2.K + 1/2Tr[shear2] * (U.shear2.U) - 1/4(Tr[shear2]^2)* U.U -Det[shear] (U.shear.U + K.shear.K * B.B+ B.shear.B * K.K))];
Simplify[a1nocs- (B.shear2.B* K.shear2.K + 1/2Tr[shear2] * (U.shear2.U) - 1/4(Tr[shear2]^2)* U.U -Det[shear] (U.shear.U + K.shear.K * B.B+ B.shear.B * K.K)- (B.K)^2 B.shear.B K.shear.K - Det[shear]^2)]

Second coefficient

Simplify[a2 - I cs^2(Det[shear]* (K.K) -1/2(K.shear.K) * Tr[shear2] + 2 (B.K)^2 K.shear.K)- (-I)( B.shear.B * (K.shear2.K) - K.shear.K * (B.shear2.B) + 1/2Tr[shear2](B.B * K.shear.K - K.K * B.shear.B)+  (B.K)^2((B . shear . B) * (K.K) -K.shear.K * (B.B) ))]

part with the speed of sound

a2cs = a2 - (a2/.cs->0);

a2csD = Expand[Simplify[a2cs / (I cs^2)]];
partnoB = a2csD /. Table[B[[i]]->0, {i,1,Length[B]}];
partB = a2csD - (a2csD /. Table[B[[i]]->0, {i,1,Length[B]}]);

decomposition of the part without B: quadratic in K and cubic in shear

kterms = Outer[Times, K, K];
l = Length[kterms] ;
Kterms = kterms[[1]];
For[i=1,i<l, i++, 
Kterms = Join[Kterms, kterms[[i+1]]]
]
Kterms = DeleteDuplicates[Kterms];

shearvec = {s1,s2};
sterms = Outer[Times,shearvec,shearvec];
l = Length[sterms];
Sterms = sterms[[1]];
For[i=1,i<l, i++, 
Sterms = Join[Sterms, sterms[[i+1]]]
]
Sterms = DeleteDuplicates[Sterms] ;

sterms = Outer[Times, Sterms,shearvec];
l = Length[sterms];
Sterms = sterms[[1]];
For[i=1,i<l, i++, 
Sterms = Join[Sterms, sterms[[i+1]]]
]
Sterms = DeleteDuplicates[Sterms] ;

totterms = Outer[Times, Kterms,Sterms];
l = Length[totterms];
TOTterms = totterms[[1]];
For[i=1,i<l, i++, 
TOTterms = Join[TOTterms, totterms[[i+1]]]
]
TOTterms =  DeleteDuplicates[TOTterms] ;

sc1 = Expand[(K.shear.K) * Tr[shear2]]  ;
sc2 = Expand[ Det[shear]* (K.K)] ;
scalars = {sc1, sc2};

alphas= {a,b};
expr = Expand[partnoB  - scalars.alphas  ];


l = Length[TOTterms];
subs=  Table[TOTterms[[i]]->0 , {i,1,l}] ;
Subs = Table[Delete[subs,i],{i,1,l}];
Expr = Table[Collect[ ReplaceAll[expr,Subs[[i]]],TOTterms,Simplify],{i,1,l}];

Expr = ReplaceAll[Expr, {k1->1,k2->1,k3->1,s1->1,s2->1}];
Expr = DeleteCases[Expr,0];
Length[Expr] ;
If[ DuplicateFreeQ[Expr] , 
Expr = Expr , 
Expr = DeleteDuplicates[Expand[Expr]] ;
]
Length[Expr];
Solve[Table[Expr[[i]] == 0 , {i,1,Length[Expr]}], alphas]

Simplify[partnoB - (Det[shear]* (K.K) -1/2(K.shear.K) * Tr[shear2] )];
Simplify[partB - 2 (B.K)^2 K.shear.K];
Simplify[a2cs - I cs^2(Det[shear]* (K.K) -1/2(K.shear.K) * Tr[shear2] + 2 (B.K)^2 K.shear.K)]

part without the speed of sound

a2nocs = a2/.cs->0;
a2nocsD =Expand[ Simplify[a2nocs / (-I)]];

partlinS =  Expand[(Expand[ Coefficient[a2nocsD,s1]] /. s2->0 )*s1+ (Expand[ Coefficient[a2nocsD,s2]] /. s1->0 )*s2] ;
partcubS = Expand[(Expand[ Coefficient[a2nocsD,s1^3]] /. s2->0 )*s1^3 +  (Expand[ Coefficient[a2nocsD,s1^2 s2]] /. s2->0 )*s1^2 s2+ (Expand[ Coefficient[a2nocsD,s1 s2^2]] /. (s1^2)->0 )*s1 s2^2  + (Expand[ Coefficient[a2nocsD,s2^3]] /. s1->0 )*s2^3] ;
Expand[a2nocsD - partlinS -partcubS ]

decomposition of the terms cubic in shear: quadratic in B and K

kterms = Outer[Times, K, K];
l = Length[kterms] ;
Kterms = kterms[[1]];
For[i=1,i<l, i++, 
Kterms = Join[Kterms, kterms[[i+1]]]
]
Kterms = DeleteDuplicates[Kterms];

bterms = Outer[Times, B, B];
l = Length[bterms] ;
Bterms = bterms[[1]];
For[i=1,i<l, i++, 
Bterms = Join[Bterms, bterms[[i+1]]]
]
Bterms = DeleteDuplicates[Bterms];

bkterms = Outer[Times,Bterms,Kterms];
l = Length[bkterms];
BKterms = bkterms[[1]];
For[i=1,i<l, i++, 
BKterms = Join[BKterms, bkterms[[i+1]]]
]
BKterms = DeleteDuplicates[BKterms];

shearvec = {s1,s2};
sterms = Outer[Times, shearvec,shearvec];
l = Length[sterms];
Sterms = sterms[[1]];
For[i=1,i<l, i++, 
Sterms = Join[Sterms, sterms[[i+1]]]
]
Sterms =  DeleteDuplicates[Sterms] ;

sterms = Outer[Times, shearvec,Sterms];
l = Length[sterms];
Sterms = sterms[[1]];
For[i=1,i<l, i++, 
Sterms = Join[Sterms, sterms[[i+1]]]
]
Sterms =  DeleteDuplicates[Sterms] ;

totterms = Outer[Times,Sterms,BKterms];
l = Length[totterms];
TOTterms = totterms[[1]];
For[i=1,i<l, i++, 
TOTterms = Join[TOTterms, totterms[[i+1]]]
]
TOTterms = DeleteDuplicates[TOTterms];

sc1 = Expand[Det[shear] * (B.K)^2 ]  ;
sc2 = Expand[Det[shear] * (U.U) ]  ;
sc3 = Expand[Det[shear] * (B.B)(K.K) ]  ;
sc4 = Expand[K.shear.K * (B.shear2.B) ]  ;
sc5 = Expand[B.shear.B * (K.shear2.K) ]  ;
sc6 = Expand[U.shear.U * Tr[shear2]];
sc7 = Expand[K.shear.K * (B.B) * Tr[shear2]];
sc8 = Expand[B.shear.B * (K.K) * Tr[shear2]];


scalars = {0, sc2,sc3,sc4,sc5,sc6,sc7,sc8};

alphas= {c1,c2,c3,c4,c5,c6,c7,c8};
expr = Expand[partcubS - scalars.alphas  ];

l = Length[TOTterms];
subs=  Table[TOTterms[[i]]->0 , {i,1,l}] ;
Subs = Table[Delete[subs,i],{i,1,l}];
Expr = Table[Collect[ ReplaceAll[expr,Subs[[i]]],TOTterms,Simplify],{i,1,l}];

Expr = ReplaceAll[Expr, {b1->1,b2->1,b3->1,k1->1,k2->1,k3->1,s1->1,s2->1}];
Expr = DeleteCases[Expr,0];
Length[Expr] ;
If[ DuplicateFreeQ[Expr] , 
Expr = Expr , 
Expr = DeleteDuplicates[Expand[Expr]] ;
]
Length[Expr];
Solve[Table[Expr[[i]] == 0 , {i,1,Length[Expr]}], alphas]

Simplify[partcubS - 1/2(2 B.shear.B * (K.shear2.K) - 2K.shear.K * (B.shear2.B) + K.shear.K * (B.B) * Tr[shear2]- B.shear.B * (K.K) * Tr[shear2])];

decomposition of the terms linear in shear: quartic in B and K

rest  = Expand[Simplify[partlinS / (B.K)^2]];

bterms = Outer[Times, B,B];
l = Length[bterms];
Bterms = bterms[[1]];
For[i=1,i<l, i++, 
Bterms = Join[Bterms, bterms[[i+1]]]
]
kterms = Outer[Times, K, K];
l = Length[kterms] ;
Kterms = kterms[[1]];
For[i=1,i<l, i++, 
Kterms = Join[Kterms, kterms[[i+1]]]
]
Kterms = DeleteDuplicates[Kterms];
Bterms = DeleteDuplicates[Bterms];
toterms =  Outer[Times,Bterms,Kterms];
l = Length[toterms];
TOTerms = toterms[[1]];
For[i=1,i<l, i++, 
TOTerms = Join[TOTerms, toterms[[i+1]]]
]
TOTerms = DeleteDuplicates[TOTerms] ;
shearvec = {s1,s2};
toterms = Outer[Times,TOTerms,shearvec];
l = Length[toterms];
TOTerms = toterms[[1]];
For[i=1,i<l, i++, 
TOTerms = Join[TOTerms, toterms[[i+1]]]
]
TOTerms = DeleteDuplicates[TOTerms] ;

sc1 = Expand[U.shear.U]  ;
sc2 = Expand[ (B . shear . K)* (B.K)] ;
sc3 = Expand[K.shear.K * (B.B)];
sc4= Expand[ (B . shear . B) * (K.K) ] ; 


scalars = {0, sc2,sc3,sc4};

alphas= {a,b,c,d};
expr = Expand[rest - scalars.alphas  ];


l = Length[TOTerms];
subs=  Table[TOTerms[[i]]->0 , {i,1,l}] ;
Subs = Table[Delete[subs,i],{i,1,l}];
Expr = Table[Collect[ ReplaceAll[expr,Subs[[i]]],TOTerms,Simplify],{i,1,l}];

Expr = ReplaceAll[Expr, {b1->1,b2->1,b3->1,k1->1,k2->1,k3->1,s1->1,s2->1}];
Expr = DeleteCases[Expr,0];
Length[Expr] ;
If[ DuplicateFreeQ[Expr] , 
Expr = Expr , 
Expr = DeleteDuplicates[Expand[Expr]] ;
]
Length[Expr];
Solve[Table[Expr[[i]] == 0 , {i,1,Length[Expr]}], alphas]

Simplify[partlinS - (B.K)^2((B . shear . B) * (K.K) -K.shear.K * (B.B) )];

Simplify[a2nocs - (-I)(1/2(2 B.shear.B * (K.shear2.K) - 2K.shear.K * (B.shear2.B) + K.shear.K * (B.B) * Tr[shear2]- B.shear.B * (K.K) * Tr[shear2])+  (B.K)^2((B . shear . B) * (K.K) -K.shear.K * (B.B) ))];

Third coefficient

Simplify[a3 - ( cs^2(Tr[shear2] * K.K  - K.shear2.K-2 (B.K)^2 K.K)+ (K . shear . K)* (B . shear . B) + U. shear2 .U + Tr[shear2]* (B.K)^2- B.B *K.K(B.K)^2 - 1/4Tr[shear2]^2)]

part with the speed of sound

a3cs = a3 - (a3 /. cs->0 );

partnoB = a3cs /. Table[ B[[i]] ->0 ,{i,1,Length[B]}];

decomposition of the terms with cs but no B: quadratic in K and shear

kterms = Outer[Times, K, K];
l = Length[kterms] ;
Kterms = kterms[[1]];
For[i=1,i<l, i++, 
Kterms = Join[Kterms, kterms[[i+1]]]
]
Kterms = DeleteDuplicates[Kterms];

shearvec = {s1,s2};
sterms = Outer[Times,shearvec,shearvec];
l = Length[sterms];
Sterms = sterms[[1]];
For[i=1,i<l, i++, 
Sterms = Join[Sterms, sterms[[i+1]]]
]
Sterms = DeleteDuplicates[Sterms] ;

totterms = Outer[Times, Kterms, Sterms];
l = Length[totterms] ;
TOTterms = totterms[[1]];
For[i=1,i<l, i++, 
TOTterms = Join[TOTterms, totterms[[i+1]]]
]
TOTterms = DeleteDuplicates[TOTterms];

sc1 = Expand[K.shear2.K];
sc2 = Expand[(K.shear.K) *( K.shear.K)];
sc3 = Expand[Tr[shear2] * K.K] ;

alphas = {a,b,c} ;
scalars = {sc1,sc2,sc3};
partnoBD = Simplify[partnoB /(cs^2)];
expr =Expand[ partnoBD- alphas.scalars] ;

l = Length[TOTterms];
subs=  Table[TOTterms[[i]]->0 , {i,1,l}] ;
Subs = Table[Delete[subs,i],{i,1,l}];
Expr = Table[Collect[ ReplaceAll[expr,Subs[[i]]],TOTterms,Simplify],{i,1,l}];

Expr = ReplaceAll[Expr, {k1->1,k2->1,k3->1,s1->1,s2->1}];
Expr = DeleteCases[Expr,0];
Length[Expr] ;
If[ DuplicateFreeQ[Expr] , 
 , 
Expr = DeleteDuplicates[Expand[Expr]] ;
]
Length[Expr]; Solve[Expr==0, alphas]

Simplify[partnoB - cs^2(Tr[shear2] * K.K  - K.shear2.K)];

partB = a3cs - (a3cs /. Table[ B[[i]] ->0 ,{i,1,Length[B]}]);
Simplify[partB  - (-2cs^2 (B.K)^2 K.K)];

Simplify[a3cs - cs^2(Tr[shear2] * K.K  - K.shear2.K-2 (B.K)^2 K.K)]

part without the speed of sound

a3nocs = a3 /. cs->0 ;

partnoB = a3nocs /. Table[ B[[i]] ->0 ,{i,1,Length[B]}];
Simplify[partnoB - (- (1/2Tr[shear2])^2)];

partB =  a3nocs  - ( a3nocs /. Table[ B[[i]] ->0 ,{i,1,Length[B]}]);
partBnoS = partB /. {s1 ->0,s2->0};
Simplify[partBnoS  - (- B.B (B.K)^2 K.K)];

partBS = partB -(partB /. {s1 ->0,s2->0} );

decomposition: it is quadratic in B and K and shear

bterms = Outer[Times, B,B];
l = Length[bterms];
Bterms = bterms[[1]];
For[i=1,i<l, i++, 
Bterms = Join[Bterms, bterms[[i+1]]]
]
kterms = Outer[Times, K, K];
l = Length[kterms] ;
Kterms = kterms[[1]];
For[i=1,i<l, i++, 
Kterms = Join[Kterms, kterms[[i+1]]]
]
shearvec = {s1,s2};
sterms = Outer[Times, shearvec,shearvec];
l = Length[sterms] ;
Sterms = sterms[[1]] ;
For[i=1,i<l, i++, 
Sterms = Join[Sterms, sterms[[i+1]]]
]
Kterms = DeleteDuplicates[Kterms];
Bterms = DeleteDuplicates[Bterms];
Sterms = DeleteDuplicates[Sterms];
toterms =  Outer[Times,Bterms,Kterms];
l = Length[toterms];
TOTerms = toterms[[1]];
For[i=1,i<l, i++, 
TOTerms = Join[TOTerms, toterms[[i+1]]]
]
TOTerms = DeleteDuplicates[TOTerms] ;
toterms = Outer[Times,TOTerms,Sterms];
l = Length[toterms];
TOTerms = toterms[[1]];
For[i=1,i<l, i++, 
TOTerms = Join[TOTerms, toterms[[i+1]]]
]
TOTerms = DeleteDuplicates[TOTerms] ;

sc1 = Expand[ (K . shear . B)^2 ] ;
sc2 =  Expand[(K . shear . K)* (B . shear . B) ] ;
trace2 = Tr[shear2];
sc3 = Expand[U. shear2 .U]  ;
sc4 = Expand[(B. shear2.K) * (B.K)] ;
sc5 = Expand[ trace2 * (B.K)^2] ;
sc6 = Expand[ trace2 * (U.U)];

levi = LeviCivitaTensor[3] ;
W =TensorContract[  Outer[Times, levi, shear] ,{{3,4}}] ;
bkW = TensorContract[ Outer[Times, W, B, K] , {{1,4},{2,5}}] ;
Wbk = TensorContract [ Outer[Times, W, B, K] , {{2,4},{3,5}}];
sc7 = Expand[ bkW . Wbk ];

scalars = {0,sc2,sc3,sc4,sc5,sc6,sc7};
alphas= {a,b,c,d,e,f,g};
scalars.alphas;
expr =  partBS- scalars.alphas  ;

l = Length[TOTerms];
subs=  Table[TOTerms[[i]]->0 , {i,1,l}] ;
Subs = Table[Delete[subs,i],{i,1,l}];
Expr = Table[Collect[ ReplaceAll[expr,Subs[[i]]],TOTerms,Simplify],{i,1,l}];

Expr = ReplaceAll[Expr, {b1->1,b2->1,b3->1,k1->1,k2->1,k3->1,s1->1,s2->1}];
Expr = DeleteCases[Expr,0];
Length[Expr] ;
If[ DuplicateFreeQ[Expr] , 
 , 
Expr = DeleteDuplicates[Expand[Expr]] ;
]
Length[Expr];
Solve[Expr==0, alphas]

Simplify[partBS - ((K . shear . K)* (B . shear . B) + U. shear2 .U + Tr[shear2]* (B.K)^2)];

Simplify[partB - ((K . shear . K)* (B . shear . B) + U. shear2 .U + Tr[shear2]* (B.K)^2- B.B (B.K)^2 K.K)];

Simplify[a3nocs - ((K . shear . K)* (B . shear . B) + U. shear2 .U + Tr[shear2] * (B.K)^2- B.B (B.K)^2 K.K- (1/2Tr[shear2])^2)];

Fourth coefficient

Simplify[a4 - (-I)(cs^2 K.shear.K + (K . shear . K)* (B.B)-(B.shear.B) * (K.K))]

ia4 = Simplify[I a4];
ia4cs = ia4 - (ia4 /. cs->0);
Simplify[ia4cs - cs^2 K.shear.K];

ia4nocs = ia4 /. cs->0;

decomposition: the part without cs is quadratic in B and K and linear in the shear

bterms = Outer[Times, B,B];
l = Length[bterms];
Bterms = bterms[[1]];
For[i=1,i<l, i++, 
Bterms = Join[Bterms, bterms[[i+1]]]
]
kterms = Outer[Times, K, K];
l = Length[kterms] ;
Kterms = kterms[[1]];
For[i=1,i<l, i++, 
Kterms = Join[Kterms, kterms[[i+1]]]
]
Kterms = DeleteDuplicates[Kterms];
Bterms = DeleteDuplicates[Bterms];
toterms =  Outer[Times,Bterms,Kterms];
l = Length[toterms];
TOTerms = toterms[[1]];
For[i=1,i<l, i++, 
TOTerms = Join[TOTerms, toterms[[i+1]]]
]
TOTerms = DeleteDuplicates[TOTerms] ;
shearvec = {s1,s2};
toterms = Outer[Times,TOTerms,shearvec];
l = Length[toterms];
TOTerms = toterms[[1]];
For[i=1,i<l, i++, 
TOTerms = Join[TOTerms, toterms[[i+1]]]
]
TOTerms = DeleteDuplicates[TOTerms] ;

sc1 = Expand[U.shear.U]  ;
sc2 = Expand[ (B . shear . K)* (B.K)] ;
sc3 = Expand[(K . shear . K)* (B . B)];
sc4= Expand[ (B . shear . B) * (K.K) ] ; 


scalars = {0, sc2,sc3,sc4};

alphas= {a,b,c,d};
expr = Expand[ia4nocs- scalars.alphas  ];


l = Length[TOTerms];
subs=  Table[TOTerms[[i]]->0 , {i,1,l}] ;
Subs = Table[Delete[subs,i],{i,1,l}];
Expr = Table[Collect[ ReplaceAll[expr,Subs[[i]]],TOTerms,Simplify],{i,1,l}];

Expr = ReplaceAll[Expr, {b1->1,b2->1,b3->1,k1->1,k2->1,k3->1,s1->1,s2->1}];
Expr = DeleteCases[Expr,0];
Length[Expr] ;
If[ DuplicateFreeQ[Expr] , 
Expr = Expr , 
Expr = DeleteDuplicates[Expand[Expr]] ;
]
Length[Expr];
Solve[Table[Expr[[i]] == 0 , {i,1,Length[Expr]}], alphas]

Simplify[ia4nocs - ((K . shear . K)* (B . B)-(B . shear . B) * (K.K))]

Fifth coefficient

Simplify[a5 - (cs^2 K.K  - Tr[shear2]+ 2 (B.K)^2 + U.U)]

a5cs = a5 - (a5/. cs->0) ;
Simplify[a5cs  - cs^2 K.K];

a5nocs = a5/. cs->0 ; 
partnoB = a5nocs /. Table[ B[[i]] ->0 ,{i,1,Length[B]}];
Simplify[partnoB - (- Tr[shear2])];

partB = a5nocs -(a5nocs /. Table[ B[[i]] ->0 ,{i,1,Length[B]}]);

decomposition of the part with B and no cs

bterms = Outer[Times, B,B];
l = Length[bterms];
Bterms = bterms[[1]];
For[i=1,i<l, i++, 
Bterms = Join[Bterms, bterms[[i+1]]]
]
Bterms = DeleteCases[Bterms,0];
Bterms = DeleteDuplicates[Bterms];

kterms = Outer[Times, K,K];
l = Length[bterms];
Kterms = kterms[[1]];
For[i=1,i<l, i++, 
Kterms = Join[Kterms, kterms[[i+1]]]
]
Kterms = DeleteDuplicates[Kterms];

bkterms = Outer[Times, Bterms,Kterms];
BKterms = bkterms[[1]];
For[ i = 1, i < Length[bkterms], i++ ,
BKterms = Join[BKterms, bkterms[[i+1]]]
]
TOTterms = DeleteDuplicates[BKterms];

sc1 = Expand[(B.K)^2];
sc2 = Expand[U.U];

alphas = {a,b} ;
scalars = {sc1,sc2};
expr =Expand[ partB - alphas.scalars] ;

l = Length[BKterms];
subs=  Table[BKterms[[i]]->0 , {i,1,l}] ;
Subs = Table[Delete[subs,i],{i,1,l}];
Expr = Table[Collect[ ReplaceAll[expr,Subs[[i]]],TOTterms,Simplify],{i,1,l}];

Expr = ReplaceAll[Expr, {bW->1,bq->1,kW->1,kq->1,ks->1}];
Expr = DeleteCases[Expr,0];
Length[Expr] ;
If[ DuplicateFreeQ[Expr] , 
 , 
Expr = DeleteDuplicates[Expand[Expr]] ;
]
Length[Expr]; Solve[Expr==0, alphas]

Simplify[partB - (2 (B.K)^2 + U.U)];

Background with vorticity

Dispersion relation

vortmat = {{-om,ro kW , ro kq, ro ks, 0,0,0},{cs^2 / ro kW, -om,0,0,-bq kq, bq kW,0},{cs^2 / ro kq ,0 , - om , I  W, bW kq,-bW kW, 0},{cs^2 / ro ks, 0, -I  W, - om, bW ks , bq ks , - bW kW - bq kq},{0,-bq kq,  bW kq , bW ks, -om , 0 ,0},{0, bq kW , - bW kW, bq ks,0, - om, - I W},{0, 0,0,- bW kW - bq kq, 0,  I W , -om}} ;
vortmat // MatrixForm 
vortdet = Simplify[Det[vortmat]];

a0 = ReplaceAll[Expand[vortdet], {om^7->0,om^6 ->0,om^5->0,om^4->0,om^3->0,om^2->0,om->0}] ;
a1 = ReplaceAll[Expand[vortdet] - a0, {om^7->0,om^6 ->0,om^5->0,om^4->0,om^3->0,om^2->0}] /. om ->1;
a2 = ReplaceAll[Expand[vortdet] - a0 - a1 om, {om^7->0,om^6 ->0,om^5->0,om^4->0,om^3->0}] /. om ->1;
a3 = ReplaceAll[Expand[vortdet] - a0 - a1 om - a2 om^2, {om^7->0,om^6 ->0,om^5->0,om^4->0}] /. om ->1;
a4 = ReplaceAll[Expand[vortdet] - a0 - a1 om - a2 om^2 - a3 om^3, {om^7->0,om^6 ->0,om^5->0}] /. om ->1;
a5 = ReplaceAll[Expand[vortdet] - a0 - a1 om - a2 om^2 - a3 om^3 - a4 om^4, {om^7->0,om^6 ->0}] /. om ->1;
a6 = ReplaceAll[Expand[vortdet] - a0 - a1 om - a2 om^2 - a3 om^3 - a4 om^4 - a5 om^5, {om^7->0}] /. om ->1;
a7 =(Expand[vortdet] - a0 - a1 om - a2 om^2 - a3 om^3 - a4 om^4 - a5 om^5 - a6 om^6)  /. om ->1;

coefflist = {a0,a1,a2,a3,a4,a5,a6,a7};
For[i=1,i < Length[coefflist]+1,i++,If[coefflist[[i]]==0, Print["0: ",i-1]] ]
For[i=1,i < Length[coefflist]+1,i++,If[coefflist[[i]]==1, Print["1: ",i-1]] ]
For[i=1,i < Length[coefflist]+1,i++,If[coefflist[[i]]==-1, Print["-1: ",i-1]] ]

All even coefficients are zero, and the seventh is -1. Save the rest.

SetDirectory[NotebookDirectory[]];
DumpSave["vortbg_coeff.mx", {a1,a3,a5}];

Coefficients decomposition

SetDirectory[NotebookDirectory[]];
<< "vortbg_coeff.mx";
B= {bW,bq,0};
K = {kW,kq,ks};
vort= {W,0,0};

First coefficient

Simplify[a1 -cs^2((B.K)^2 + vort.vort)*((K.vort)^2 + (B.K)^2 K.K) - (Cross[B,vort] .Cross[K,vort])^2]

a1nocs = a1 /. cs->0;
Simplify[a1nocs - (Cross[B,vort] .Cross[K,vort])^2];

a1cs = a1 - (a1 /. cs->0);

rest = Expand[Simplify[a1cs / (cs^2((B.K)^2 + vort.vort))]];

restnoW = rest /. W->0;
Simplify[restnoW  - ((B.K)^2 K.K)];

restW = rest - (rest /. W->0);
Simplify[restW - (K.vort)^2];

Simplify[a1cs- cs^2((B.K)^2 + vort.vort)((K.vort)^2 + (B.K)^2 K.K)];

Third coefficient

Simplify[a3 - (-cs^2(K.K vort.vort + (K.vort)^2+2(B.K)^2 K.K)-(vort.vort)^2- Cross[Cross[K,vort],B].Cross[Cross[K,vort],B] -Cross[Cross[vort,B],K].Cross[Cross[vort,B],K]-2(B.K) * (B.vort)* (vort.K)- B.B * (B.K)^2 * K.K)]

part with the speed of sound

a3cs = a3 - (a3/. cs->0);

partW = a3cs - (a3cs /. W->0);
Simplify[partW - (-cs^2 (K.K vort.vort + (K.vort)^2))];

partnoW = a3cs /. W->0 ;
Simplify[partnoW - (-2cs^2 ((B.K)^2 K.K))];

Simplify[a3cs - (-cs^2)(K.K vort.vort + (K.vort)^2+2(B.K)^2 K.K)]

part without the speed of sound

a3nocs = a3 /. cs->0 ;
Expand[%]

partnoW = a3nocs /. W->0 ;
Simplify[partnoW - (- B.B * (B.K)^2 * K.K)]

partW = a3nocs - (a3nocs /. W->0 );
coeff = CoefficientList[partW,W] ;
partWquad = coeff[[3]] W^2;
partWquart = coeff[[5]] W^4;

Simplify[partWquart - (-(vort.vort)^2)];

decomposition of the part quadratic in the vorticity: quadratic in B and K

bterms = Outer[Times, B,B];
l = Length[bterms];
Bterms = bterms[[1]];
For[i=1,i<l, i++, 
Bterms = Join[Bterms, bterms[[i+1]]]
]
Bterms = DeleteCases[Bterms,0];
Bterms = DeleteDuplicates[Bterms];

kterms = Outer[Times, K,K];
l = Length[bterms];
Kterms = kterms[[1]];
For[i=1,i<l, i++, 
Kterms = Join[Kterms, kterms[[i+1]]]
]
Kterms = DeleteDuplicates[Kterms];

bkterms = Outer[Times, Bterms,Kterms];
BKterms = bkterms[[1]];
For[ i = 1, i < Length[bkterms], i++ ,
BKterms = Join[BKterms, bkterms[[i+1]]]
]
BKterms = DeleteDuplicates[BKterms];

wterms = Outer[Times, vort,vort];
l = Length[wterms];
Wterms = wterms[[1]];
For[i=1,i<l, i++, 
Wterms = Join[Wterms, wterms[[i+1]]]
]
Wterms = DeleteCases[Wterms,0];
Wterms = DeleteDuplicates[Wterms];


bkwterms = Outer[Times, BKterms,Wterms];
l = Length[bkwterms];
BKWterms = bkwterms[[1]];
For[i=1,i<l, i++, 
BKWterms = Join[BKWterms, bkwterms[[i+1]]]
]
BKWWterms = DeleteCases[BKWterms,0];
TOTterms = DeleteDuplicates[BKWterms];

sc1 = Expand[Cross[Cross[B,K],vort].Cross[Cross[B,K],vort]];
sc2 = Expand[Cross[Cross[K,vort],B].Cross[Cross[K,vort],B]];
sc3 = Expand[Cross[Cross[vort,B],K].Cross[Cross[vort,B],K]];
sc4 = Expand[(B.K) * (B.vort)* (vort.K)];

alphas = {a,b,c,d} ;
scalars = {sc1,sc2,sc3,sc4};
expr =Expand[ partWquad - alphas.scalars] ;

l = Length[BKterms];
subs=  Table[BKterms[[i]]->0 , {i,1,l}] ;
Subs = Table[Delete[subs,i],{i,1,l}];
Expr = Table[Collect[ ReplaceAll[expr,Subs[[i]]],TOTterms,Simplify],{i,1,l}];

Expr = ReplaceAll[Expr, {bW->1,bq->1,kW->1,kq->1,ks->1,W->1}];
Expr = DeleteCases[Expr,0];
Length[Expr] ;
If[ DuplicateFreeQ[Expr] , 
 , 
Expr = DeleteDuplicates[Expand[Expr]] ;
]
Length[Expr]; Solve[Expr==0, alphas]

Simplify[partWquad - (- Cross[Cross[K,vort],B].Cross[Cross[K,vort],B] -Cross[Cross[vort,B],K].Cross[Cross[vort,B],K]-2(B.K) * (B.vort)* (vort.K))];

Simplify[a3nocs - (-(vort.vort)^2- Cross[Cross[K,vort],B].Cross[Cross[K,vort],B] -Cross[Cross[vort,B],K].Cross[Cross[vort,B],K]-2(B.K) * (B.vort)* (vort.K)- B.B * (B.K)^2 * K.K)];

Fifth coefficient

Simplify[a5 - (cs^2 K.K + 2(B.K)^2 + Cross[B,K].Cross[B,K] + 2 vort.vort )]

a5cs = a5 - (a5 /. cs->0);
Simplify[a5cs - cs^2 K.K];

a5nocs = a5 /. cs->0;
partnoW = a5nocs  /. W->0;

decomposition of the terms without vorticity

bterms = Outer[Times, B,B];
l = Length[bterms];
Bterms = bterms[[1]];
For[i=1,i<l, i++, 
Bterms = Join[Bterms, bterms[[i+1]]]
]
Bterms = DeleteCases[Bterms,0];
Bterms = DeleteDuplicates[Bterms];

kterms = Outer[Times, K,K];
l = Length[bterms];
Kterms = kterms[[1]];
For[i=1,i<l, i++, 
Kterms = Join[Kterms, kterms[[i+1]]]
]
Kterms = DeleteDuplicates[Kterms];

bkterms = Outer[Times, Bterms,Kterms];
BKterms = bkterms[[1]];
For[ i = 1, i < Length[bkterms], i++ ,
BKterms = Join[BKterms, bkterms[[i+1]]]
]
TOTterms = DeleteDuplicates[BKterms];

sc1 = Expand[(B.K)^2];
sc2 = Expand[Cross[B,K].Cross[B,K]];

alphas = {a,b} ;
scalars = {sc1,sc2};
expr =Expand[ partnoW - alphas.scalars] ;

l = Length[BKterms];
subs=  Table[BKterms[[i]]->0 , {i,1,l}] ;
Subs = Table[Delete[subs,i],{i,1,l}];
Expr = Table[Collect[ ReplaceAll[expr,Subs[[i]]],TOTterms,Simplify],{i,1,l}];

Expr = ReplaceAll[Expr, {bW->1,bq->1,kW->1,kq->1,ks->1}];
Expr = DeleteCases[Expr,0];
Length[Expr] ;
If[ DuplicateFreeQ[Expr] , 
 , 
Expr = DeleteDuplicates[Expand[Expr]] ;
]
Length[Expr]; Solve[Expr==0, alphas]

Simplify[a5nocs- (2(B.K)^2 + Cross[B,K].Cross[B,K] + 2 vort.vort)];