@Muhammad Usman Your MatrixEquation1 is a 20x20 matrix of integers and involves an unknown parameter lambda.  You want to determine lambda so that the determinant of the matrix is zero.  So all you need to do is:

Z := Determinant(evalm(MatrixEquation1)):
fsolve(Z);

That produces the 8 roots that you have already obtained.

For some reason that I don't understand, you pollute your equations with several unnecessary evalfs. That introduces numerical truncation errors that lead to the issue that you have observed.  Just don't do evalf and you will be okay.  You don't even need to change Digits.  The default Digits:=10 works just fine.

@apronya Calculate the real and imaginary parts of the solution separately.

Download draft-2.mw

@apronya There are two distinct aspects in what I have suggested: (a) we solve the differential equation for r(q) instead of q(r), and (b) we solve it symbolically rather than numerically.  If your more complicated differential equation does not admit a symbolic solution, then (b) is out but (a) may still apply.  To proceed further, we need to see what your differential equation looks like.

@raj2018 See if this answers your questions.

Download mw.mw

@jalal The animations, as you have envisioned them, do convery some information in the cases of the sin and tan functions, but they are not practical in the cases of the cos and cot functions.  That's because sin and tan measure vertical distances in the trigonometric circle, but cos and cot correspond to horizontal distances and cannot be animated in the same way.

What your students should learn (and remember forever!) is where sin, cos, tan, and cot appear relative to the trigonometric circle.  They already know about the sin and cos.  Here is what tan and cot look like:

I don't see a pedagogical value to the aninated versions.  They may look cute, but they tend to obfuscate rather than illuminate.  I would avoid them.

Download tan-and-cot-diagrams.mw

I don't see an option for that.  I have often wished that there was one.

@acer Thanks for the further elaboration.  It looks like InertForm:-Display works better than Typesetting:-Typesett in the following example.

@acer Thanks for pointing to the Typeset and EV exports of Typesetting (both undocumented).  These can come handy in some situations and it's good to have them in one's toolbox.

@acer You wrote "I prefer the technique of using, Typsetting:-Typeset".  I am unable to find a help page on Typsetting:-Typeset.  Can you tell me where it is documented?

I see that you have a system of three first order differential equations in three unknowns.  The system involves eight unspecified coefficients aS, mu1, mu4, etc.  No initial conditions are specified.

You specify data which is a collection of four lists of four numbers each.  It's not clear how data is meant to be related to the system of differential equations.

I cannot make sense of the rest of the worksheet.  To get any help, you should explain in words, as I have done above, what it is that you want to achieve.

@Dkunb The next step still depends on what you want to do with rho_1 and rho_2.  Suppose you have calculated rho_1 and rho_2.  How are you going to use them?

By the way, I don't know the meanings of the superscripts T, *, and †.  I suppose T stands for transpose.  I don't know what the other two mean.

@Dkunb When the size of a result exceeds a set limit, Maple calculates but does not display it because displaying such large expressions would be generally useless to anyone.  It is possible to raise or completely eliminate the limiting ceiling through Maple's configuration options, but I don't see a point in doing that under the circumstances.

Based on what we see in the toy example, we expect your dsol to be a list of the form
     [ z[0,1(t)=...,  z[0,2](t)=...,  z[1,1](t)=..., z[1,2](t)=... ]
To see what your z[0,1(t) looks like, do
subs(dsol, z[0,1](t));

You will see that z[0,1(t) is a HUGE expression involving things like RootOf(_Z^4, ...).  That's because you are solving a linear system of four first order differential equations.  The characteristic polynomial of such a system is a 4th degree polynomial whose roots do not have nicely presentable forms.  (It was just lucky that things worked out neatly for the toy problem.) Things get worse as you increase the size of n.   The prudent thing to do then is to solve the system numerically.  Maple can do that, but how you would organize the code depends on what you want to do with the result.

So, to go any further, let's ask the question:  Suppose you have obtained the solution Z(t).  What are you going to do with it?

@acer Thanks for checking.  Mine was done in Maple 2023. 

One way of making the code work in Maple 2021 is to load the Physics package as you have noted.  A better solution is to change the definition of the IC to

IC := eval(Z, t=0) =
	Array(0..n-1, 0..n-1, (i,j) -> `if`(i=j,1,0));

which will work in Maple 2021 without loading the Physics package.

@Dkunb What kind of object is Z(t) in dz?  Is it a scalar, a vector, a matrix?  The expression you have for dz makes no sense under any of those interpreations.  Perhaps you have something else in mind?

I am fond of following Polya's advice which in effect says if you don't know how to solve a problem, do a simpler one first.  In your case, can you just drop all those messy preliminaries?  Just make up a simple example, such as:

I am not saying that this is a suitable example since I have no idea what you have in mind for your dz.  Change that example to make it similar, in essence, to your dz.  Add extra complications only after you get that cartoon model to work.

How do you expect to solve a differential equation numerically when some of the parameters have no numerical values?  Do:

indets(eq1, name);

to see what eq1 depends on.  Maple responds with:
        {M, N, m, mu, n, nu, t, KroneckerDelta[m, n], KroneckerDelta[mu, nu]}

Furthermore, it looks like the purpose of your zint array is to supply initial conditions to the differential equation.  But the initial conditions should like like z(m,nu,0)=something.  Yours doesn't do that.

I suggest forgetting about the array for the moment. That's the least of your problems.  Just solve one of your differential equations.  Once you figure out how that works, only then you will worry about the array.