Lets suppose I have $\mathcal{\hat{H}}\psi=(\mathcal{\hat{H}}_0+\lambda \hat{V})\psi=0$ where $\hat{V}$ is some -non linear - perturbation to $\mathcal{\hat{H}}_0$ controlled by the parameter $\lambda$ and $\mathcal{\hat{H}}_0 $ is a linear operator. Moreover, $\mathcal{\hat{H}}_0 $ is the heat operator:
\begin{equation*} \mathcal{\hat{H}}_0 = \frac{\partial}{\partial t}-\frac{\partial^2}{\partial x^2 } \end{equation*} and the perturbation is: \begin{equation} \hat{V}\psi=\left(\frac{\partial^2\psi }{\partial x^2}\right)^2 \end{equation}
I would like to study this problem with the following initial-boundary conditions: \begin{equation*} \displaystyle \left\{\begin{array}{l} \psi(x,t) \rightarrow e^x\quad\textrm{if}\quad x\rightarrow \infty \\ \\ \displaystyle \psi(x,t) \rightarrow 0\quad\textrm{if}\quad x\rightarrow 0 \\ \\ \displaystyle \psi(x,t)= e^{max(x,0)}-1 \quad\textrm{if}\quad t= 0 \end{array}\right. \end{equation*}
If $\psi(x,t)$ is an analytical function of $\lambda$ then it can be written as: \begin{equation} \psi=\psi_0+\lambda\psi_1+...+\lambda^n\psi_n\,. \label{psi-lambda-expansion} \end{equation}
This summation will be studied up to $n=N$. Using the above eq. we can write: \begin{equation*} \mathcal{\hat{H}}\psi=\mathcal{\hat{H}}_0\psi_0 +\lambda(\mathcal{\hat{H}}_0\psi_1+\hat{V}\psi_0)+O(\lambda^2)\,. \end{equation*} From above it follows that $\mathcal{\hat{H}}\psi=0$, $\forall\lambda$ if and only if every coefficient in the polynomial expansion vanishes. Note that now the term $\hat{V}\psi_0$ could be analytically known if I define proper initial-boundary conditions to $\mathcal{\hat{H}_0}$. Indeed, I can easily solve the equation $\mathcal{\hat{H}_0}\psi_0=0$.
My problem is how to introduce consistent initial-boundary conditions to $\psi_i$. In my first try I though that maybe if I go with: \begin{equation*} \displaystyle \left\{\begin{array}{l} \psi_i(x,t) \rightarrow \frac{1}{(N+1)\lambda^i}e^x\quad\textrm{if}\quad x\rightarrow \infty \\ \\ \displaystyle \psi_i(x,t) \rightarrow 0\quad\textrm{if}\quad x\rightarrow 0 \\ \\ \displaystyle \psi_i(x,t)=\frac{1}{(N+1)\lambda^i} e^{max(x,0)}-1 \quad\textrm{if}\quad t= 0 \end{array}\right. \end{equation*}
the summation of the boundary conditions for each $\psi_i$ will sum up to the required conditions to $\psi$ but this will not work because then the solution for each $\psi_i$ will depend upon $1/\lambda^i$ and $\psi$ will not have a Taylor expansion on $\lambda$. How can I define consistent boundary-initial conditions effectively solve this problem at least at order $\lambda$?.