Computation of numerical integral involving convolution

Question

I have to solve the following convolution related numerical integration problem in R or perhaps computer algebra system like Maxima.
Integral[({k(y)-l(y)}^2)dy]
where
k(.) is the pdf of a standard normal distribution
l(y)=integral[k(z)*k(z+y)dz] (standard convolution)
z and y are scalars
The domain of y is -inf to +inf.
The integral in function l(.) is an indefinite integral. Do I need to add any additional assumption on z to obtain this?
Thank you.

Answer 1

R does not do symbolic integration, just numerical integration. There is the Ryacas package which intefaces with Yacas, a symbolic math program that may help.

See the distr package for possible help with the convolution parts (it will do the convolutions, I just don't know if the result will be integrable symbolicly).

You can numerically integrate the convolutions from distr using the integrate function, but all the parameters need to be specified as numbers not variables.

Answer 2

Here is a symbolic solution from Mathematica: