Exponent of $GL(n,q)$.

Let us write $G=GL(n,q)$ for your group.

By the theorem that states that the so called Block Jordan Canonical Form of matrices over finite fields exist, we know that every matrix in $G$ is conjugate to one which is a block matrix with blocks of the form $$\begin{pmatrix}C_f\\I&C_f\\&I&C_f\\&&\ddots&\ddots\end{pmatrix} \tag{$\star$}$$ with $f$ an irreducible polynomial over the field $\mathbb F_q$ of some degree $d$, $C_d$ its companion matrix and $I$ the identity matrix of size $d$; see, for example, this notes.

It follows that the exponent of $G$ is the LCM of the orders of all matrices of this form with total size at most $n$.

So let $A$ be a matrix of the form $(\star)$, with $r\times r$ blocks, $f\in\mathbb F_q[X]$ an irreducible polynomial of degree $d$, $C_f$ its companion matrix and $I=I_d$ the $d\times d$ identity matrix. If we let $D$ be the «diagonal part» and $N$ be the «lower triangular part» (so that $D$ is a block diagonal matrix with $r$ blocks equal to $C_f$ along the diagonal and $N$ is a block matrix with non-zero blocks equal to $I_d$ exactly where these appear in $A$), we have $A=D+N$, $DN=DN$ and $N^r=0$. It follows immediately from this and the binomial formula that $$A^s=D^s+\binom{s}{1}D^{s-1}N+\binom{s}{2}D^{s-2}N^2+\cdots+\binom{s}{r-1}D^{s-r+1}N^{r-1}$$ for all $s\geq0$. It is easy to see what each term in this sum looks like, and using that we can see that $A^s=I_{dr}$ if and only if $D^s=I_{dr}$ and $\binom{s}{i}D^{s-i}=0$ for all $i\in\{1,\dots,r-1\}$. This happens iff $$\text{$C^s=I_d$ and $\binom{s}{i}=0$ for all $i\in\{1,\dots,r-1\}$.}$$ Of course, the order of $A$ is the least positive number $s$ which satisfies this condition. Plainly, it depends only on $f$ and on $r$, so we can call it $o(f,r)$. Moreover, it is clear that $o(f,r)\mid o(f,r+1)$, so the exponent of $G$ is $$LCM\{o(f,\lfloor n/\deg(f)\rfloor:f\in\mathbb F_q[X],\text{$f$ irreducible}, \deg(f)\leq n\}$$

If $n=2$, it is not very difficult to compute this: I get $LCM(p(p^{2n}-1),p(p^n-1))$, which happily equals with the formula I gave in a comment to the question.

One can in fact give a quite explicit formula for the exponent of $GL(n,\Bbb F_q)$, and a fairly explicit description of the set of orders of elements that can occur.

First, the order of $A\in GL(n,\Bbb F_q)$ only depends on the minimal polynomial $f$ of $A$, and is equal to the minimal exponent $e>0$ such that $X^e-1$ is divisible by $f$. Since $A$ is invertible the constant term of $f$ is nonzero, but otherwise $f$ can be essentially made to be any polynomial of degree at most $n$ in $\Bbb F_q[X]$ using companion matrices ("essentially" because extending a matrix with an identity matrix block may change its minimal polynomial: it becomes the $\def\lcm{\operatorname{lcm}}\lcm$ of $X-1$ and the old minimal polynomial; this does not however affect the associated value $e$). So we get that the set of orders of elements $A$ is equal to the set of exponents $e$ associated to polynomials of degree at most $n$ in $\Bbb F_q[X]$, and the exponent of $GL(n,\Bbb F_q)$ is the least $e>0$ for which $X^e-1$ is divisible by all such polynomials.

However not all such polynomials need to be considered. If we decompose a polynomial into pairwise relatively prime factors, then the associated $e$ will be the least common multiple of those numbers for the factors, by the Chinese remainder theorem for $K[X]$. Also if $f$ is irreducible of degree $d$ then the associated exponent $e_f$ will be a divisor of $q^d-1$ (because $\Bbb F_q[X]/(f)$ is a finite field of order $q^d$), and $e_f=q^d-1$ will occur for certain such polynomials, called primitive polynomials (because the multiplicative group of that finite field is cyclic). We must also consider the exponents associated to powers $f^i$ of irreducible polynomials (as long as their degree is at most $n$); here too the interesting case will be when $f$ is a primitive polynomial.

Since $X^{q^d-1}-1$ is the product of all monic irreducible polynomials with degree dividing $d$, except $X$, it is divisible by $f$ but not by $f^2$; it then follows that $q^{e_f}-1$, which is divisible by $f$ by assumption, is not divisible by $f^2$ either. Say $X^{e_f}=1+af$ with $a$ not divisible by $f$. Then by the binomial formula $X^{e_fk}\not\equiv 1\pmod{f^2}$ for $0<k<p$, where $p$ is the characteristic of $\Bbb F_q$, but on the other hand $X^{e_fp}=1+a^pf^p$ by the Frobenius endomorphism of $\Bbb F_q[X]$. Similarly $X^{e_fpk}\not\equiv 1\pmod{f^{p+1}}$ for $0<k<p$, but $X^{e_fp^2}=1+a^{p^2}f^{p^2}$. Continuing like this, it follows that the smallest exponent $e>0$ with $X^e\equiv 1\pmod{f^i}$ is given by $$ e=e_fp^{\lceil\log_p i\rceil} $$ which is maximal (in the multiplicative sense) among irreducible polynomials of degree $d$ when $f$ is primitive, and then gives $$ e=(q^d-1)p^{\lceil\log_p i\rceil} $$

The problem of finding all orders of elements of $GL(n,\Bbb F_q)$ becomes that of finding the numbers expressible as the least common multiple of a collection of such numbers $e$, each for some $(d,i)$, under the constraint that the sum of the values $di=\deg(f^i)$ be at most $n$ (of course all divisors of such numbers also occur as order of an elements of $GL(n,\Bbb F_q)$). I don't know of any really transparent way of describing this set of orders. In any case, as the link in the question indicates, the fact that $\Bbb F_q[A]$ has only $q^d-1$ nonzero elements where $d$ is the degree of the minimal polynomial of $A$ shows that the maximum (in the additive sense) possible order $q^n-1$ is attained for a single pair $(d,i)=(n,1)$.

However the exponent of $GL(n,\Bbb F_q)$ can be explicitly given. The multiplicity of $p$ in it is $\lceil\log_p n\rceil$, obtained for $(d,i)=(1,n)$. The factor relatively prime to $p$ is $\lcm(q-1,q^2-1,\ldots,q^n-1)$. It is known (see for instance this question) that $\gcd(q^k-1,q^l-1)=q^{\gcd(k,l)}-1$, so the only common factors to consider here are the ones that occur in the list themselves, and indeed given the preceding numbers, the number $q^k-1$ contributes as fresh factor the value $\Phi_k(q)$, where $\Phi_k$ is the $k$-th cyclotomic polynomial. Thus finally, the exponent of $GL(n,\Bbb F_q)$ is equal to $$ p^{\lceil\log_p n\rceil} \prod_{k=1}^n\Phi_k(q) \qquad\text{where }p\text{ is the characteristic of }\Bbb F_q. $$

As a sanity check, I'll plug in $n=2$ to give the exponent $p(q-1)(q+1)=p(q^2-1)$. For $n=3$ one gets exponent $2^2(q^2-1)(q^2+q+1)$ in characteristic $2$ (indeed the exponent of $GL(3,\Bbb F_2)$, which has order $168$, is $4\times3\times7=84$) and $p(q^2-1)(q^2+q+1)$ in any other characteristic $p$.