I can't answer your second question, but I have a pretty fast recurrence for computing a single Littlewood-Richardson coefficient for the complete flag variety. It also works for equivariant coefficients, but unlike other equivariant recurrences I've seen you can just set the equivariant part to 0 and get a recurrence for the ordinary coefficients. It does not restrict to a recurrence for Grassmannian coefficients only, so we have to expand our attention to the entire symmetric group. It is not published anywhere so I will prove it here.

The main result of my dissertation is a Leibniz formula for divided difference operators. If $p$ and $q$ are two polynomials in variables $x_1,x_2,\ldots,x_n$ and $w$ is a permutation, then
$$\partial_w(pq)=\sum_{u,v}{c_{u,v}^w(x_1,\ldots,x_n)\partial_u(p)\partial_v(q)}$$
where $c_{u,v}^w(x_1,\ldots,x_n)$ is the equivariant Littlewood-Richardson coefficient, which is the coefficient appearing in the expansion
$$S_u(x;y)S_v(x;y)=\sum_{w}{c_{u,v}^w(y_1,\ldots,y_n)S_w(x;y)}$$
of a product of two double Schubert polynomials. When $\ell(u)+\ell(v)=\ell(w)$, these coefficients are the ordinary Littlewood-Richardson coefficients (by which I mean for the complete flag variety).

For an integer $k$ I denote by $s_k$ the simple transposition $(k,k+1)$. Now, let $s_i$ be a right descent of $w$, so that
$$\partial_w=\partial_{ws_i}\partial_{s_i}$$
Then since
$$\partial_{s_i}(pq)=\partial_{s_i}(p)q+p\partial_{s_i}(q)-(x_i-x_{i+1})\partial_{s_i}(p)\partial_{s_i}(q)$$
we can write
$$\partial_w(pq)=\partial_{ws_i}(\partial_{s_i}(p)q+p\partial_{s_i}(q)-(x_i-x_{i+1})\partial_{s_i}(p)\partial_{s_i}(q))=\sum_{u,v}{c_{us_i,v}^{ws_i}\partial_{us_i}\partial_{s_i}(p)\partial_v(q)}+\sum_{u,v}{c_{u,vs_i}^{ws_i}\partial_u(p)\partial_{vs_i}\partial_{s_i}(q)}+\sum_{u',u,v}{c_{u',u,v}^{ws_i}\partial_{u'}(x_i-x_{i+1})\partial_{us_i}\partial_{s_i}(p)\partial_{vs_i}\partial_{s_i}(q)}$$
For the last sum,
$$\sum_{u',u,v}{c_{u',u,v}^{ws_i}\partial_{u'}(x_i-x_{i+1})\partial_{us_i}\partial_{s_i}(p)\partial_{vs_i}\partial_{s_i}(q)}=\sum_{u',v',u,v}{c_{u',v'}^{ws_i}c_{us_i,vs_i}^{v'}\partial_{u'}(x_i-x_{i+1})\partial_{us_i}\partial_{s_i}(p)\partial_{vs_i}\partial_{s_i}(q)}$$
If $v'=ws_i$, then we get a contribution of
$$\sum_{u',u,v}{c_{us_i,vs_i}^{ws_i}c_{u',ws_i}^{ws_i}\partial_{u'}(x_i-x_{i+1})\partial_{us_i}\partial_{s_i}(p)\partial_{vs_i}\partial_{s_i}(q)}=\sum_{u,v}{c_{us_i,vs_i}^{ws_i}(x_{w(i+1)}-x_{w(i)})\partial_{us_i}\partial_{s_i}(p)\partial_{vs_i}\partial_{s_i}(q)}$$
(I snuck in that last equality, it's true because acting with $ws_i$ is equivalent to summing over all $u'$ what you get by acting with $c_{u',ws_i}^{ws_i}\partial_{u'}$.) If $v'<ws_i$ (in Bruhat order) and the term is nonzero, then $\ell(u')=1$, so $u'=s_j$ for some $j$. Indeed, either $j=i-1$, $j=i$, or $j=i+1$. Thus we get
$$\sum_{s_j,v',u,v}{c_{us_i,vs_i}^{v'}c_{s_j,v'}^{ws_i}\partial_{s_j}(x_i-x_{i+1})\partial_{us_i}\partial_{s_i}(p)\partial_{vs_i}\partial_{s_i}(q)}$$
By Monk's formula, in order for this to be nonzero we must have $\ell(v')=\ell(ws_i)-1=\ell(w)-2$ and there exist $a\leq j<b$ such that
$$ws_i=v'(a,b)$$
where $(a,b)$ denotes a transposition.

Long story short, picking $i$ so that $w(i)>w(i+1)$ we have

(1) If $u(i)<u(i+1)$ and $v(i)<v(i+1)$, then
$$c_{u,v}^w=0$$

(2) If $u(i)>u(i+1)$, and $v(i)<v(i+1)$, then
$$c_{u,v}^w=c_{us_i,v}^{ws_i}$$

(3) If $u(i)<u(i+1)$, and $v(i)>v(i+1)$, then
$$c_{u,v}^w=c_{u,vs_i}^{ws_i}$$

Now for the fun part.

(4) If $u(i)>u(i+1)$ and $v(i)>v(i+1)$, then
$$c_{u,v}^w=c_{us_i,v}^{ws_i}+c_{u,vs_i}^{ws_i}+\sum_{\stackrel{a\leq i-1<b}{\ell(ws_i(a,b))=\ell(w)-2}}{c_{us_i,vs_i}^{ws_i(a,b)}}+\sum_{\stackrel{a'\leq i+1<b'}{\ell(ws_i(a',b'))=\ell(w)-2}}{c_{us_i,vs_i}^{ws_i(a',b')}}-2\sum_{\stackrel{a''\leq i<b''}{\ell(ws_i(a'',b''))=\ell(w)-2}}{c_{us_i,vs_i}^{ws_i(a'',b'')}}-(x_{w(i+1)}-x_{w(i)})c_{us_i,vs_i}^{ws_i}$$
giving a recurrence for the equivariant coefficients in terms of $x_1,\ldots,x_n$. The formatting is a bit off, so note that the first two summations have nothing in front of them (so coefficient $+1$), and the third summation has a coefficient of $-2$. If you're not interested in the positive degree coefficients, ignore the last term (setting $x=0$). This is not so much fun for a human to use but on a computer it gets you the single coefficient fast (provided you either memoize or traverse the Bruhat graph intelligently).

The base case(s) for the recurrence is
$$c_{u,1}^u=c_{1,v}^v=1$$
or, if you're really a minimalist, it's sufficient to only assume
$$c_{1,1}^1=1$$

allthe terms in the expansion of $s_\lambda s_\mu$, one method is to work with a sufficiently large but finite number of variables, express the product as a linear combination of Schur functions with unknown coefficients, write the Schur functions as bialternants (quotients of determinants), and specialize the variables to real numbers sufficiently generically to be able solve the resulting system of linear equations for the unknown coefficients. I believe that John Stembridge uses this technique for some of his SF computations. $\endgroup$