Fourier method
We can attempt to find a periodic solution by decomposing the function $p(x)$ according to its Fourier series; solving the PDE for each term in the Fourier series independently; and then summing these solutions. The result isn't pretty, but it does seem to work.
a = v;
b = 1;
g = -Subscript[\[Lambda], i];
c[n_] = Subscript[\[Beta], i] Subscript[\[Psi], 0] FourierCoefficient[Abs[Sin[\[Pi] x]], x, n, FourierParameters -> {1, -2 \[Pi]}]
Note that FourierCoefficients expects the function to be defined over a region symmetric about the origin. This necessitates the use of $|\sin(\pi x)|$ rather than just $\sin(\pi x)$, so that the function is "correct" over the range $[-\frac12,\frac12]$.
ic = w[x, 0] == 0;
fourierpde[n_] = a*D[w[x, t], x] + b*D[w[x, t], t] == g*w[x, t] + c[n] Exp[2 \[Pi] I n x]
fouriersoln[n_] = FullSimplify[DSolve[{fourierpde[n], ic}, w[x, t], {x, t}]]
Well, that seems nice enough, at least. Now we just need to add these up.
analyticsoln[x_, t_] = Sum[w[x, t] /. fouriersoln[n], {n, -\[Infinity], \[Infinity]}];
Oy. But remarkably, it does give us a solution that agrees well with the numerical solution:
(Left to right: analytic Fourier series solution, numerical solution, difference between solutions. Note the difference in scale for the third graph.)
A few notes on this method:
When I plotted the analytic solution, I did have to wrap analyticsolution[x,t] in Re[] to get a better-quality plot. Without this, the plot has some small "gaps" that I believe are due to rounding error leaving analyticsolution[x,t] with a non-negligible imaginary part.
These small imaginary parts could presumably be eliminated by taking FourierCosCoefficient and FourierSinCoefficient separately, solving the PDE for both $p(x) = \sin(2 \pi n x)$ and $\cos(2 \pi n x)$, and then summing. However, I would not be surprised if under this method, Mathematica takes much longer to sum the series (or is unable to find a closed form for the solution at all.)
I tried applying FullSimplify to analyticsolution[x,t], but it did not return a simplified result in any reasonable amount of time.
An approximate solution, which would be easier for Mathematica to plot and manipulate, could be obtained by truncating the final resummation (i.e. drop all terms with $|n|$ above some threshold. In addition, if you change $p(x)$ and Mathematica is unable to resum the resulting infinite series, then an approximate series solution might be the best you can do.
If you wish to solve this equation for inhomogeneous initial conditions $w(x,0) = f(x) \neq 0$, I believe that you could to so as follows: Solve the PDE for $p(x) = 0$ and $w(x,0) = f(x)$; and then add the resulting inhomogeneous source-free solution to the solution obtained via the method above. Note that you will need to periodically extend $f(x)$ over the entire real line in the same way to make this work. You could even, in principle, decompose it in a Fourier series and obtain the solution for inhomogeneous ICs in that way.
Original method (flawed)
My original answer follows. However, the resulting solutions were discontinuous across the characteristics passing through the points $x =$ integer.
Mathematica can be coaxed into providing a solution by extending the source function $p(x)$ to a periodic version that covers the entire real line (with $p(x) = p(x-1)$ for all $x$), and then solving the PDE over the entire real line.
Clear["Global`*"]
a = 1;
b = v;
g = -Subscript[\[Lambda], i];
p = Subscript[\[Beta], i] Subscript[\[Psi], 0] Sin[Pi x] ( 2 HeavisideTheta[Sin[Pi x]] - 1)
ic = w[x, 0] == 0;
bc = w[0, t] == w[1, t];
pde = a*D[w[x, t], x] + b*D[w[x, t], t] == g*w[x, t] + p
DSolve[{pde, ic}, w[x, t], {x, t}]
FullSimplify[%]
The function Sin[Pi x] ( 2 HeavisideTheta[Sin[Pi x]] - 1) gives a rectified sine wave, which has the property that $p(x) = p(x-1)$ for all $x$. A similar type of result (not identical in form, but presumably functionally equivalent) can be found by using
p = Subscript[\[Beta], i] Subscript[\[Psi], 0]*Piecewise[{{Sin[Pi x], Sin[Pi x] >=0}, {- Sin[Pi x], Sin[Pi x] < 0}}]
instead.
However, this solution appears to contain unnatural "shocks" (discontinuities) along the characteristics passing through the "matching point" $x = 0/1$.
Looking at the difference between this solution and the numerical solution (calculated above), the difference appears to be constant in the regions between these characteristics. This suggests that there might be some way to salvage this method by more carefully defining the functions (and/or that there might be a subtle error in my method, or in Mathematica).
