solving/approximating integral of lognormal cdf I am trying to solve the integral,
\begin{equation}
\int_{0}^b\Phi_{LN}(x, \mu,\sigma)dx
\end{equation}
where $\Phi_{LN}(x, \mu,\sigma)$ is the lognormal cdf, evaluated at $x$.  i.e., $$\int^x_{0}\frac{1}{\sqrt{2\pi\sigma}y}\exp\left({-\frac{1}{2}\left(\frac{\log y-\mu}{\sigma}\right)^2}\right)dy$$.
Being able to efficiently compute the outer integral (top equation) would be great. This would ideally involve expressing it in a form that most software could efficiently solve/approximate.
I thought I could have solved this problem by solving the integral to a standard normal cdf, and simply exchanging $x$ with $\frac{log(x)-\mu}{\sigma}$, however this does not seemed to have worked.
In case someone wants to use the integral of the standard normal cdf to answer this question it is shown following,
$$
\int\Phi_{N(0,1)}(x) = x\Phi_{N(0,1)}(x) - \int x\Phi^{'}_{N(0,1)}(x)dx
$$
and here is my failed attempt to convert it to work as the lognormal equivalent,
> # define function to integrate N(0,1) CDF
> int_LN_CDF = function(x, mu, sigma){
+   # transform to lognormal RV
+   y = (log(x) - mu)/sigma
+   # integral of standard normal cdf
+   y*pnorm(y) + 1/sqrt(2*pi)*exp(-(y^2)/2)
+ }
> # set arbitary values of mu and sigma 
> mu = log(2); sigma = 0.4
> # test function against numerical approach
> int_LN_CDF(2, mu, sigma) - int_LN_CDF(1, mu, sigma)
[1] 0.3820694
> integrate(function(y) pnorm(y, mu, sigma), lower = 1, upper = 2)
[1] 0.9491263 with absolute error < 1.1e-14

Any help or guidance would be appreciated.
Thanks
 A: Don't reinvent wheels that involve tedious numerical analysis.
NIST tells me the CDF of the lognormal distribution with $x \ge 0$ and $\sigma > 0$ is 
$$ \Phi_N \left(\frac{\ln(x)} {\sigma} \right) \hspace{.2in} x \ge 0; \sigma > 0  $$
where $\Phi_N$ is the CDF of the standard normal distribution 
$$ \Phi_N(x) = \int_{-\infty}^{x} \frac{e^{-u^{2}/2}} {\sqrt{2\pi}} \,\mathrm{d} u  \text{.}  $$
They also provide Dataplot and R codes on this page
The Gnu Scientific Library contains an implementation of the CDF of the lognormal distribution.
Boost has an implementation.  (This page also better explains the $P$ and $Q$ function references on the Gnu Scientific Library page.)
SciPy has an implementation.
A: You can use a similar approach as for the standard normal distribution used here solving/approximating integral of standard normal cdf.
I'll denote by $\Phi$ and $\varphi$ the distribution and density function of the standard normal distribution - the general case should be feasible from there.
The distribution function of the log-normal distribution is then $F(x) = \Phi(\log(x))$.
Starting with
\begin{align}
\frac{d}{dx} xF(x) = F(x) + \varphi(\log(x))
\end{align}
we need to find a function $L(x)$ such that $\frac{d}{dx} L(x) = -\varphi(\log(x))$.
Then we can finally define $I(x) = xF(x) + L(x)$ with $\frac{d}{dx} I(x) = F(x)$ which is what you are looking for.
Finding $L$ is not straight forward - so I asked Wolfram Alpha for help.
With
$L(x) = \frac{1}{2 e^{-0.5}} \mbox{erf} ((1-\log(x)) / \sqrt{2})$, where $\mbox{erf}$ is the error function
\begin{align}
\mbox{erf}(x) = \frac{2}{\sqrt{\pi}} \int_0^x e^{-t^2} dt
\end{align}
you should get your result.

Edit: to make my suggestion in the comment concerning the general case more explicit I suggest to look at
\begin{align}
L(x) = \frac{1}{2}e^{(\mu+\sigma^2/2)}\mbox{erf}\left(\left(\sigma - \frac{\log(x)-\mu}{\sigma}\right)/\sqrt{2}\right).
\end{align}
With $\frac{d}{dx}\mbox{erf}(x) = \frac{2}{\sqrt{\pi}}e^{-x²}$ we then get
\begin{align}
\frac{d}{dx} L(x) &= -\frac{1}{2}e^{(\mu+\sigma^2/2)}\frac{2}{\sqrt{\pi}}e^{-\frac{1}{2}\left(\sigma - \frac{\log(x)-\mu}{\sigma}\right)^2}\frac{1}{x\sigma\sqrt{2}}\\
&= -\frac{1}{x\sigma\sqrt{2\pi}}e^{(\mu+\sigma^2/2)}e^{-\frac{1}{2}\left(\sigma^2 + \left(\frac{\log(x)-\mu}{\sigma}\right)^2 - 2(\log(x)-\mu)\right)}\\
&= -\frac{1}{x\sigma\sqrt{2\pi}}e^{(\mu+\sigma^2/2)}e^{-\frac{1}{2}\sigma^2}e^{-\frac{1}{2}\left(\frac{\log(x)-\mu}{\sigma}\right)^2}e^{\log(x)}e^{-\mu}\\
&= -\frac{1}{\sigma\sqrt{2\pi}}e^{-\frac{1}{2}\left(\frac{\log(x)-\mu}{\sigma}\right)^2}\\
&= -\frac{1}{\sigma}\varphi\left(\frac{\log(x)-\mu}{\sigma}\right).
\end{align}
Turning back to the general problem with $F_{\mu,\sigma}(x) = \Phi\left(\frac{\log(x)-\mu}{\sigma}\right)$ we now go on defining $I(x) = xF_{\mu,\sigma}(x) + L(x)$ to get
\begin{align}
\frac{d}{dx}I(x) &= F_{\mu,\sigma}(x) + xf_{\mu,\sigma}(x) - \frac{1}{\sigma }\varphi\left(\frac{\log(x)-\mu}{\sigma}\right)\\
&= F_{\mu,\sigma}(x),
\end{align}
i.e. $I(x)$ is the anti-derivative of $F_{\mu,\sigma}(x)$.
