with the matter bounce scenario with the background (\[background\]), this expression is equal to $v_k^{(0)}= -\frac{B_1(k)}{4\sqrt{3}}\rho_c\tau^2+\frac{4B_2(k)}{\sqrt{3}\rho_c}\frac{1}{\tau}$, which one has to match with (\[xxx\]) so as to obtain $$\begin{aligned}
B_1(k)=\sqrt{\frac{8}{3}}\frac{ k^{3/2}}{\rho_c} \quad \mbox{and} \quad B_2(k)=-i\sqrt{\frac{3}{8}}\frac{\rho_c}{2 k^{3/2}}.\end{aligned}$$ When one has the expression of $v_k^{(0)}$ one can use its definition to write $\Phi_k^{(0)}$ as a function of $\delta\sigma_k^{(0)}$. Then, inserting this expression into Eq. $i-0$, one obtains a solvable first-order differential equation in $\delta\sigma_k^{(0)}$. To obtain the next iteration one has to solve, using the method of variation of constants for second order differential equations, the equation $$\begin{aligned}
v_k''+ \left(k^2-\frac{z''}{z}\right)v_k=-\frac{k^2}{z}\left(z\bar{\phi}'\left( \frac{\bar{f}_{\chi\chi}\delta\sigma_k^{(0)}}{({\bar{\phi}')}^3 }  \right)'\right)'.\end{aligned}$$ Once $v_k^{(1)}$ has been obtained, one has to use the definition of $v$ and Eq. $i-0$ to find $\delta\sigma_k^{(1)}$, and so continue successively to obtain the next iteration.

For example, we take the background (\[background\]). Since from the equation (\[z3\]) one can see that the function $\delta\sigma_k^{(0)}(t)$ has two simple poles, one in the expanding and the other in the contracting phase, namely $t_+>0$ and $t_-=-t_+$, to regularize the integral one has to add an imaginary part to the poles obtaining $t_++i\epsilon$ and $t_--i\epsilon$, then after an integration by parts one will have $$\begin{aligned}
\label{zzz}
v_k^{(1)}=B_1(k)z+B_2(k)z\int_{-\infty}^{\tau}\frac{1}{z^2}d\bar\tau
-k^2z \int_{-\infty}^{\tau}\frac{\bar{\phi}'}{z}\left( \frac{\bar{f}_{\chi\chi}\delta\sigma_k^{(0)}}{({\bar{\phi}')}^3 }  