NSolve[Rationalize[f1[x] == f2[x], 0], x, Reals], 100]
yields three solutions (with or without
N), which is the minimum number of solutions if
Zl ρ is positive.
The following, which sets the precision of the input to match the working precision,
NSolve[SetPrecision[f1[x] == f2[x], 100], x, Reals, WorkingPrecision -> 100]
also yields three solutions, after a while, which agree to 26 digits or so with the first code. (Using
SetPrecision[f1[x] == f2[x], Infinity] gives an equivalent result much faster.)
Precision: The input is basically
MachinePrecision. The OP's code,
NSolve[f1[x] == f2[x], x, Reals, WorkingPrecision -> 100] (with
Reals added per OP's comment), is asking for a working precision that is higher than the input precision. The output consists of three accurate solutions at
MachinePrecision, which is what I would expect. But I also expected a warning about the precision mis-match, like the one you get with
FindRoot. The solutions are exactly the same as those obtained without
WorkingPrecision -> 100. It seems the option is ignored, but without a warning.
I do not see anything wrong with the results. Generally in Mathematica to get a high-precision output, you have to tell it explicitly that the input is high-precision. I think this is a good approach, since if the error in the input is not as small as the precision claims, misleading results occur.
One might wonder which of the two results above is more accurate. That depends on which of the codes
Rationalize[f1[x] == f2[x], 0]
SetPrecision[f1[x] == f2[x], 100]
have coefficients that are closer to the true values of the coefficients. There is no way to know that from what is given in the question. If the digits of each parameter reflect the precision to which they are known, then the machine precision results are already more precise than the input and an error analysis would be appropriate.