I know nothing about the range of the data for your problem. Assuming that the values for Rho, Mu and z are real numbers you can gain insight into your problem by combining bbgodfrey's comment with a plot using Manipulate.
For example, if Rho and Mu are known parameters you can see how the solution varies as you change the value of z.
Manipulate[
...

Although Belisarius' creative solution is entirely satisfactory, a solution symbolic at every step may be useful. To begin, define
x[t_] := -9 Sin[2 t] - 5 Sin[3 t]
y[t_] := 9 Cos[2 t] - 5 Cos[3 t]
and note that t = π corresponds to the uppermost point in the star in the question, {0, 14}}. From there, the point {0, -5} can be reached by increasing or ...

If you think you will run into that ambiguity a lot you might want to define a function to carry out your calculation:
Clear[matrixmult]
matrixmult[m_?MatrixQ, v_?VectorQ] := v.m.v
matrixmult[m_?NumberQ, v_?VectorQ] := v.(m IdentityMatrix[3]).v
I used two conditional definitions: the correct definitions will be picked depending on the type of the ...

I can't open the 9MB file as it is a .7z extension and I don't know what that is.
However, if your expressions has sets of Denominators that are the same then the following will spread the Simplify over your kernels by ParallelMaping it.
tmp is the first 5 summands you shared in PasteBin
Simplify@Total@ParallelMap[Simplify[Total[#]] &, GatherBy[List ...

Mathematica is a registered trademark of Wolfram Research, Inc. While the mark is used herein with the limited permission of Wolfram Research, Stack Exchange and this site disclaim all affiliation therewith.