1 Introduction

Problem of studying the stability of functional equations has begun with a question posed by S. Ulam (see, e.g., [17]) and an answer given by D.H. Hyers [13]. Since then a number of papers investigating the so called now Hyers–Ulam stability have appeared. The results concern also various generalizations of the problem and these kind of research have their origins in the papers by T. Aoki [1], D.G. Bourgin [7], Th. Rassias [16], P. Gavruta [11].

Let X and Y be linear spaces over the same field \(\mathbb {F}\in \{\mathbb {R},\mathbb {C}\}\), \(a_1,a_2,b_1,b_2\in \mathbb {F}\setminus \{0\},\) \(C_{1},C_{2},C_{3},C_{4}\in \mathbb {F}\) and \(f: X^2\rightarrow Y\). In [10], K. Ciepliński starting with a bilinear mapping, i.e., linear in each of its arguments, considered the following functional equation

$$\begin{aligned}&f(a_1x_1+a_2x_2,b_1y_1+b_2y_2)=C_{1}f(x_1,y_1)\nonumber \\&\quad +\,C_{2}f(x_1,y_2)+C_{3}f(x_2,y_1)+C_{4}f(x_2,y_2) \end{aligned}$$
(1)

for all \(x_1,x_2,y_1,y_2 \in X, \) and investigated, among others, its Hyers–Ulam stability in Banach spaces. In fact, he proved the stability without knowing the general solution of (1) and under some additional assumptions. In [6], the authors described the form of solutions of (1). They were also studying relations between (1) and bilinear mappings.

In the present paper, firstly knowing already the form of solutions of (1) we prove its Hyers–Ulam stability, also in the cases excluded in [10]. Secondly, applying the fixed point method, we study the generalized stability of (1) for the same classes of control functions.

Particular cases of (1) are, among others, the following three functional equations:

$$\begin{aligned}&f(x_1+x_2,y_1+y_2)=f(x_1,y_1)+f(x_1,y_2)+f(x_2,y_1)+f(x_2,y_2), \\&4f\Big (\frac{x_1+x_2}{2},\frac{y_1+y_2}{2}\Big )=f(x_1,y_1)+f(x_1,y_2)+f(x_2,y_1)+f(x_2,y_2), \\&2f\Big (x_1+x_2,\frac{y_1+y_2}{2}\Big )=f(x_1,y_1)+f(x_1,y_2)+f(x_2,y_1)+f(x_2,y_2), \end{aligned}$$

that is, the equations which characterize biadditive, bi-Jensen and Cauchy-Jensen mappings, respectively. Therefore, our results generalize stability outcomes for these equations (see, e.g., [2,3,4,5, 9, 14, 15]).

Define \(C:=C_1+C_2+C_3+C_4.\)

For the convenience of the reader we recall here a result describing the solutions of (1) (see [6], and also [12], where Y is an arbitrary field of characteristic different from two).

Theorem 1

If \(f :X^2 \rightarrow Y\) satisfies (1), then there exist a biadditive function \(g:X^2\rightarrow Y,\) additive functions \(\varphi , \psi :X \rightarrow Y\) and a constant \(\delta \in Y\) such that

$$\begin{aligned} f(x,y)=g(x,y)+\varphi (x)+\psi (y)+\delta , \end{aligned}$$
(2)

and

$$\begin{aligned}&\begin{array}{l} g(a_1x, b_1y)=C_{1}g(x,y),\\ g(a_1x, b_2y)=C_{2}g(x,y),\\ g(a_2x, b_1y)=C_{3}g(x,y),\\ g(a_2x, b_2y)=C_{4}g(x,y), \end{array} \end{aligned}$$
(3)
$$\begin{aligned}&\begin{array}{l} \varphi (a_1x)=(C_{1}+C_{2})\varphi (x),\\ \varphi (a_2x)=(C_{3}+C_{4})\varphi (x), \end{array} \end{aligned}$$
(4)
$$\begin{aligned}&\begin{array}{l} \psi (b_1y)=(C_{1}+C_{3})\psi (y),\\ \psi (b_2y)=(C_{2}+C_{4})\psi (y), \end{array} \end{aligned}$$
(5)

for all \(x,y \in X\), and

$$\begin{aligned} \delta ( C- 1)=0. \end{aligned}$$
(6)

Conversely, each function f of the form (2) with g biadditive, \(\varphi , \psi \) additive, and such that conditions (3), (4), (5), (6) are satisfied, is a solution of (1).

Throughout this paper we keep the standard notation: \(\mathbb {N}, \mathbb {R}\) and \(\mathbb {C}\) stand for the sets of all positive integers, all real numbers and all complex numbers, respectively. Moreover, we denote \(\mathbb {R}_+ :=[0, \infty ),\) \(\mathbb {N}_0:=\mathbb {N}\cup \{0\}\) and we adopt the convention \(0^0=1\).

2 Hyers–Ulam Stability of (1)

We start the section with recalling two stability results: for the Cauchy equation (see [13]) and for the biadditivity equation (see, e.g., [5, 9]).

Lemma 1

Let \((H,+)\) be an abelian group and \((Y,\Vert \cdot \Vert )\) be a Banach space. Given \(\varepsilon >0\) assume that \(f:H \rightarrow Y\) satisfies

$$\begin{aligned} \Vert f(x+y)-f(x)-f(y)\Vert \le \varepsilon , \qquad x,y\in H. \end{aligned}$$

Then there exists an additive function \(F:H \rightarrow Y\) such that

$$\begin{aligned} \Vert f(x)-F(x)\Vert \le \varepsilon , \qquad x \in H. \end{aligned}$$

Moreover, F is a unique function satisfying the above condition and it is of the form \(F(x)=\lim \limits _{n \rightarrow \infty } \frac{1}{2^n}f(2^n x)\) for all \(x \in H\).

Lemma 2

Let \((H,+)\) be an abelian group and \((Y,\Vert \cdot \Vert )\) be a Banach space. Given \(\varepsilon >0\) assume that \(g:H^2 \rightarrow Y\) satisfies

$$\begin{aligned}&\Vert g(x_1 +x_2, y_1+y_2) -g(x_1,y_1) - g(x_1,y_2) \\&\quad - g(x_2,y_1) - g(x_2,y_2)\Vert \le \varepsilon ,\qquad x_1,x_2,y_1,y_2\in H. \end{aligned}$$

Then there exists an additive function \(G:H^2 \rightarrow Y\) such that

$$\begin{aligned} \Vert g(x,y)-G(x,y)\Vert \le \frac{1}{3}\varepsilon , \qquad x,y \in H. \end{aligned}$$

Moreover, G is a unique function satisfying the above condition and it is of the form \(G(x,y)=\lim \limits _{n \rightarrow \infty } \frac{1}{4^n}g(2^n x, 2^n y)\) for all \(x,y \in H\).

Now we are able to present the main result of this section.

Theorem 2

Let \((Y,\Vert \cdot \Vert )\) be a Banach space and \(\varepsilon > 0\). Assume that \( f:X^2 \rightarrow Y\) is a mapping such that

$$\begin{aligned} \Vert f(a_1x_1+a_2x_2,b_1y_1+b_2y_2)&-C_{1}f(x_1,y_1)-C_{2}f(x_1,y_2)\nonumber \\&-C_{3}f(x_2,y_1)-C_{4}f(x_2,y_2)\Vert \le \varepsilon \end{aligned}$$
(7)

for \(x_1,x_2,y_1,y_2 \in X.\) Then there exists a solution \(F:X^2 \rightarrow Y\) of (1) such that

$$\begin{aligned} \Vert f(x,y)-f(0,0)-F(x,y)\Vert \le 14 \varepsilon , \qquad x,y \in X. \end{aligned}$$
(8)

Moreover, if \(C\ne 1\) then F is a unique solution of (1) such that (8) holds.

Proof

Immediately from (7) we obtain the following inequalities

$$\begin{aligned}&\Vert f(a_1x_1+a_2x_2,0)-C_{1}f(x_1,0)-C_{2}f(x_1,0)-C_{3}f(x_2,0)-C_{4}f(x_2,0)\Vert \le \varepsilon ,\nonumber \\ \end{aligned}$$
(9)
$$\begin{aligned}&\Vert f(0,b_1y_1+b_2y_2)-C_{1}f(0,y_1)-C_{2}f(0,y_2)-C_{3}f(0,y_1)-C_{4}f(0,y_2)\Vert \le \varepsilon , \nonumber \\ \end{aligned}$$
(10)

for all \(x_1,x_2,y_1,y_2 \in X\), and

$$\begin{aligned} \Vert (1-C)f(0,0)\Vert \le \varepsilon . \end{aligned}$$
(11)

Therefore, the functions \(\varphi (x):=f(x,0)-f(0,0)\) for \(x \in X,\) and \(\psi (y):=f(0,y)-f(0,0)\) for \(y \in X,\) satisfy the conditions

$$\begin{aligned} \Vert \varphi (a_1x_1+a_2x_2)-C_{1}\varphi (x_1)-C_{2}\varphi (x_1)-C_{3}\varphi (x_2)-C_{4}\varphi (x_2)\Vert \le 2\varepsilon \end{aligned}$$
(12)

and

$$\begin{aligned} \Vert \psi (b_1y_1+b_2y_2)-C_{1}\psi (y_1)-C_{2}\psi (y_2)-C_{3}\psi (y_1)-C_{4}\psi (y_2)\Vert \le 2\varepsilon , \end{aligned}$$
(13)

respectively.

By (7) we also have

$$\begin{aligned} \begin{array}{l} \Vert f(a_1x_1,b_1y_1)-C_{1}f(x_1,y_1)-C_{2}f(x_1,0)-C_{3}f(0,y_1)-C_4f(0,0)\Vert \le \varepsilon , \\ \Vert f(a_1x_1,b_2y_2)-C_{1}f(x_1,0)-C_{2}f(x_1,y_2)-C_3f(0,0)-C_{4}f(0,y_2)\Vert \le \varepsilon , \\ \Vert f(a_2x_2,b_1y_1)-C_{1}f(0,y_1)-C_2f(0,0)-C_{3}f(x_2,y_1)-C_{4}f(x_2,0)\Vert \le \varepsilon ,\\ \Vert f(a_2x_2,b_2y_2)-C_1f(0,0)-C_{2}f(0,y_2)-C_{3}f(x_2,0)-C_{4}f(x_2,y_2)\Vert \le \varepsilon \end{array}\nonumber \\ \end{aligned}$$
(14)

and, moreover,

$$\begin{aligned} \begin{array}{l} \Vert f(a_1x_1,0)-(C_{1}+C_{2})f(x_1,0)-(C_{3}+C_{4})f(0,0)\Vert \le \varepsilon , \\ \Vert f(a_2x_2,0)-(C_{3}+C_{4})f(x_2,0)-(C_{1}+C_{2})f(0,0)\Vert \le \varepsilon , \\ \Vert f(0,b_1y_1)-(C_{1}+C_{3})f(0,y_1)-(C_{2}+C_{4})f(0,0)\Vert \le \varepsilon ,\\ \Vert f(0,b_2y_2)-(C_{2}+C_{4})f(0,y_2)-(C_{1}+C_{3})f(0,0)\Vert \le \varepsilon . \end{array} \end{aligned}$$
(15)

From (7), (11), (14) and (15) it follows that

$$\begin{aligned}&\Vert f(a_1x_1 + a_2x_2, b_1y_1 +b_2y_2) - f(a_1x_1,b_1y_1) - f(a_1x_1,b_2y_2) \\&\quad - f(a_2x_2,b_1y_1) - f(a_2x_2,b_2y_2)+f(a_1x_1,0)+ f(a_2x_2,0)\\&\quad + f(0,b_1y_1)+ f(0,b_2y_2)-f(0,0)\Vert \le 10\varepsilon , \end{aligned}$$

and, since \(a_1a_2b_1b_2\ne 0\),

$$\begin{aligned} \begin{array}{ll} &{} \Vert f(x_1 +x_2, y_1+y_2) -\, f(x_1,y_1) - f(x_1,y_2) - f(x_2,y_1) - f(x_2,y_2)\\ &{} \quad +\,f(x_1,0)+ f(x_2,0)+ f(0,y_1)+ f(0,y_2)-f(0,0)\Vert \le 10\varepsilon . \end{array} \end{aligned}$$
(16)

From (9), (11) and (15) we obtain

$$\begin{aligned} \Vert f(x_1+x_2,0)-f(x_1,0)-f(x_2,0)+f(0,0)\Vert \le 4\varepsilon , \end{aligned}$$
(17)

and by (10), (11) and (15) we have

$$\begin{aligned} \Vert f(0,y_1+y_2)-f(0,y_1)-f(0,y_2)+f(0,0)\Vert \le 4\varepsilon , \end{aligned}$$
(18)

so, for all \(x_1,x_2,y_1,y_2 \in X\),

$$\begin{aligned} \Vert \varphi (x_1+x_2)-\varphi (x_1)-\varphi (x_2)\Vert \le 4\varepsilon ~\text{ and }~ \Vert \psi (y_1+y_2)-\psi (y_1)-\psi (y_2)\Vert \le 4\varepsilon . \end{aligned}$$

On account of Lemma 1, there exist a unique additive function \(\Phi \) and a unique additive function \(\Psi \) such that

$$\begin{aligned} \Vert \varphi (x)-\Phi (x)\Vert \le 4\varepsilon , \quad \; \Vert \psi (x)-\Psi (x)\Vert \le 4 \varepsilon , \qquad x \in X, \end{aligned}$$
(19)

with

$$\begin{aligned} \Phi (x)=\lim _{n \rightarrow \infty } \frac{1}{2^n}\varphi \left( 2^n x\right) ,\quad \; \Psi (x)= \lim _{n \rightarrow \infty } \frac{1}{2^n}\psi \left( 2^n x\right) ,\qquad x \in X. \end{aligned}$$

Therefore using (12) and (13), we derive that the functions \(F_1(x,y):=\Phi (x)\) and \(F_2(x,y):=\Psi (y)\), for \(x,y \in X\), satisfy (1).

Let us define \(g :X^2 \rightarrow Y\) by

$$\begin{aligned} g(x,y):=f(x,y)-f(x,0)-f(0,y)+f(0,0), \quad x,y \in X. \end{aligned}$$
(20)

Then

$$\begin{aligned} f(x,y)=g(x,y)+f(x,0)+f(0,y)-f(0,0)=g(x,y)+\varphi (x)+\psi (y)+f(0,0) \end{aligned}$$

and by (7), (9), (10) and (11), we get

$$\begin{aligned} \begin{array}{ll} \Vert g(a_1x_1 + a_2x_2,b_1y_1+b_2y_2)\!\!\! &{} -\;C_{1}g(x_1,y_1)-C_{2}g(x_1,y_2)\\ &{} - \;C_{3}g(x_2,y_1)-C_{4}g(x_2,y_2)\Vert \le 4\varepsilon . \end{array} \end{aligned}$$
(21)

On account of (16), (17), (18) and (20), we obtain

$$\begin{aligned} \Vert g(x_1 +x_2, y_1+y_2) -g(x_1,y_1) - g(x_1,y_2) - g(x_2,y_1) - g(x_2,y_2)\Vert \le 18\varepsilon . \end{aligned}$$

By Lemma 2, there exists a unique biadditive function G such that

$$\begin{aligned} \Vert g(x,y)-G(x,y)\Vert \le 6 \varepsilon , \end{aligned}$$
(22)

and, moreover, \(G(x,y)=\lim \limits _{n \rightarrow \infty } \frac{1}{4^n} g(2^n x,2^n y)\). Using (21), we obtain that G satisfies (1).

Let us define

$$\begin{aligned} F(x,y):=G(x,y)+\Phi (x)+\Psi (y), \quad x,y \in X. \end{aligned}$$
(23)

Function F satisfies (1) and from (19) and (22) we get

$$\begin{aligned}&\Vert f(x,y)-f(0,0)-F(x,y)\Vert \le \Vert g(x,y)-G(x,y) \Vert +\Vert \varphi (x)-\Phi (x)\Vert \\&\quad +\Vert \psi (y)-\Psi (y)\Vert \le 14\varepsilon . \end{aligned}$$

For the proof of the uniqueness, assume that \(C\ne 1\) and let \(F'\) be another function satisfying (1) and inequality (8). Therefore, \(F'\) is of the form (cf., Theorem 1)

$$\begin{aligned} F'(x,y)=G'(x,y)+\Phi '(x)+\Psi '(y)+\delta ', \quad x,y \in X, \end{aligned}$$

with biadditive \(G'\), additive \(\Phi '\) and \(\Psi '\), satisfying (3), (4) and (5), respectively, and with \(\delta '=0\) in the case \(C\ne 1\).

We have for all \(x,y \in X\), \(n\in \mathbb {N},\)

$$\begin{aligned}&\Vert F(nx,ny)-F'(nx,ny)\Vert \le 28 \varepsilon ,\\&\Vert G(nx,ny)+\Phi (nx)+\Psi (ny) - G'(nx,ny)-\Phi '(nx)-\Psi '(ny)\Vert \le 28 \varepsilon ,\\&\Vert n^2\big (G(x,y)-G'(x,y)\big ) + n\big (\Phi (x)+\Psi (y)-\Phi '(x)-\Psi '(y)\big )\Vert \le 28\varepsilon . \end{aligned}$$

Dividing the above inequality by \(n^2\) side by side and letting n tend to infinity we obtain \(G=G'\), and consequently,

$$\begin{aligned} \Phi (x)+\Psi (y)=\Phi '(x)+\Psi '(y)\quad x,y \in X. \end{aligned}$$

It is now enough to set \(y=0\) and then \(x=0\) in order to obtain \(\Phi =\Phi '\) and \(\Psi =\Psi '\), respectively. \(\square \)

Remark 1

A thorough inspection of the proof of Theorem 2 shows that in the case \(C=1\) we are able to obtain a better approximation. Namely, if \( f:X^2 \rightarrow Y\) is a mapping satisfying (7) for \(x_1,x_2,y_1,y_2 \in X\) and \(C=1\), then there exists a solution \(F:X^2 \rightarrow Y\) of (1) such that

$$\begin{aligned} \Vert f(x,y)-f(0,0)-F(x,y)\Vert \le 11\varepsilon , \qquad x,y \in X. \end{aligned}$$
(24)

Remark 2

It is also easy to observe that in the case \(C=1\) we do not have the uniqueness of function F in (8). Indeed, each function \({\overline{F}}:X^2\rightarrow Y\),

$$\begin{aligned} \overline{F}:=G+\Phi +\Psi +\delta ' \end{aligned}$$

with \(G,\Phi , \Psi \) defined as in the proof of Theorem 2, and with \(\delta '\in Y\) such that

$$\begin{aligned} \Vert \delta '\Vert \le 3 \varepsilon \end{aligned}$$

satisfies, on account of Remark 1, conditions (1) and (8).

3 Generalized Stability of (1)

In this section we provide some results concerning generalized stability with various approximation functions. In what follows we will use a notation

$$\begin{aligned}&(\Phi f)(x_1,y_1,x_2,y_2) := \; f(a_1x_1+a_2x_2,b_1y_1+b_2y_2)\\&\quad -C_1f(x_1,y_1)-C_2f(x_1,y_2)-C_3f(x_2,y_1)-C_4f(x_2,y_2) \end{aligned}$$

for \(x_1,x_2,y_1,y_2 \in X\). Let us also denote \(a:=a_1+a_2\) and \(b:=b_1+b_2\).

Our first result reads as follows.

Theorem 3

Suppose that \((Y,\Vert \cdot \Vert )\) is a Banach space, \(C\ne 0\), \(a\ne 0,\) \(b\ne 0\). Let \(f:X^2\rightarrow Y\) and \(\theta :X^4\rightarrow \mathbb {R}_+\) be mappings satisfying the inequality

$$\begin{aligned}&\Vert (\Phi f)(x_1,y_1,x_2,y_2)\Vert \le \theta (x_1,y_1,x_2,y_2),\qquad x_1,x_2,y_1,y_2\in X. \end{aligned}$$
(25)

Assume, further, that for an \(s\in \{-1,1\}\) (depending on abC) we have

$$\begin{aligned} \varepsilon ^*(x,y)\!:=\! \sum _{n=0}^{\infty }\frac{\theta \big (a^{sn+\frac{s-1}{2}}x,b^{sn+\frac{s-1}{2}}y,a^{sn+\frac{s-1}{2}}x,b^{sn+\frac{s-1}{2}}y\big )}{|C|^{sn+\frac{s+1}{2}}} \!<\!\infty , \quad x,\,y\!\in \! X,\nonumber \\ \end{aligned}$$
(26)

and

$$\begin{aligned} \lim _{n\rightarrow \infty }\frac{\theta \big (a^{sn}x_1,b^{sn}y_1,a^{sn}x_2,b^{sn}y_2\big )}{|C|^{sn}} =0,\qquad x_1,x_2,y_1,y_2\in X. \end{aligned}$$
(27)

Then there exists a unique solution \(F:X^2\rightarrow Y\) of (1) such that

$$\begin{aligned} \Vert f(x,y)-F(x,y)\Vert \le \varepsilon ^*(x,y),\qquad x,y \in X. \end{aligned}$$
(28)

Proof

Putting \(x_1=x_2=x\) and \(y_1=y_2=y\) in (25) we get

$$\begin{aligned} \Vert f\big (ax,by\big )-Cf(x,y)\Vert \le \theta (x,y,x,y), \qquad x,y\in X, \end{aligned}$$

whence,

$$\begin{aligned} \bigg \Vert \frac{f\big (ax,by\big )}{C} - f(x,y)\bigg \Vert \le \frac{1}{|C|}\theta (x,y,x,y),\qquad x,y\in X. \end{aligned}$$
(29)

Similarly, putting \(x_1=x_2=\frac{x}{a}\) and \(y_1=y_2=\frac{y}{b}\) in (25) we obtain

$$\begin{aligned} \Big \Vert f(x,y) - Cf\Big (\frac{x}{a},\frac{y}{b}\Big )\Big \Vert \le \theta \Big (\frac{x}{a},\frac{y}{b},\frac{x}{a},\frac{y}{b}\Big ),\qquad x,y\in X. \end{aligned}$$
(30)

Define

$$\begin{aligned} (\mathcal {T}\xi )(x,y):=\frac{1}{C^s}\xi (a^sx,b^sy),\qquad \xi \in Y^{X^2},\;x,y \in X, \end{aligned}$$
(31)

and

$$\begin{aligned} \varepsilon (x,y):=\left\{ \begin{array}{ll} \dfrac{1}{|C|}\theta (x,y,x,y), &{} \text{ for } \; s=1,\\ &{}\\ \displaystyle \theta \Big (\frac{x}{a},\frac{y}{b},\frac{x}{a},\frac{y}{b}\Big ), &{} \text{ for } \; s=-1, \end{array} \right. \end{aligned}$$

for all \(x,y \in X\). Then, for any \(\xi ,\mu :{X^2}\rightarrow Y,\) \(x,y\in X\) we have

$$\begin{aligned} \Vert (\mathcal {T}\xi )(x,y)-(\mathcal {T}\mu )(x,y)\Vert =\frac{1}{|C|^s}\big \Vert \xi \big ( a^sx,b^sy \big )- \mu \big (a^sx,b^sy\big )\big \Vert \end{aligned}$$

and by (29) and (30),

$$\begin{aligned} \Vert (\mathcal {T}f)(x,y)-f(x,y)\Vert \le \varepsilon (x,y),\qquad x,y\in X. \end{aligned}$$

Next, put

$$\begin{aligned} (\Lambda \eta )(x,y):=\frac{1}{|C|^s}\eta \big (a^sx,b^sy\big ),\qquad \eta \in \mathbb {R}_+ ^{X^2}\!\!,\; x,y\in X. \end{aligned}$$

As one can check,

$$\begin{aligned} (\Lambda ^n\varepsilon )(x,y){=}\frac{\varepsilon \big (a^{sn}x,b^{sn}y\big )}{|C|^{sn}}{=} \left\{ \begin{array}{ll} \displaystyle \frac{\theta (a^{n}x,b^{n}y,a^{n}x,b^{n}y)}{|C|^{n+1}}, &{} \text{ for } \; s=1,\\ \displaystyle |C|^n\theta \Big (\frac{x}{a^{n+1}},\frac{y}{b^{n+1}},\frac{x}{a^{n+1}},\frac{y}{b^{n+1}}\Big ), &{} \text{ for } \; s=-1, \end{array} \right. \end{aligned}$$

for all \( x,y\in X,\; n\in \mathbb {N}_0\).

The operators \(\mathcal {T}:Y^{X^2}\rightarrow Y^{X^2}\) and \(\Lambda :\mathbb {R}_+ ^{X^2}\rightarrow \mathbb {R}_+ ^{X^2}\) satisfy the assumptions of Theorem 1 in [8], therefore, there exists a unique fixed point \(F:X^2\rightarrow Y\) of \(\mathcal {T}\) such that (28) holds. Moreover,

$$\begin{aligned} F(x,y)=\lim _{n\rightarrow \infty }(\mathcal {T}^nf)(x,y), \qquad x,y\in X. \end{aligned}$$
(32)

Now, we prove that for any \(x_1,x_2,y_1,y_2\in X\) and \(n \in \mathbb {N}_0\) we have

$$\begin{aligned} \big \Vert \big (\Phi (\mathcal {T}^n f)\big ) (x_1,y_1,x_2,y_2) \big \Vert \le \frac{\theta \big (a^{sn}x_1,b^{sn}y_1,a^{sn}x_2,b^{sn}y_2\big )}{|C|^{sn}}. \end{aligned}$$
(33)

Since the case \(n=0\) is just (25), fix an \(n \in \mathbb {N}_0\) and assume that (33) holds for any \(x_1,x_2,y_1,y_2\in X\). Then for any \(x_1,x_2,y_1,y_2\in X\) we get

$$\begin{aligned}&\big \Vert \big (\Phi (\mathcal {T}^{n+1} f)\big )(x_1,y_1,x_2,y_2) \big \Vert \\&=\big \Vert \big (\mathcal {T}(\mathcal {T}^{n}f)\big )(a_1x_1+a_2x_2,b_1y_1+b_2y_2)\\&\quad -C_1\big (\mathcal {T}(\mathcal {T}^{n}f)\big )(x_1,y_1)-C_2\big (\mathcal {T}(\mathcal {T}^{n}f)\big )(x_1,y_2)\\&\quad -C_3\big (\mathcal {T}(\mathcal {T}^{n}f)\big )(x_2,y_1) -C_4\big (\mathcal {T}(\mathcal {T}^{n}f)\big )(x_2,y_2) \big \Vert \\&=\Big \Vert \frac{1}{C^{s}} (\mathcal {T}^{n}f)\big (a^s(a_1x_1+a_2x_2),b^s(b_1y_1+b_2y_2)\big )\\&\quad -\,C_1\frac{1}{C^{s}} (\mathcal {T}^{n}f)(a^sx_1,b^sy_1) - C_2\frac{1}{C^{s}}(\mathcal {T}^{n}f)(a^sx_1,b^sy_2)\\&\quad - C_3\frac{1}{C^{s}}(\mathcal {T}^{n}f)(a^sx_2,b^sy_1) -\,C_4\frac{1}{C^{s}}(\mathcal {T}^{n}f)(a^sx_2,b^sy_2)\Big \Vert \\&=\frac{1}{|C|^{s}}\big \Vert \big (\Phi (\mathcal {T}^{n} f)\big )(a^{s}x_1,b^{s}y_1,a^{s}x_2,b^{s}y_2) \big \Vert \\&\le \frac{\theta \big (a^{s(n+1)}x_1,b^{s(n+1)}y_1,a^{s(n+1)}x_2,b^{s(n+1)}y_2\big )}{|C|^{s(n+1)}}, \end{aligned}$$

and thus, (33) holds for any \(x_1,x_2,y_1,y_2\in X\) and \(n \in \mathbb {N}_0\).

Letting \(n\rightarrow \infty \) in (33) and using (27) we finally obtain

$$\begin{aligned} (\Phi F)(x_1,y_1,x_2,y_2)=0,\qquad x_1,x_2,y_1,y_2\in X, \end{aligned}$$

which means that function F satisfies (1).

For the proof of uniqueness, assume that \(F'\) is another function satisfying (1) and (28). We have for all \(x,y \in X\), \(l\in \mathbb {N}_0\)

$$\begin{aligned}&\Vert F(x,y)-F'(x,y)\Vert \\ =&\;\Big \Vert \frac{1}{C^{sl}}F(a^{sl}x,b^{sl}y)-\frac{1}{C^{sl}}F'(a^{sl}x,b^{sl}y) \Big \Vert \\ \le&\; \frac{1}{|C|^{sl}}\big (\Vert F(a^{sl}x,b^{sl}y)- f(a^{sl}x,b^{sl}y)\Vert \\&+\Vert F'(a^{sl}x,b^{sl}y)-f(a^{sl}x,b^{sl}y) \Vert \big )\\ \le&\; 2\sum _{n=0}^{\infty }\frac{\theta \big (a^{s(l+n)+\frac{s-1}{2}}x,b^{s(l+n)+\frac{s-1}{2}}y,a^{s(l+n)+\frac{s-1}{2}}x,b^{s(l+n)+\frac{s-1}{2}}y\big )}{|C|^{s(n+l)+\frac{s+1}{2}}}\\ =&\;2\sum _{n=l}^{\infty }\frac{\theta \big (a^{sn+\frac{s-1}{2}}x,b^{sn+\frac{s-1}{2}}y,a^{sn+\frac{s-1}{2}}x,b^{sn+\frac{s-1}{2}}y\big )}{|C|^{sn+\frac{s+1}{2}}}, \end{aligned}$$

whence letting \(l \rightarrow \infty \) and using (26) we obtain \(F(x,y)=F'(x,y)\) for all \(x,y\in X\), which finishes the proof. \(\square \)

Theorem 3 with \(\theta (x_1,y_1,x_2,y_2):=\varepsilon > 0\) gives immediately the classical Hyers–Ulam stability result for (1). Namely, we have the following corollary.

Corollary 1

Let \((Y,\Vert \cdot \Vert )\) be a Banach space, \(\varepsilon >0\), \(C\ne 0\), \(|C|\ne 1\), \(a\ne 0\) and \(b \ne 0\). If \(f:X^2\rightarrow Y\) satisfies the inequality

$$\begin{aligned} \Vert (\Phi f)(x_1,y_1,x_2,y_2)\Vert \le \varepsilon ,\qquad x_1,x_2,y_1,y_2\in X, \end{aligned}$$

then there exists a unique solution \(F:X^2\rightarrow Y\) of (1) such that

$$\begin{aligned} \Vert f(x,y)-F(x,y)\Vert \le \frac{\varepsilon }{|1-|C||} ,\qquad x,y \in X. \end{aligned}$$

Proof

From (26) we have

$$\begin{aligned} \varepsilon ^*(x,y)&=\left\{ \begin{array}{ll} \displaystyle \sum _{n=0}^{\infty }\frac{\varepsilon }{|C|^{n+1}} , &{} \text{ for } \; |C|>1 \\ \displaystyle \sum _{n=0}^{\infty }|C|^n\varepsilon , &{} \text{ for } \; 0<|C|<1 \end{array} \right. =\;\left\{ \begin{array}{ll} \displaystyle \frac{\varepsilon }{|C|-1} , &{} \text{ for } \; |C|>1\\ \dfrac{\varepsilon }{1-|C|}, &{} \text{ for } \; 0<|C|<1 \end{array} \right. \\&=\frac{\varepsilon }{|1-|C||}, \quad \text{ for }\; \; C\in \mathbb {R}\setminus \{-1,0,1\}. \end{aligned}$$

\(\square \)

Remark 3

Studying the proof of Theorem 3 one can make several observations:

  • We do not demand that the coefficients \(a_1,a_2,b_1,b_2\) are non-zero.

  • If \(C=0\) then for \(\varepsilon ^*\) in (26) to be well defined we take \(s=-1\). If also \(a\ne 0,\,b \ne 0\), then in Theorem 3, f satisfies the condition

    $$\begin{aligned} \Vert f(x,y)\Vert \le \theta \Big ( \frac{x}{a}, \frac{y}{b}, \frac{x}{a}, \frac{y}{b}\Big ), \qquad x,y \in X, \end{aligned}$$

    and in Corollary 1, f is bounded by \(\varepsilon \). Both, in the theorem and in the corollary, we have then

    $$\begin{aligned} F(x,y)=\lim _{n\rightarrow \infty }(\mathcal {T}^nf)(x,y)= \lim _{n\rightarrow \infty } C^{n}f\Big (\frac{x}{a^n},\frac{y}{b^n}\Big ) =0, \qquad x,y\in X. \end{aligned}$$
  • If \(a=0=b\) (and \(|C|>1\), for (26) to be satisfied), we take \(s=1\), and we have

    $$\begin{aligned} \Big \Vert f(x,y)-\frac{f(0,0)}{C}\Big \Vert \le \frac{1}{|C|}\theta (x,y,x,y), \;\;x,y \in X,\;\; \end{aligned}$$
    (34)

    in Theorem 3, and with \(\theta (x,y,x,y)=\varepsilon \), in Corollary 1.

    Then

    $$\begin{aligned} F(x,y)=\lim _{n\rightarrow \infty }(\mathcal {T}^nf)(x,y)= \lim _{n\rightarrow \infty } \frac{1}{C^{n}}f(0,0) =0. \end{aligned}$$

    From (34), it follows that in Theorem 3, f is majorized by the function

    $$\begin{aligned} X^2\ni (x,y)\mapsto \frac{1}{|C|}\theta (x,y,x,y)+\frac{\theta (0,0,0,0)}{|C-1||C|}, \end{aligned}$$

    and in Corollary 1, it is simply bounded.

  • If \(a=0\) and \(b\ne 0\) (and \(|C|>1\)) then \(s=1\) and the approximating function F depends only on one variable

    $$\begin{aligned} F(x,y)=\lim _{n\rightarrow \infty }(\mathcal {T}^nf)(x,y)= \lim _{n\rightarrow \infty }\frac{1}{C^ n}f(0,b^ny), \qquad x,y\in X . \end{aligned}$$

    Analogous approach we have for \(a\ne 0\) and \(b=0\).

  • If \(|C|>1\) then \(s=1\), and Corollary 1 coincides with the result of Ciepliñski from [10].

Theorem 4

Let \((Y,\Vert \cdot \Vert )\) be a Banach space. Assume that \(f:X^2\rightarrow Y\) and \(\theta :X^4\rightarrow \mathbb {R}_+\) are mappings satisfying inequality (25) and the conditions

$$\begin{aligned} \varepsilon ^*(x,y):=\sum _{n=0}^{\infty }\sum _{i+j+k+l=n}{n \atopwithdelims (){i,j,k,l}}|C_1|^i|C_2|^j|C_3|^k|C_4|^{l} \delta _1^{(i,j,k,l)}(x,y,x,y) <\infty ,\quad \end{aligned}$$
(35)

for \(x,y\in X\) and

$$\begin{aligned} \lim _{n\rightarrow \infty }\sum _{i+j+k+l=n}{n \atopwithdelims (){i,j,k,l}}|C_1|^i|C_2|^j|C_3|^k|C_4|^{l}\delta _0^{(i,j,k,l)}(x_1,y_1,x_2,y_2) =0,\qquad \end{aligned}$$
(36)

for \(x_1,x_2,y_1,y_2\in X,\) where

\(\delta _m^{(i,j,k,l)}(x_1,y_1,x_2,y_2):=\)

$$\begin{aligned}&\theta \Big (\frac{x_1}{(2a_1)^{i+j+m}(2a_2)^{k+l}}, \frac{y_1}{(2b_1)^{i+k+m}(2b_2)^{j+l}},\\&\quad \frac{x_2}{(2a_1)^{i+j}(2a_2)^{k+l+m}}, \frac{y_2}{(2b_1)^{i+k}(2b_2)^{j+l+m}} \Big ). \end{aligned}$$

Then there exists a unique solution \(F:X^2\rightarrow Y\) of (1) such that condition (28) holds.

Proof

Putting \(x_1=\frac{x}{2a_1},\) \(x_2=\frac{x}{2a_2}\), \(y_1=\frac{y}{2b_1}\) and \(y_2=\frac{y}{2b_2}\) in (25) (with \(x,y \in X\)) we get

$$\begin{aligned}&\Big \Vert f(x,y)- C_1f\Big (\frac{x}{2a_1},\frac{y}{2b_1}\Big )- C_2f\Big (\frac{x}{2a_1},\frac{y}{2b_2}\Big )\nonumber \\&\quad \quad - C_3f\Big (\frac{x}{2a_2},\frac{y}{2b_1}\Big )- C_4f\Big (\frac{x}{2a_2},\frac{y}{2b_2}\Big )\Big \Vert \nonumber \\&\quad \le \theta \Big (\frac{x}{2a_1},\frac{y}{2b_1},\frac{x}{2a_2},\frac{y}{2b_2}\Big ), \qquad x,y\in X. \end{aligned}$$
(37)

Define

$$\begin{aligned} (\mathcal {T}\xi )(x,y):=C_1\xi \Big (\frac{x}{2a_1},\frac{y}{2b_1}\Big )+ C_2\xi \Big (\frac{x}{2a_1},\frac{y}{2b_2}\Big )+C_3\xi \Big (\frac{x}{2a_2},\frac{y}{2b_1}\Big )\\ +C_4\xi \Big (\frac{x}{2a_2},\frac{y}{2b_2}\Big ), \qquad \xi \in Y^{X^2},\; x,y\in X, \end{aligned}$$

and

$$\begin{aligned} \varepsilon (x,y):=\theta \Big (\frac{x}{2a_1},\frac{y}{2b_1},\frac{x}{2a_2},\frac{y}{2b_2}\Big ),\qquad x,y\in X. \end{aligned}$$

Then, by (37), we obtain

$$\begin{aligned} \Vert (\mathcal {T}f)(x,y)-f(x,y)\Vert \le \varepsilon (x,y),\qquad x,y\in X. \end{aligned}$$

Put also

$$\begin{aligned} (\Lambda \eta )(x,y):=|C_1|\eta \Big (\frac{x}{2a_1},\frac{y}{2b_1}\Big )+ |C_2|\eta \Big (\frac{x}{2a_1},\frac{y}{2b_2}\Big )+ |C_3|\eta \Big (\frac{x}{2a_2},\frac{y}{2b_1}\Big )\\ +|C_4|\eta \Big (\frac{x}{2a_2},\frac{y}{2b_2}\Big ),\qquad \eta \in \mathbb {R}_+ ^{X^2},\; x,y\in X. \end{aligned}$$

Now, using induction, we show that for any \(n \in \mathbb {N}_0\) and \(x,y\in X\) we have

$$\begin{aligned} \begin{array}{ll} (\Lambda ^n\varepsilon )(x,y) = &{} \displaystyle \sum _{i+j+k+l=n}{n \atopwithdelims (){i,j,k,l}}|C_1|^{i}|C_2|^{j}|C_3|^{k}|C_4|^{l}\\ &{} \displaystyle \times \,\varepsilon \Big (\Big (\frac{1}{2a_1}\Big )^{i+j}\Big (\frac{1}{2a_2}\Big )^{k+l}\!\!\!x, \Big (\frac{1}{2b_1}\Big )^{i+k}\Big (\frac{1}{2b_2}\Big )^{j+l}\!\!\!y\Big ). \end{array} \end{aligned}$$
(38)

Fix \(x,y\in X\). Clearly, (38) is true for \(n=0\). Next, fix an \(n \in \mathbb {N}_0\) and assume that (38) holds. Then

$$\begin{aligned}&(\Lambda ^{n+1}\varepsilon )(x,y)=\big (\Lambda (\Lambda ^n\varepsilon )\big )(x,y)\\&\quad =\;|C_1|(\Lambda ^n\varepsilon )\Big (\frac{x}{2a_1},\frac{y}{2b_1}\Big )+ |C_2|(\Lambda ^n\varepsilon )\Big (\frac{x}{2a_1},\frac{y}{2b_2}\Big )\\&\quad \quad + |C_3|(\Lambda ^n\varepsilon )\Big (\frac{x}{2a_2},\frac{y}{2b_1}\Big )+ |C_4|(\Lambda ^n\varepsilon )\Big (\frac{x}{2a_2},\frac{y}{2b_2}\Big )\\&\quad = \sum _{i+j+k+l=n}\!{n \atopwithdelims (){i,j,k,l}}|C_1|^{i+1}|C_2|^j|C_3|^k|C_4|^{l} \\&\quad \quad \times \varepsilon \bigg (\!\frac{x}{(2a_1)^{i+j+1}(2a_2)^{k+l}}, \frac{y}{(2b_1)^{i+k+1}(2b_2)^{j+l}}\!\bigg )\\&\quad \quad +\sum _{i+j+k+l=n}\!{n \atopwithdelims (){i,j,k,l}}|C_1|^{i}|C_2|^{j+1}|C_3|^k|C_4|^{l} \\&\quad \quad \times \varepsilon \bigg (\!\frac{x}{(2a_1)^{i+j+1}(2a_2)^{k+l}}, \frac{y}{(2b_1)^{i+k}(2b_2)^{j+l+1}}\!\bigg )\\&\quad \quad +\sum _{i+j+k+l=n}\!{n \atopwithdelims (){i,j,k,l}}|C_1|^{i}|C_2|^{j}|C_3|^{k+1}|C_4|^{l}\\&\quad \quad \times \varepsilon \bigg (\!\frac{x}{(2a_1)^{i+j}(2a_2)^{k+l+1}}, \frac{y}{(2b_1)^{i+k+1}(2b_2)^{j+l}}\!\bigg )\\&\quad \quad +\sum _{i+j+k+l=n}\!{n \atopwithdelims (){i,j,k,l}}|C_1|^{i}|C_2|^{j}|C_3|^{k}|C_4|^{l+1}\\&\quad \quad \times \varepsilon \bigg (\!\frac{x}{(2a_1)^{i+j}(2a_2)^{k+l+1}}, \frac{y}{(2b_1)^{i+k}(2b_2)^{j+l+1}}\!\bigg )\\&\quad = \sum _{i+j+k+l=n+1}\!{n +1\atopwithdelims (){i,j,k,l}}|C_1|^i|C_2|^j|C_3|^k|C_4|^{l}\\&\quad \quad \times \varepsilon \bigg (\!\frac{x}{(2a_1)^{i+j}(2a_2)^{k+l}}, \frac{y}{(2b_1)^{i+k}(2b_2)^{j+l}}\!\bigg ), \end{aligned}$$

and thus (38) is true for any \(n \in \mathbb {N}_0\) and \(x,y\in X\).

One can now show that the operators \(\mathcal {T}:Y^{X^2}\rightarrow Y^{X^2}\) and \(\Lambda :{\mathbb {R}_+}^{X^2}\rightarrow {\mathbb {R}_+}^{X^2}\) satisfy the assumptions of Theorem 1 in [8] and therefore there exists a unique fixed point \(F:X^2\rightarrow Y\) of \(\mathcal {T}\) such that (28) holds. Moreover, F is given by (G).

We prove that for any \(x_1,x_2,y_1,y_2\in X\) and \(n \in \mathbb {N}_0\) we have

$$\begin{aligned}&\big \Vert \big (\Phi (\mathcal {T}^n f)\big ) (x_1,y_1,x_2,y_2) \big \Vert \nonumber \\&\quad \le \sum _{i+j+k+l=n}{n \atopwithdelims (){i,j,k,l}}|C_1|^i|C_2|^j|C_3|^k|C_4|^{l}\delta _0^{(i,j,k,l)}(x_1,y_1,x_2,y_2). \end{aligned}$$
(39)

Since the case \(n=0\) is just (25), fix an \(n \in \mathbb {N}_0\) and assume that (39) holds for any \(x_1,x_2,y_1,y_2\in X\). Then for any \(x_1,x_2,y_1,y_2\in X\) we get

$$\begin{aligned}&\big \Vert \big (\Phi (\mathcal {T}^{n+1} f)\big )(x_1,y_1,x_2,y_2) \big \Vert =\big \Vert \big (\mathcal {T}(\mathcal {T}^{n}f)\big )(a_1x_1+a_2x_2,b_1y_1+b_2y_2)\\&\quad \quad -C_1\big (\mathcal {T}(\mathcal {T}^{n}f)\big )(x_1,y_1)-C_2\big (\mathcal {T}(\mathcal {T}^{n}f)\big )(x_1,y_2)\\&\quad \quad -C_3\big (\mathcal {T}(\mathcal {T}^{n}f)\big )(x_2,y_1)-C_4\big (\mathcal {T}(\mathcal {T}^{n}f)\big )(x_2,y_2)\big \Vert \\&\quad =\Big \Vert C_1(\mathcal {T}^{n}f)\Big (\frac{a_1x_1+a_2x_2}{2a_1},\frac{b_1y_1+b_2y_2}{2b_1}\Big )\\&\quad \quad +\, C_2(\mathcal {T}^{n}f)\Big (\frac{a_1x_1+a_2x_2}{2a_1},\frac{b_1y_1+b_2y_2}{2b_2}\Big )\\&\quad \quad +\, C_3(\mathcal {T}^{n}f)\Big (\frac{a_1x_1+a_2x_2}{2a_2},\frac{b_1y_1+b_2y_2}{2b_1}\Big )\\&\quad \quad +\, C_4(\mathcal {T}^{n}f)\Big (\frac{a_1x_1+a_2x_2}{2a_2},\frac{b_1y_1+b_2y_2}{2b_2}\Big )\\&\quad \quad -\, C_1\Big (C_1(\mathcal {T}^{n} f)\Big (\frac{x_1}{2a_1},\frac{y_1}{2b_1}\Big )+ C_2 (\mathcal {T}^{n}f)\Big (\frac{x_1}{2a_1},\frac{y_1}{2b_2}\Big )\\&\quad \quad +\, C_3 (\mathcal {T}^{n}f)\Big (\frac{x_1}{2a_2},\frac{y_1}{2b_1}\Big ) + \, C_4 (\mathcal {T}^{n}f)\Big (\frac{x_1}{2a_2},\frac{y_1}{2b_2}\Big )\Big )\\&\quad \quad - \, C_2\Big (C_1(\mathcal {T}^{n} f)\Big (\frac{x_1}{2a_1},\frac{y_2}{2b_1}\Big ) +\, C_2 (\mathcal {T}^{n}f)\Big (\frac{x_1}{2a_1},\frac{y_2}{2b_2}\Big )\\&\quad \quad +\, C_3 (\mathcal {T}^{n}f)\Big (\frac{x_1}{2a_2},\frac{y_2}{2b_1}\Big )+ C_4 (\mathcal {T}^{n}f)\Big (\frac{x_1}{2a_2},\frac{y_2}{2b_2}\Big )\Big )\\&\quad \quad -\, C_3\Big (C_1(\mathcal {T}^{n} f)\Big (\frac{x_2}{2a_1},\frac{y_1}{2b_1}\Big )+\, C_2 (\mathcal {T}^{n}f)\Big (\frac{x_2}{2a_1},\frac{y_1}{2b_2}\Big )\\&\quad \quad +\, C_3 (\mathcal {T}^{n}f)\Big (\frac{x_2}{2a_2},\frac{y_1}{2b_1}\Big ) +\, C_4 (\mathcal {T}^{n}f)\Big (\frac{x_2}{2a_2},\frac{y_1}{2b_2}\Big )\Big )\\&\quad \quad -\, C_4\Big (C_1(\mathcal {T}^{n} f)\Big (\frac{x_2}{2a_1},\frac{y_2}{2b_1}\Big )+ C_2 (\mathcal {T}^{n}f)\Big (\frac{x_2}{2a_1},\frac{y_2}{2b_2}\Big )\\&\quad \quad +\, C_3 (\mathcal {T}^{n}f)\Big (\frac{x_2}{2a_2},\frac{y_2}{2b_1}\Big )+ C_4 (\mathcal {T}^{n}f)\Big (\frac{x_2}{2a_2},\frac{y_2}{2b_2}\Big )\Big )\Big \Vert \\&\quad \le |C_1|\Big \Vert \big (\Phi (\mathcal {T}^{n}f)\big )\Big (\frac{x_1}{2a_1},\frac{y_1}{2b_1},\frac{x_2}{2a_1},\frac{y_2}{2b_1}\Big ) \Big \Vert \\&\quad \quad +\, |C_2|\Big \Vert \big (\Phi (\mathcal {T}^{n}f)\big )\Big (\frac{x_1}{2a_1},\frac{y_1}{2b_2},\frac{x_2}{2a_1},\frac{y_2}{2b_2}\Big ) \Big \Vert \\&\quad \quad +\, |C_3|\Big \Vert \big (\Phi (\mathcal {T}^{n}f)\big )\Big (\frac{x_1}{2a_2},\frac{y_1}{2b_1},\frac{x_2}{2a_2},\frac{y_2}{2b_1}\Big ) \Big \Vert \\&\quad \quad +\, |C_4|\Big \Vert \big (\Phi (\mathcal {T}^{n}f)\big )\Big (\frac{x_1}{2a_2},\frac{y_1}{2b_2},\frac{x_2}{2a_2},\frac{y_2}{2b_2}\Big ) \Big \Vert \\&\quad \le \sum _{i+j+k+l=n}{n \atopwithdelims (){i,j,k,l}}|C_1|^{i+1}|C_2|^j|C_3|^k|C_4|^{l} \, \\&\quad \quad \times \delta _0^{(i,j,k,l)} \Big (\frac{x_1}{2a_1},\frac{y_1}{2b_1},\frac{x_2}{2a_1},\frac{y_2}{2b_1}\Big )\\&\quad \quad +\,\sum _{i+j+k+l=n}{n \atopwithdelims (){i,j,k,l}}|C_1|^{i}|C_2|^{j+1}|C_3|^k|C_4|^{l}\, \\&\quad \quad \times \delta _0^{(i,j,k,l)} \Big (\frac{x_1}{2a_1},\frac{y_1}{2b_2},\frac{x_2}{2a_1},\\&\quad \quad +\, \sum _{i+j+k+l=n}{n \atopwithdelims (){i,j,k,l}}|C_1|^{i}|C_2|^j|C_3|^{k+1}|C_4|^{l} \, \delta _0^{(i,j,k,l)}\\&\quad \quad \times \delta _0^{(i,j,k,l)}\Big (\frac{x_1}{2a_2},\frac{y_1}{2b_1},\frac{x_2}{2a_2},\frac{y_2}{2b_1}\Big )\\&\quad \quad +\, \sum _{i+j+k+l=n}{n \atopwithdelims (){i,j,k,l}}|C_1|^{i}|C_2|^j|C_3|^k|C_4|^{l+1} \, \\&\quad \quad \times \delta _0^{(i,j,k,l)} \Big (\frac{x_1}{2a_2},\frac{y_1}{2b_2},\frac{x_2}{2a_2},\frac{y_2}{2b_2}\Big )\\&\quad =\sum _{i+j+k+l=n+1}{n+1 \atopwithdelims (){i,j,k,l}}|C_1|^i|C_2|^j|C_3|^k|C_4|^{l} \, \delta _0^{(i,j,k,l)}(x_1,y_1,x_2,y_2). \end{aligned}$$

We have thus shown that (39) holds for \(x_1,x_2,y_1,y_2\in X\) and \(n \in \mathbb {N}_0\). Letting \(n\rightarrow \infty \) in (39) and using (36) we see that

$$\begin{aligned} (\Phi F)(x_1,y_1,x_2,y_2)=0,\qquad x_1,x_2,y_1,y_2\in X, \end{aligned}$$

which means that function F satisfies (1).

For the proof of uniqueness, assume that \(F'\) is another function satisfying (1) and (28). Then, for any \(m \in \mathbb {N}\) we have

$$\begin{aligned}&\!\!\!\!\!\!\Vert F(x,y)-F'(x,y)\Vert \\ =&\Big \Vert \sum _{i+j+k+l=m}{m \atopwithdelims (){i,j,k,l}}C_1^{i}C_2^{j}C_3^kC_4^{l} \\&\quad \times \Big [ F \Big (\frac{x}{(2a_1)^{i+j}(2a_2)^{k+l}}, \frac{y}{(2b_1)^{i+k}(2b_2)^{j+l}}\Big ) \\&\quad - F' \Big (\frac{x}{(2a_1)^{i+j}(2a_2)^{k+l}}, \frac{y}{(2b_1)^{i+k}(2b_2)^{j+l}}\Big )\Big ]\Big \Vert \\ \le&\sum _{i+j+k+l=m}{m \atopwithdelims (){i,j,k,l}}|C_1|^{i}|C_2|^{j}|C_3|^k|C_4|^{l} \\&\quad \times \Big \Vert F \Big (\frac{x}{(2a_1)^{i+j}(2a_2)^{k+l}}, \frac{y}{(2b_1)^{i+k}(2b_2)^{j+l}}\Big )\\&\quad - F' \Big (\frac{x}{(2a_1)^{i+j}(2a_2)^{k+l}}, \frac{y}{(2b_1)^{i+k}(2b_2)^{j+l}}\Big )\Big \Vert \\ \le&\sum _{i+j+k+l=m}{m \atopwithdelims (){i,j,k,l}}|C_1|^{i}|C_2|^{j}|C_3|^k|C_4|^{l} \\&\quad \times \bigg (\Big \Vert F \Big (\frac{x}{(2a_1)^{i+j}(2a_2)^{k+l}}, \frac{y}{(2b_1)^{i+k}(2b_2)^{j+l}}\Big )\\&\quad - f \Big (\frac{x}{(2a_1)^{i+j}(2a_2)^{k+l}}, \frac{y}{(2b_1)^{i+k}(2b_2)^{j+l}}\Big )\Big \Vert \\&\quad +\Big \Vert f \Big (\frac{x}{(2a_1)^{i+j}(2a_2)^{k+l}}, \frac{y}{(2b_1)^{i+k}(2b_2)^{j+l}}\Big )\\&\quad \quad - F' \Big (\frac{x}{(2a_1)^{i+j}(2a_2)^{k+l}}, \frac{y}{(2b_1)^{i+k}(2b_2)^{j+l}}\Big )\Big \Vert \bigg )\\ =&\;2 \sum _{n=0}^{\infty } \sum _{i+j+k+l=m}\sum _{\overline{i}+\overline{j}+\overline{k}+\overline{l}=n}{m \atopwithdelims (){i,j,k,l}}{n \atopwithdelims (){\overline{i},\overline{j},\overline{k},\overline{l}}}|C_1|^{i+\overline{i}}|C_2|^{j+\overline{j}}|C_3|^{k+\overline{k}}|C_4|^{l+\overline{l}} \\&\quad \times \theta \Big (\frac{x}{(2a_1)^{i+\overline{i}+j+\overline{j}+1}(2a_2)^{k+\overline{k}+l+\overline{l}}}, \frac{y}{(2b_1)^{i+\overline{i}+k+\overline{k}+1}(2b_2)^{j+\overline{j}+l+\overline{l}}}, \\&\qquad \qquad \qquad \frac{x}{(2a_1)^{i+\overline{i}+j+\overline{j}+1}(2a_2)^{k+\overline{k}+l+\overline{l}}}, \frac{y}{(2b_1)^{i+\overline{i}+k+\overline{k}+1}(2b_2)^{j+\overline{j}+l+\overline{l}}}\Big )\\ =&\; 2 \sum _{n=0}^{\infty } \;\sum _{i+j+k+l=n+m}{n+m \atopwithdelims (){i,j,k,l}}|C_1|^{i}|C_2|^{j}|C_3|^{k}|C_4|^{l} \\&\quad \times \theta \Big (\frac{x}{(2a_1)^{i+j+1}(2a_2)^{k+l}}, \frac{y}{(2b_1)^{i+k+1}(2b_2)^{j+l}}, \\&\qquad \qquad \qquad \frac{x}{(2a_1)^{i+j+1}(2a_2)^{k+l}}, \frac{y}{(2b_1)^{i+k+1}(2b_2)^{j+l}}\Big )\\ =&\; 2 \sum _{n=m}^{\infty } \;\sum _{i+j+k+l=n}{n \atopwithdelims (){i,j,k,l}}|C_1|^{i}|C_2|^{j}|C_3|^{k}|C_4|^{l} \\&\quad \times \theta \Big (\frac{x}{(2a_1)^{i+j+1}(2a_2)^{k+l}}, \frac{y}{(2b_1)^{i+k+1}(2b_2)^{j+l}},\\&\qquad \qquad \qquad \frac{x}{(2a_1)^{i+j+1}(2a_2)^{k+l}}, \frac{y}{(2b_1)^{i+k+1}(2b_2)^{j+l}}\Big ). \end{aligned}$$

Tending now with m to infinity, on the account of the assumption, it follows that \(F=F'\), which completes the proof. \(\square \)

Theorem 4 with \(\theta (x_1,y_1,x_2,y_2):=\varepsilon >0\) gives immediately the following corollary on the classical Hyers–Ulam stability of (1).

Corollary 2

Let \((Y,\Vert \cdot \Vert )\) be a Banach space, \(\varepsilon > 0\) and \(|C_1|+|C_2|+|C_3|+|C_4| <1\). If \(f:X^2\rightarrow Y\) satisfies the inequality

$$\begin{aligned} \Vert (\Phi f)(x_1,y_1,x_2,y_2)\Vert \le \varepsilon ,\qquad x_1,x_2,y_1,y_2\in X, \end{aligned}$$

then there exists a solution \(F:X^2\rightarrow Y\) of (1) such that

$$\begin{aligned} \Vert f(x,y)-F(x,y)\Vert \le \frac{\varepsilon }{1-(|C_1|+|C_2|+|C_3|+|C_4|)} ,\qquad x,y \in X. \end{aligned}$$