Here's a "typical" set of charts:
>Huh? S(t) starts out at a value of 1?
We can consider the numbers as fractions of the total population. Then S(0) = 1 just means that the entire population is susceptible.
>So it looks like I(t) has a maximum of about 9%.
Yes, but that'll depend upon your choices of parameters a and b. In the above example I chose A = 0.5 and B = 0.3.
>I assume that these differential equations are hard to solve.
Dividing [2] by [1] gives:
[4] dI/dS =  1 + k/S ... where k = b/a.
The solution is:
[5] I(t) = I(0)  (S(t)  S(0)) + k log(S(t)/S(0))
If we assume that S(0) = 1 (meaning that 100% of the population is Susceptible), then we have:
[A]
I(t) = I(0)  (S(t)  1) + k log(S(t)) ... where k = b/a < 1

That'll give a graph like this (for a = 0.5, b = 0.3 and I(0) = 0.0000001%):
It shows S(t) decreasing from 100% to 32% and I(t) reaching a maximum of 9.3% then decreasing to 0.
>I(0) is pretty small don't you think?
If there are 6 billion people on the planet and 6 are currently infected, then that'd give I(0) = 0.0000001%.
>Then you're saying that, eventually, 9.3% of the world population will become infected?
Did I mention that it'll depend upon your choices of parameters a and b?
Further, I(t) gives the number of people "currently" infected.
Some may have recovered or died ... and gone to the Rpopulation.
More importantly, it'll depend upon how well the SIR mathematical description matches the actual evolution of the disease.
  
From [A], if we put I(0) ≈ 0 and S(t) = 1  x with x ≈ 0 (so we're talking about the initial spreading of the disease), then we get:
[6] I(t) ≈ x + k log(1x) or I(t) ≈ (1k)x +k x^{2}/2
For our example above, k = b/a = 0.3/0.5 = 0.6.
The approximation to the differential equation [6] compares favourably with the actual solution
[A]:
>For the start of the pandemic, right?
Right ... and we might use that to gauge the values of a and b which give the "best" fit to the actual data.
Of course, we just have a theoretical expression relating I to R ... not I to t.
To get I(t) vs t we need to solve [2], namely dI/dt = a S(t) I(t)  b I(t), using [A].
We could also try to solve [1] using [A]. That'd give:
[7] dS/dt =  a S(t) { I(0)  (S(t)  1) + k log(S(t)) } ... where k = b/a.
Alas, I can't do that, so I just did a numerical integration of [1] and [2] (using a RunkeKutta routine).
  
Here's a spreadsheet to play with. Click on the picture to download the spreadsheet.
Notes:
 The parameter a is a time scale. If the time variable t is replaced by θ = at, then the differential equations becomes:
[1A] dS/dθ =  S(t) I(t)
[2A] dI/dθ = S(t) I(t)  k I(t)
[3A] dR/dθ = k I(t)
 Increasing a then makes the disease evolution speed up. Maxima are reached soooner. The disease terminates sooner.
 Since the maximum of I(t) occurs when S(t) has fallen to the value k = b/a, changing k then varies when that occurs.
 I can't find any "good" parameter values !! (YOU try.)
Maybe I'll play with the SEIR model ...