How to use Mathematica in my Problem as I would use it in Matlab

Question

I would like to calculate the convolution of a function with itself for $n$ times. Let $n$ be $10$. For $n=1$, I wont make any calculation and $f0$ is my function, when $n=2$, I will convolve $f0$ with itself and get $g=f0*f0$ when $n=3$ I will get $gg=g*f0$ etc..

As it can be seen the total number of convolutions that I need to make is $9$. At each iteration I will compute some other functions. In Matlab I do it as follows

h=f0;
for n=1:9
some computations related to "f0"
h=conv(h,f0);
some computations related to "h"
end

I would like to do the same in Mathematica and below is the code that I typed:

f0[x_] := PDF[NormalDistribution[-2, 2], x]
f1[x_] := PDF[NormalDistribution[2, 2], x]


q[1, B_, A_, f_] := f

q[n_ /; n > 1, B_, A_, f_] :=Block[{\[Omega]},q[n, B, A, f] =Function[x,Evaluate[Integrate[q[n - 1, B, A, f][\[Omega]] f[x - \[Omega]],{\[Omega], B, A}]]]]

p[n_, B_, A_, f_] :=Integrate[q[n, B, A, f][x], {x, -\[Infinity], B}] +Integrate[q[n, B, A, f][x], {x, A, \[Infinity]}]

In short $q$ is calculating the convolution and $p$ some integrals over a given $q$.

If I want to calculate for example

Sum[p[n, B, A, f0]*(k3*n + k2*p[n, B, A, f1]), {n, 1, 10}]

where $k2$ and $k3$ are some constants. Then Mathematica will make 108 convolutions and the following calculations. Instead I can make the same thing in Matlab with 9 convolutions and the following calculations.

As a result of this I decided to change my code parts as follows

q0[1] := PDF[NormalDistribution[-2, 2], x]
For[n = 2, n < 10, n++,q0[n][B_, A_] :=Convolve[q0[n - 1][B, A], f0(UnitStep[x - B] - UnitStep[x - A]), x,y]]

this code shows what I want to do but It is incorrect. I want to calculate q0[1][x], q0[2][x], q0[3][x],.. which are simply the one time, 2 times, and 3 times convolutions of $f0$ with itself (in the range from $B$ to $A$ ).

$$q_n(x)=\int_{B}^A q_{n-1} (\omega)f(x-\omega)\mathrm{d} \omega,\quad q_1=f,\quad n\geq 1.$$

I would like to do some other operations for example

Sum[p[n, B, A, f0]*(k3*n + k2*p[n, B, A, f1]), {n, 1, 10}]

Thank you very much. Please let me know if something is unclear in the same for loop. How I can manage it?

Answer 1

This appears to be a follow up of this question.

As suggested altering using list will provide the intermediate results. Using:

qlist[f_, g_, x_, b_, a_, n_] := 
 Nest[Convolve[#, f (UnitStep[x - b] - UnitStep[x - a]), x, y] &, g, 
  n] 

To achieve you goal and assuming k3,k2 are or will be defined your goal can be achieved:

list1=qlist[f0[x],f0[x],x,b,a,9];
list2=qlist[f1[x],f1[x],x,b,a,9];
Total@MapThread[#1(k3 #3+k2 #2)&,{list1,list2,Range[10]}]

18 convolutions (Matlab wil have to do 9 for each initial each function). NestList will produce 10 elements (starting condition and 9 convolutions).

I am certain there may be more efficient methods but I believe this is consistent with your aims.

Answer 2

Does this do what you want:

cf[f_, n_] := NestList[Convolve[#1, f, x, y] /.y -> x &, f, n]

Usage:

f0[x_] := PDF[NormalDistribution[-2, 2], x]

Then:

cf[f0[x], 3]

Gives:

 {E^(-(1/8) (2 + x)^2)/(2 Sqrt[2 \[Pi]]), E^(-(1/16) (4 + x)^2)/(
 4 Sqrt[\[Pi]]), E^(-(1/24) (6 + x)^2)/(2 Sqrt[6 \[Pi]]), E^(-(1/32) (8 + x)^2)/(4 Sqrt[2 \[Pi]])}

If you want cf[f0[x],1] + cf[f0[x],2] + ... + cf[f0[x],10] i.e. the sum of each elememt, then just do:

Total[cf[f0[x], 10]]