I am a bit late to the party, but I have had similar questions to yours in the past. I will summarise below what I know, which has been able to "quell" my dissatisfactions about integrable systems for the time being. Maybe it's just a placebo though...
Definition (at least one that I like)
A system with $n$ degrees of freedom, and $n$ constants of motion,
such that the Poisson bracket of any pair of constants of motion
vanishes, is known as a completely integrable system. Such a
collection of constants of motion are said to be in involution with
each other.
Completely integrable is opposed to just "partially" integrable, meaning that you cannot get a full analytical solution, and to superintegrable which are systems with $>n$ constants of motion. For example orbital motion is superintegrable (when neglecting inter-planetary interactions, or 'perturbations') because you not only have energy and angular momentum, but also the Runge-Lenz vector as conserved quantities..
For $F$ to be a constant of motion, you need $\{F, H\} = 0$, while for two constants of motion to be in involution you need $\{F_1, F_2\} = 0$. This essentially means that any two constants of motion are "compatible" with one another and can be employed at the same time. An example where this is not true is, jumping to quantum mechanics and hence replacing Poisson brackets with commutation relations, angular momentum. While each component of the angular momentum operator commutes with the Hamiltonian $[ \hat L_i , \hat H] = 0$, they are not in involution with each other as $[ \hat L_i , \hat L_j] = \epsilon_{ijk} \hat L_k$. Hence you cannot "use" all three of them at the same time, but usually opt for the conventional $\hat{\mathbf{L}}^2$ and $\hat L_z$.
Practical definition
Integrable systems are nonlinear differential equations which ‘in
principle’ can be solved analytically. This means that the solution
can be reduced to a finite number of algebraic operations and
integrations.
This definition is taken from here which also makes an interesting distinction between ODEs and PDEs, arguing that integrability is not fully defined for the latter.
An example would be a $1D$ system with Hamiltonian $H(p,q) = p^2/2 + V(q)$, obeying the usual Hamilton's equations $\dot q = p, \dot p = -\mathrm{d}V/\mathrm{d}q$. Using the conserved quantity (energy) $E = p^2/2 + V(q)$ coming out of $\{H,H\} =0$, you can write $p = \pm \sqrt{2(E-V(q))}$ and hence:
$$ t = \pm \frac{\mathrm{d}q}{\sqrt{2(E-V(q))}},$$
which then you can invert to find $q(t)$. I only have one integral relating $q$ and $t$. In an $n$-dimensional system $q_1, q_2,... q_i$, if one can still break down the problem to $n$ integrals each involving a specific $q_i$ and $t$ only, then one would have fully integrated the system without needing to solve coupled differential equations simultaneously (and very probably numerically).
And I think this qualifies as an answer to why these systems are associated with being nice.
Closed form solution
Apart from having an 'analytical' solution, meaning you can write $N$ equations each with only one $q_i$ as outlined above (as opposed to $N$ coupled equations), the term 'closed form solution' may have a geometrical interpretation.
Given a system with variables $q_1, q_2, ... q_i, ... q_n$ and $p_1, p_2, ... p_i, ... p_n$, the parameter space is $\mathbb{R}^{2n}$.
The actual trajectory that solutions to the equations of motion $(q_s,p_s)$ will follow live in the phase space $U \subseteq \mathbb{R}^{2n}$. For instance a 2D harmonic oscillator may give you closed Lissajous figures when $\omega_1/\omega_2 \in \mathbb{Q}$ (hence $U \subset\mathbb{R}^{2n}$), but these curves will fill the whole space densely for $\omega_1/\omega_2 \not \in \mathbb{Q}$ (hence $U = \mathbb{R}^{2n}$).
A more topological explanation for this involves lines on invariant tori and Hopf fibration, which I don't know enough about.
Infinite conserved quantities
Take $N$ free particles in a closed container that are not interacting among each other, but only bounce back from the container walls. Each particle conserves momentum (and energy, though this is trivial in the absence of a potential). Make $N\rightarrow \infty$, and you have infinitely many conserved quantities.
The system might look like a mess since you have a zillion particles going all over the place, but each particle is doing its own thing, following an equation of motion that is independent of what the other particles are doing.
Given the initial conditions, + the conservation of each particle's momentum, the system is completely integrable.
Thermalisation
Continuing the $N$ particle example from above. If the system is integrable, then it cannot thermalise.
Thermalisation means that, eventually, the velocity/momentum/energy distribution tends to a Maxwell-Boltzmann profile (or whatever for quantum degenerate gases). This can only be true if particles are allowed to interact (at least a bit) so as to "redistribute" the momentum. Total energy (and total momentum) will still be conserved, but by allowing inter-particle interactions you are not enforcing the $N$ ($\rightarrow \infty$) conservations of each momentum.
You might have heard of a phenomenon called Many-Body Localisation (MBL), where a (quantum) system in the presence of weak disorder seems to remain localised despite the presence of interactions, thereby not reaching thermal equilibrium. This is connected ot the Eigenstate Thermalisation Hypothesis (mentioned in the other answer) which seems to tbe one of the few criteria to classify this MBL phase, as you can't use any symmetries and stuff since it's not an equilibrium state.
Well, a mathematical model to explain how this might be possible (some maths for this in presented in this review) assumes a set of localised conserved charges that are constants of motion for the system, sometimes referred to as $\ell$-bits, making the system (“locally”) integrable.
Relation to chaos
This question here is interesting about this, but to be honest I have yet to understand the precise definition of chaos. Sometimes it seems to be a deterministic system that heavily depends on initial conditions, while some other times it means that the 'approximate present cannot predict the future'.
