\def\papertitle{Existence  of random coupled system of fractional differential equations in generalized Banach space with retarded and advanced arguments}
\def\paperauthorA{\underline{Fouad Fredj$^1$}, Hadda Hammouche$^1$}
\def\paperauthorB{Hadda Hammouche}
% you can add more authors using alphabets after \paperauthor
% 
%------------------------------------------------------------------------------------------
% Do not delete this part starting from here and up to .....
\documentclass[oneside]{article}
\usepackage{icma2021,amssymb,amsmath,mathrsfs} 
%\setcounter{page}{1}
%\input{../../expages.tex}\setpagenumber{1}
%\input{expages.tex}\setpagenumber{1}
\ninept

\usepackage{mathptmx} 
% Saves a lot of ouptut space in PDF... after conversion with the distiller	  
% Delete if you cannot get PS fonts working on your system.

% pdf-tex settings: detect automatically if run by latex or pdflatex
\newif\ifpdf
\ifx\pdfoutput\relax
\else
\ifcase\pdfoutput
\pdffalse
\else
\pdftrue
\fi

\ifpdf % compiling with pdflatex
\usepackage[pdftex,
pdftitle={\papertitle},
%pdfauthor={\paperauthorA},
% for more than one author use pdfauthor={\paperauthorA,\paperauthorB, .....},
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%ù
colorlinks=false, % links are activated as colror boxes instead of color text
bookmarksnumbered, %use section numbers with bookmarks
pdfstartview= XYZ	% starts with zoom=100% instead of full screen; especially useful if working with a big screen :-)
]{hyperref}
\pdfcompresslevel=9
\usepackage[pdftex]{graphicx}
\usepackage[figure,table]{hypcap}
\else % compiling with latex
\usepackage[dvips]{graphicx}
\usepackage[dvips,
colorlinks=false,           % no color links
bookmarksnumbered, %use section numbers with bookmarks
pdfstartview= XYZ	% starts with zoom=100% instead of full screen
]{hyperref}
% hyperrefs are active in the pdf file after conversion
\usepackage[figure,table]{hypcap}
\fi

\newcommand\dsum{\displaystyle\sum}
\newcommand\dprod{\displaystyle\prod}

\newcommand{\ww}{\widetilde}
\newcommand{\hh}{\widehat}
\newcommand{\bb}{\overline}
\newcommand{\cc}{\cdot}

\newtheorem{Theorem}{Theorem}[section]
\newtheorem{lem}[Theorem]{Lemma}
%\theoremstyle{plain}
%\theoremstyle{definition}
\newtheorem{defi}[Theorem]{Definition}
\newtheorem{ex}{Example}
%\newtheorem*{pro}{Property}

%------------------------------------------------------------------------------------------
\title{\papertitle}
\affiliation{\paperauthorA}	% This command replaces \name{The DAFx Crew} 
{$1-$Mathematics and Applied Sciences Laboratory,\\ Ghardaia University, Ghardaia 47000, Algeria.}

%\usepackage{layout}
%------------------------------------------------------------------------------------------
%%%%%%%%%%%%%%%%%%%%%%% to here %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{document}
	
	%\layout
	
	% more pdf-tex settings:
	\ifpdf % used graphic file format for pdflatex
	\DeclareGraphicsExtensions{.png,.jpg,.pdf}
	\else  % used graphic file format for latex
	\DeclareGraphicsExtensions{.eps}
	\fi
	
	\maketitle

	
\begin{abstract}
This paper aims to investigate  a random coupled system with multiple fractional derivatives of $\psi-$Caputo with retarded and advanced arguments, and proves the uniqueness of random solution by applying random versions of the Pervo fixed point theorem, while the existence of  solutions is derived by a random version of a Krasnoselskii-type fixed point theorem. Obtained results are supported by examples and illustrated in the last section
\end{abstract}
{\textbf {Keywords:}}$\psi-$ fractional derivative, random coupled system, generalized Banach spaces, retarded arguments, advanced arguments. \\
{\bf {Mathematics Subject Classification:}} 26A33; 45D05.

\section{Introduction}
In this paper, we investigate the existence and uniqueness  of the following  nonlinear random coupled system of  $\psi$-Caputo fractional integro-differential equations
\begin{equation} \label{pb}
\left\{ \begin{array}{ll}

^cD^{\alpha_1;\psi} x(t,\upsilon) + \sum_{i=1}^{m}I^{\gamma_{1,i};\psi}p_{1,i}(t,x^t(\upsilon),y^t(\upsilon),\upsilon)=q_1(t,x^t(\upsilon),y^t(\upsilon),\upsilon)\\
\\
^cD^{\alpha_2;\psi} y(t,\upsilon) + \sum_{i=1}^{m}I^{\gamma_{2,i};\psi}p_{2,i}(t,x^t(\upsilon),y^t(\upsilon),\upsilon)=q_2(t,x^t(\upsilon),y^t(\upsilon),\upsilon)\\
\end{array} \right.; t\in J,\upsilon\in\Omega.
\end{equation}

\begin{equation} \label{cond}
\left\{ \begin{array}{ll}

(x(t,\upsilon), y(t,\upsilon))=(\eta_1(t,\upsilon),\eta_2(t,\upsilon)),\quad t\in [a-r,a],r>0,\upsilon\in \Omega \\
\\
(x(t,\upsilon), y(t,\upsilon))=(\xi_1(t,\upsilon),\xi_2(t,\upsilon)),\quad t\in [T,T+l],l>0,\upsilon\in \Omega \\
\end{array} \right. \quad ;\upsilon\in\Omega,
\end{equation}
where $J=[a,T]$, $D^{\alpha_j;\psi}$ denotes the $\psi$-Caputo fractional derivative of order $1< \alpha_j\leq 2 $, $I^{\gamma_{i,j};\psi}$  is  the $\psi$-Riemann-Liouville fractional integral of orders $\gamma_{i,j} > 0$, $(\Omega,\mathcal{A})$ is measurable space, $p_i,q_{i,j}: J\times C([-r,l],\mathbb{R}^n)\times C([-r,l],\mathbb{R}^n)\times\Omega\rightarrow\mathbb{R}^m$ are given functions,$\eta_j\in C([a-r,a],\mathbb{R}^n) $ with $\eta_j(a,\upsilon)=0$ and  $\xi\in C([T,T+l],\mathbb{R}^n)$ with $\xi(T,\upsilon)=0$; $j=1,2,i=1\cdots m$.
\section{Preliminaries}
In this section, we introduce some preliminaries and lemmas that will be used throughout this paper and we will prove an auxiliary lemma, which plays a key role in defining a fixed point problem associated with the given problem. By $C([-r,l],\mathbb{R}^n)$ we denote the Banach space of all continuous functions from $[-r,l]$ into $\mathbb{R}^n$ provided with the norm

\[\|x\|_{[-r,l]}=\sup\lbrace\|x(t)\|:-r\leq t\leq l\rbrace.\]
and $C([a,T],\mathbb{R}^n)$ is Banach space equipped with norm 
\[\|x\|_{[a,T]}=\sup\lbrace\|x(t)\|:a\leq t\leq T\rbrace.\]
$AC(J,\mathbb{R}^n)$ is the space of absolutely continuous functions on $J$, and we denote by $AC^1(J,\mathbb{R}^n)$ the space of functions $x(t)$ which have continuous derivatives on $J$:
\[AC^1(J,\mathbb{R}^n)=\lbrace\|x\rightarrow\mathbb{R}^n:x'\in AC^1(J,\mathbb{R}^n)\rbrace .\]
where \[x'(t)=t\dfrac{d}{dt}g(t),\quad t\in J.\]
we denote the space $X$ by
\[X=\lbrace x:[a-r,T+l]\rightarrow \mathbb{R}^n:x_{|[a-r,a]}\in C([a-r,a]),x_{|[a,T]}\in AC([a,T]) \text{ and }x_{|[T,T+l]}\in C([T,T+l])\rbrace\]
where the aforementioned space are supplemented with the following norms
\[\|x\|_{[a-r,a]}=\sup\lbrace\|x(t)\|:a-r\leq t\leq a\rbrace,\]
\[\|x\|_{[T,T+l]}=\sup\lbrace\|x(t)\|:T\leq t\leq T+l\rbrace,\]
\[\|x\|_X=\sup\lbrace\|x(t)\|:a-r\leq t\leq T+l\rbrace.\]
The product space $X\times X$ is provided 
with the norm
\[\|(x,y)\|_{X\times X}:= \|x\|_X + \|y\|_X.\]
\begin{defi}\cite{Almeida}
Let $\alpha>0$ and an increasing function $\psi:J\rightarrow \mathbb{R}$ satisfy $\psi'(t)\neq 0$ for
all $t \in J$. The $\psi$-Riemann–Liouville fractional integral of a function $x:[a, T] \rightarrow \mathbb{R}$ is defined by

\[I^{\alpha;\psi}x(t)=\dfrac{1}{\Gamma(\alpha)}\int^t_a \psi'(s)(\psi(t) - \psi(s))^{\alpha} x(s)ds,\qquad 0<a<s<t.\]

\end{defi}


\begin{defi}\cite{Almeida}
%Let $\psi'(x)\neq 0 (- \infty\leq a<t<b \leq\infty)$ and 
The $\psi$-Caputo fractional derivative of order $\alpha > 0$ for a function $x \in C^n[0, \infty)$ is defined by
\[ ^c D^{\alpha;\psi}x(t)=\dfrac{1}{\Gamma(n-\alpha)}\int^t_a \psi'(s)(\psi(t) - \psi(s))^{n-\alpha-1}D_{\psi}^n x(s)ds,\qquad 0<a<s<t.\]

where $n=[p]+1$ and $D_{\psi}^n = \bigg(\dfrac{1}{\psi'(t)}\dfrac{d}{dt}\bigg)^n$.
\end{defi}

\begin{lem}\label{lm:1}\cite{Almeida}
Let $p>0$. The following holds 

\begin{itemize} 
\item If $x\in C(J,\mathbb{R})$, then

\[^cD^{\alpha;\psi} I^{\alpha;\psi} x(t)= x(t),\quad t\in J,\]

\item If $x\in C^n(J,\mathbb{R})$, $n-1<p<n$, then

\[I^{\alpha;\psi} {^cD^{\alpha;\psi}} x(t)= x(t)-\sum^{n-1}_{k=0} c_k \Psi^k_0(t),\quad t\in J.\]

where, $\Psi^{k}_0(t) = (\psi(t) - \psi(0))^{k},\Psi_0(t) = (\psi(t) - \psi(0))$ and $c_k = \dfrac{D_{\psi}^{k}x(0)}{k!}$.
\end{itemize} 

\end{lem}            

\begin{lem}
Let $ 1 < \alpha < 2 $. For any functions $  F , P_i \in C(J,\mathbb{R})$,  $\eta\in C([a-r,a],\mathbb{R}) $ with $\eta(a)=0$ and  $\xi\in C([T,T+l],\mathbb{R}) $ with $\xi(T)=0$. Then the following linear problem

\begin{align}\label{lipb}
\left\lbrace \begin{array}{lllll}
\begin{split}
^c\mathit{D}^{\alpha;\psi}x(t) +\sum_{i=1}^m I^{\gamma;\psi}P_i(t) = Q(t) , \quad t \in J ,\\
u(t)=\eta(t), \qquad\qquad t\in[a-r,a], r>0,\\
u(t)=\xi(t), \qquad\qquad t\in [T,T+l], l>0.
\end{split}
\end{array} \right.
\end{align}
\end{lem}
has a unique solution, which is given by
%\begin{equation}\label{sol}
%u(t)=\dfrac{1}{\Gamma(\alpha)}\left\lbrace\begin{array}{lll}
%\eta(t),\text{ if } t\in [a-r,a],\\
%%\displaystyle{\frac{1}{\Gamma(\alpha+\gamma)}\sum_{i=1}^m\int^T_a G_{\gamma+\alpha}(t,s)P_i(s)ds - \frac{1}{\Gamma(\alpha)}\int^T_a G_{\alpha}(t,s)F(s)ds},\text{ if } t\in J,\\
%\xi(t),\text{ if } t\in [T,T+l],
%\end{array}\right.
%\end{equation}

\begin{align}\label{sol}
\begin{split}
x(t)=\dfrac{1}{\Gamma(\alpha)}\left\lbrace\begin{array}{lll}
 & \eta(t),\quad\text{ if } t\in [a-r,a],\\
 \\
 & \displaystyle{I^{\alpha;\psi}Q(t)- \sum_{i=1}^m I^{\gamma+\alpha;\psi} P_i(t)}\\ 
 & \displaystyle{ + \dfrac{\Psi_a(t)}{\Psi_a(T)}\bigg(\sum_{i=1}^m I^{\gamma+\alpha;\psi} P_i(T) - I^{\alpha;\psi}Q(T)}\bigg),\quad\text{ if } t\in J,\\
 \\
 & \xi(t),\quad\text{ if } t\in [T,T+l].
\end{array}\right.
\end{split}
\end{align}
%where
%\begin{align}\label{grn1}
%\begin{split}
%  G_{\alpha}&(t,s)=\\
% &\dfrac{1}{\Gamma(\alpha)}\left\lbrace\begin{array}{lll}
%\frac{\big(\psi(t)-\psi(a)\big)\psi'(s)\big(\psi(T)-\psi(a)\big)^{\alpha-1}}{\psi(T)-\psi(a)}-\psi'(s)\big(\psi(t)-\psi(a)\big)^{\alpha-1}, a \leq s \leq t \leq T,\\
%%\frac{\big(\psi(t)-\psi(a)\big)\psi'(s)\big(\psi(T)-\psi(a)\big)^{\alpha-1}}{\psi(T)-\psi(a)},\quad a \leq t \leq s \leq T.
%\end{array}\right.
%\end{split}
%\end{align}
%and
%\begin{equation}\label{grn2}
%G_{\alpha}(t,s)=\dfrac{1}{\Gamma(\gamma_i+\alpha)}\left\lbrace\begin{array}{lll}
%%\frac{\big(\psi(t)-\psi(a)\big)\psi'(s)\big(\psi(T)-\psi(a)\big)^{\gamma_i+\alpha-1}}{\psi(T)-\psi(a)}-\psi'(s)\big(\psi(t)-\psi(a)\big)^{\gamma_i+\alpha-1}, a \leq s \leq t \leq T,\\
%%\frac{\big(\psi(t)-\psi(a)\big)\psi'(s)\big(\psi(T)-\psi(a)\big)^{\gamma_i+\alpha-1}}{\psi(T)-\psi(a)},\quad a \leq t \leq s \leq T.
%\end{array}\right.
%\end{equation}
\begin{proof}
Applying the$\psi$-Riemann-Liouville fractional integral of order $\alpha$ to both side of the equation in \eqref{lipb}, and using Lemma \ref{lm:1}, we get


\begin{align}\label{Eq:I}
\begin{split}
&  x(t)= I^{\alpha;\psi}Q(t)- \sum_{i=1}^m I^{\gamma+\alpha;\psi} P_i(t) +c_0 + c_1\Psi_a(t) ;\qquad c_0,c_1\in\mathbb{R}.
 \end{split}
\end{align}
Using the fact that $x(a)=\eta(a)=0$ , $x(T)=\xi(T)=0$, and from \eqref{Eq:I}, we find

\[x(a)=c_0=0\]
and
\[x(T)=I^{\alpha;\psi}Q(T)- \sum_{i=1}^m I^{\gamma+\alpha;\psi} P_i(T) + c_1\Psi_a(t)=0 \]

Some simple computations give us

\[c_1=\dfrac{1}{\Psi_a(T)}\bigg(\sum_{i=1}^m I^{\gamma+\alpha;\psi} P_i(T) - I^{\alpha;\psi}Q(T)\bigg)\]
Inserting $c_0$ and $c_1$ in eqution \eqref{Eq:I}, which leads to solution \eqref{sol}.

\end{proof}

\par Let $\mathcal{A}$ and $\mathcal{B}$ be two $\sigma$-algebra of Borel subsets of $\mathbb{R}^m$.

\begin{defi}
Let $(\mathbb{R}^m,\mathcal{A})$ and $(\mathbb{R}^m,\mathcal{B})$ be two measurable spaces. A mapping $F:(\mathbb{R}^m,\mathcal{A})\rightarrow(\mathbb{R}^m,\mathcal{B})$  is said to be measurable if the $\sigma$-algebra of Borel $F^{-1}(\mathcal{B})\subset\mathcal{A}$ i.e.
\[ F^{-1}(A)\subset\mathcal{A}\qquad\text{ for all } A\subset\mathcal{B}.\]
\end{defi}

\begin{defi}\label{def1}\cite{Sincer}
A function $g:\Omega\times \mathbb{R}^m\rightarrow \mathbb{R}^m$ is called jointly measurable if $g(\cc, y)$ is measurable for all $y\in\mathbb{R}^m$ and $g(z,\cc)$ is continuous for all $z\in \Omega$.
\end{defi}

\begin{defi}\cite{Sincer}\label{def2}
A function $g:J\times \Omega\times\mathbb{R}^m\rightarrow \mathbb{R}^m$ is  Carathéodory if the following conditions are satisfied:

 \begin{itemize}

   \item[(i)] $(t,y)\rightarrow g(t,y,z)$ is jointly measurable for any $z\in\mathbb{R}^m$ and

   \item[(ii)] $z\rightarrow g(t,y,z)$ is continuous for any $t\in J$ and $y\in\Omega$.

 \end{itemize}

\end{defi}

\begin{defi}\cite{Graef}
A mapping $Q:\Omega\times Y\rightarrow Y$ is Said to be a random operator if, for any $y\in Y$, $Q(.,y)$ is measurable.
\end{defi}

\begin{defi}\cite{Sincer}
A random fixed point of a random operator $Q$ is a measurable function $g:\Omega\rightarrow Y$ such that

\[g(y)=Q(y,g(y)),\qquad\text{ for all }y\in\Omega.\]\\

\end{defi}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
\ \ \ Let $y,z\in\mathbb{R}^n$ with $y=(y_1,\cdots,y_n)$ and $z=(z_1,\cdots,z_n)$, by $y\leq z$ we mean $y_i\leq z_i$ for all $i=1,\cdots,n.$ Also we set $|y|=(|y_1|,\cdots,|y_n|)$, $\max(y,z)=(\max(y_1,z_1),\cdots,\max(y_n,z_n))$ and $\mathbb{R}_+^n=\lbrace y\in\mathbb{R}^n:y_i>0, \forall i, 1\leq i\leq n\rbrace$. If $c\in\mathbb{R}$, then $y\leq c$ means $y_i\leq c$ for each $i=1,\cdots,n.$

\begin{defi}\cite{Sincer}
Let $Y$ be a nonempty set. By a vector-valued metric on $Y$, we mean a map $d:Y\times Y\rightarrow\mathbb{R}^n$ satisfying

\begin{itemize}
\item[(i)] $d(y,z)\geq 0$ for all $y,z\in Y$ if $d(y,z)=0$, then $y=z$;\\

\item[(ii)] $d(y,z)=d(z,y)$ for all $y,z\in Y$;\\

\item[(iii)] $d(y,z)\leq d(y,w)+d(w,z)$ for all $y,z,w\in Y$.\\

\end{itemize}


\end{defi}

\ \ \ We call the pair $(Y,d)$ a generalized metric space with
\[\left(\begin{array}{llllll}
d_1(y,z)\\
\quad.\\
\quad.\\
\quad.\\
d_n(y,z)
\end{array}\right)\]

Notice that $d$ is a generalized metric space on $Y$ if and only if $d_i;$ $i=1,...,n,$ are metrics on $Y$.\\

for $r=(r_1,\ldots,r_n)\in\mathbb{R}_+^n$, we let

\[B(y_0,r)=\lbrace y\in Y:d(y_0,y)<r\rbrace,\]

denote the open ball centered at $y_0$ with radius $r$, and

\[\bb{B(y_0,r)}=\lbrace y\in Y:d(y_0,y)<r\rbrace.\]
be the closed ball centered at $y_0$ with radius $r$. We mention that for a generalized set metric space the notions of the open set, closed set, convergence, Cauchy sequence and completeness are similar to those in the usual metric space.

\begin{defi}\cite{Varga}
A matrix $M=(a_{ij})_{1\leq i,j\leq n}\in M_{n\times n}(\mathbb{R})$ is said to be convergent to zero if and only if its spectral radius $\rho(A)$ is strictly less then one. In other words this means that all the  eigenvalues of $A$ are in the open unit disc, i.e., $|\lambda|<1$; for every $\lambda\in\mathbb{C}$ with $det(A-\lambda I)=0$; where $I$ denotes the unit matrix of $M_{n\times n}(\mathbb{R})$.

\end{defi}

\begin{defi}\cite{Graef}
Let $(Y,d)$ be a generalized metric space. An operator $Q:X\rightarrow X$ is said to be contractive if there exist a matrix $M$ convergent to zero such that

\[d(Q(y),Q(z))\leq M d(y,z)\qquad\text{ for all }y,z\in Y.\]
\end{defi}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{Theorem}\label{Banach}\cite{Sincer}
Let $(\Omega,\mathcal{B})$ be a measurable space, let $Y$ be a real separable generalized Banach space and let $Q :\Omega \times Y \rightarrow Y$ be a continuous random operator. Let $ M(\vartheta)\in M_{n\times n}(\mathbb{R_+})$ be a random variable matrix such that for every $\vartheta\in\Omega$  the matrix $M(\vartheta)$ converges to 0 and 
\[d(Q(\vartheta,y_1),Q(\vartheta,y_2)) \leq M(\vartheta)d(y_1,y_2)\] 
for each $y_1,y_2\in Y$, $\vartheta\in\Omega$.\\
Then there exists a random variable $y :\Omega\rightarrow Y$ which is the unique random fixed point of $Q$.
\end{Theorem}


\begin{lem}\label{kras}\cite{Sincer}

Let $ Y $ be a separable generalized Banach space $ Y $. Suppose that $ A , B : \Omega\times Y \rightarrow Y$ are random operators such that:
\begin{enumerate}
\item A is a continuous random and $M(\vartheta)$-contraction operator,
\item B is a completely continuous random operator,
\item the matrix $I-M$ has the absolute value property.
if
\[ \Lambda=\bigg\lbrace y:\Omega\rightarrow Y \text{ is measurable }\big |\zeta(\vartheta)A(y,\vartheta)+\zeta(\vartheta)B\big(\frac{y}{\zeta(\vartheta)},\vartheta\big)=y\bigg\rbrace \]
\end{enumerate}
is bounded for all measurable mappings $\zeta :\Omega\rightarrow\mathbb{R}$ with $0<\zeta(\vartheta)<1$ on $\Omega$, then the random equation
\[y=A(y,\vartheta)+B(y,\vartheta),\qquad y\in Y,\]
has at least one solution.
\end{lem}

\section{Main results}
                                
\begin{lem}\label{EQI}
For given functions $ p_{j,i},q_j \in C(J,\mathbb{R}^n) $, $j = 1,2$ and $i = 1,\cdots,m$. A functions $x,y \in C^2$ is a random solution of systems \eqref{pb}-\eqref{cond} if and only if $x,y$ satisfies the following random coupled system integral equations

\begin{align}\label{xt}
\begin{split}
& x(t,\upsilon)=\\
 & \dfrac{1}{\Gamma(\alpha_j)}\left\lbrace\begin{array}{lll}
 & \eta_j(t,\upsilon),\quad\text{ if } t\in [a-r,a],\\
 \\
 & \displaystyle{I^{\alpha_j;\psi}q_j(t,x^t(\upsilon),y^t(\upsilon),\upsilon)- \sum_{i=1}^m I^{\gamma_{j,i}+\alpha;\psi} p_{j,i}(t,x^t(\upsilon),y^t(\upsilon),\upsilon) + \dfrac{\Psi_a(t)}{\Psi_a(T)}}\\ 
 & \displaystyle{\times\bigg(\sum_{i=1}^m I^{\gamma_{j,i}+\alpha;\psi} p_{j,i}(T,x^T(\upsilon),y^T(\upsilon),\upsilon) - I^{\alpha_j;\psi}q_j(T,x^T(\upsilon),y^T(\upsilon),\upsilon)}\bigg),\text{ if } t\in J,\\
 \\
 & \xi_j(t,\upsilon),\quad\text{ if } t\in [T,T+l].
\end{array}\right.
\end{split}
\end{align}
%and 
%\begin{align}\label{yt}
%\begin{split}
%& y(t,\upsilon)=\\
% & \dfrac{1}{\Gamma(\alpha)}\left\lbrace\begin{array}{lll}
% & \eta_2(t,\upsilon),\quad\text{ if } t\in [a-r,a],\\
% \\
% & \displaystyle{I^{\beta;\psi}q_2(t,x^t(\upsilon),y^t(\upsilon),\upsilon)- \sum_{i=1}^m I^{\zeta_i+\beta;\psi} p_{2,i}(t,x^t(\upsilon),y^t(\upsilon),\upsilon) + \dfrac{\psi(t)-\psi(a)}{\psi(T)-\psi(a)}}\\ 
% & \displaystyle{\times\bigg(\sum_{i=1}^m I^{\zeta_i+\beta;\psi} p_{2,i}(T,x^T(\upsilon),y^T(\upsilon),\upsilon) - I^{\beta;\psi}q_2(T,x^T(\upsilon),y^T(\upsilon),\upsilon)}\bigg),\text{ if } t\in J,\\
% \\
% & \xi_2(t,\upsilon),\quad\text{ if } t\in [T,T+l].
%\end{array}\right.
%\end{split}
%\end{align}

Let $X $ be the Banach space of continuous real-valued functions defined on $J$
\[ X = \lbrace x(t,\upsilon): x(t,\upsilon) \in C(J,\mathbb{R}^n)\rbrace, \]
 endowed with the norm
 \[ \|x(\cc,\upsilon)\|_X=\sup_{t\in J}\big\|x(t,\upsilon)\big\|\] 

To define a fixed point problem equivalent to the system \eqref{pb}-\eqref{cond}, we introduce the operator

\[T : J \times X \times X \times \Omega \rightarrow X \times X \]
defined by 

\[\big(T(x,y)\big)(t,\upsilon)=\left(\begin{array}{lll}
\big(T_1(x,y)\big)(t,\upsilon)\\
\\
\big(T_2(x,y)\big)(t,\upsilon)\\ 
\end{array} \right) \]
where


\begin{align}\label{Tj}
\begin{split}
& \big(T_j(x,y)\big)(t,\upsilon)=\\
 & \dfrac{1}{\Gamma(\alpha_j)}\left\lbrace\begin{array}{lll}
 & \eta_j(t,\upsilon),\quad\text{ if } t\in [a-r,a],\\
 \\
 & \displaystyle{I^{\alpha_j;\psi}q_j(t,x^t(\upsilon),y^t(\upsilon),\upsilon)- \sum_{i=1}^m I^{\gamma_{j,i}+\alpha_j;\psi} p_{j,i}(t,x^t(\upsilon),y^t(\upsilon),\upsilon) + \dfrac{\Psi_a(t)}{\Psi_a(T)}}\\ 
 & \displaystyle{\times\bigg(\sum_{i=1}^m I^{\gamma_{j,i}+\alpha_j;\psi} p_{j,i}(T,x^T(\upsilon),y^T(\upsilon),\upsilon) - I^{\alpha_j;\psi}q_j(T,x^T(\upsilon),y^T(\upsilon),\upsilon)}\bigg),\text{ if } t\in J,\\
 \\
 & \xi_j(t,\upsilon),\quad\text{ if } t\in [T,T+l].
\end{array}\right.
\end{split}
\end{align}
%and 
%\begin{align}\label{yt}
%\begin{split}
%& \big(T_2(x,y)\big)(t,\upsilon)=\\
% & \dfrac{1}{\Gamma(\alpha)}\left\lbrace\begin{array}{lll}
% & \eta_2(t,\upsilon),\quad\text{ if } t\in [a-r,a],\\
% \\
% & \displaystyle{I^{\beta;\psi}q_2(t,x^t(\upsilon),y^t(\upsilon),\upsilon)- \sum_{i=1}^m I^{\zeta_i+\beta;\psi} p_{2,i}(t,x^t(\upsilon),y^t(\upsilon),\upsilon) + \dfrac{\psi(t)-\psi(a)}{\psi(T)-\psi(a)}}\\ 
% & \displaystyle{\times\bigg(\sum_{i=1}^m I^{\zeta_i+\beta;\psi} p_{2,i}(T,x^T(\upsilon),y^T(\upsilon),\upsilon) - I^{\beta;\psi}q_2(T,x^T(\upsilon),y^T(\upsilon),\upsilon)}\bigg),\text{ if } t\in J,\\
% \\
% & \xi_2(t,\upsilon),\quad\text{ if } t\in [T,T+l].
%\end{array}\right.
%\end{split}
%\end{align}
\end{lem}

\subsection{ Uniqueness of solutions}
\par The first result is concerned with the  uniqueness of random solution for the system \eqref{pb}-\eqref{cond} and is based on random versions of the Pervo fixed point theorem.
\begin{Theorem}\label{unic}
We assume that the following hypotheses holds
\begin{itemize}

\item[(H1)] 
 	The functions $p_{j,i}$ and $q_j$ are Carathéodory; $j = 1,2$ and $i = 1,\cdot,m$.

\item[(H2)] There exist measurable functions $ \mathcal{K}_j,\mathcal{L}_j,\mathcal{M}_{j,i},\mathcal{N}_{j,i} : \Omega \rightarrow (0,\infty) ; j=1,2$ and $i = 1,\cdot,m$ such that:
	\[\|q_j(t,x,y,\upsilon) - q_j(t,\bb x,\bb y,\upsilon)\| \leq \mathcal{K}_j(t,\upsilon)\| x - \bb x \|_{[-r,l]} + \mathcal{L}_j(t,\upsilon)\| y - \bb y \|_{[-r,l]},\]
and

\[\|p_{j,i}(t,x,y,\upsilon) - p_{j,i}(t,\bb x,\bb y,\upsilon)\| \leq \mathcal{M}_{j,i}(t,\upsilon)\| x - \bb x \|_{[-r,l]} + \mathcal{N}_{j,i}(t,\upsilon)\| y - \bb y \|_{[-r,l]},\]

for a.e.$t\in J$, and each $x, y, \bb x, \bb y \in C([-r,l])$.

\item [(H3)] $ M(\upsilon) \in M_{n \times n}(\mathbb{R_+})$ is random variable matrix, such that for every $\upsilon \in \Omega$, the matrix
\begin{align*}
 M(\upsilon)= 2\left( \begin{array}{lll}
C_{\alpha_1}\mathcal{K}^{*}_1(\upsilon) + C_{\gamma_1+\alpha}\mathcal{M}^{*}_1(\upsilon) \qquad 
C_{\alpha_1}\mathcal{L}^{*}_1(\upsilon) + C_{\gamma_1+\alpha_1}\mathcal{N}^{*}_1(\upsilon)\\
\\
C_{\alpha_2}\mathcal{K}^{*}_2(\upsilon) + C_{\gamma_2+\alpha_2}\mathcal{M}^{*}_2(\upsilon)\qquad 
C_{\alpha_2}\mathcal{L}^{*}_2(\upsilon) + C_{\gamma_2+\alpha_2}\mathcal{N}^{*}_2(\upsilon)
\end{array} \right).
\end{align*}
converges to zero.
with $C_{\alpha_j}=\dfrac{\Psi_{a}^{\alpha}(T)}{\Gamma(\alpha_j+1)}, C_{\alpha_j+\gamma_j}=\dfrac{\Psi_{a}^{\alpha}(T)}{\Gamma(\gamma_j+\alpha_j+1)}$; $j=1,2,$
 and

\[\mathcal{K}^{*}_j(\upsilon)=\|\mathcal{K}_j(\cc,\upsilon)\|_{[a,T]},\mathcal{L}^{*}_j(\upsilon)=\|\mathcal{L}_j(\cc,\upsilon)\|_{[a,T]} ,\mathcal{M}^{*}_j(\upsilon)=\sum^m_{i=1}\|\mathcal{M}_{j,i}(\cc,\upsilon)\|_{[a,T]}\],\[\mathcal{N}^{*}_j(\upsilon)=\sum^m_{i=1}\|\mathcal{N}_{j,i}(\cc,\upsilon)\|_{[a,T]},\gamma_j=\sup_i\lbrace\gamma_{j,i}:i=1,\cdots m,\rbrace\]
\end{itemize}

 Then the coupled system \eqref{pb}-\eqref{cond} has a unique random solution.\\

\end{Theorem}


\begin{proof}

First, we need to show that the operator $T$ is a random operator on $X\times X$.\\ 
From (H1) and Definition \ref{def2} the maps $p_{j,i}$ and $q_j$  are measurable with respect to the variable $\upsilon$. In the view of the Definition \ref{def1}, we conclude that the maps \[\upsilon \rightarrow T_1(x, y)(t, \upsilon) \text{ and } \vartheta \rightarrow T_2(x, y)(t, \upsilon)\]
are measurable. As a result, the operator $T$ is a random operator on $X \times X \times \Omega$ into $X \times X$.\\

 Next, we prove that the operator $ T $ is contractive.\\ 
For all $ \upsilon \in \Omega$,  $(x,y),(\bb x,\bb y) \in X \times X$, and $t\in J$, we have
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\begin{align*}
 & \big\| \big(T_1(x,y)\big)(t,\upsilon) - \big(T_1(\bb x,\bb y)\big)(t,\upsilon)\big\|\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 & = I^{\alpha_1;\psi}\big\|q_1(s,x^s(\upsilon),y^s(\upsilon),\upsilon)-q_1(s,\bb x^s(\upsilon),\bb y^s(\upsilon),\upsilon)\big\|(t)\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 & + \sum_{i=1}^m I^{\gamma_{1,i}+\alpha_1;\psi} \big\|p_{1,i}(s,x^s(\upsilon),y^s(\upsilon),\upsilon) - p_{1,i}(s,\bb x^s(\upsilon),\bb y^s(\upsilon),\upsilon)\big\|(t) \\ 
 & + \dfrac{\Psi_a(t)}{\Psi_a(T)}\bigg(\sum_{i=1}^m I^{\gamma_{1,i}+\alpha_1;\psi} \big\|p_{1,i}(s,x^s(\upsilon),y^s(\upsilon),\upsilon) - p_{1,i}(s,\bb x^s(\upsilon),\bb y^s(\upsilon),\upsilon)\big\|(T)\\ 
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
 & + I^{\alpha_1;\psi}\big\|q_1(s,x^s(\upsilon),y^s(\upsilon),\upsilon) - q_1(s,\bb x^s(\upsilon),\bb y^s(\upsilon),\upsilon)\big\|(T)\bigg)\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%                                                   
 & \leq 2\bigg(\dfrac{\Psi_{a}^{\alpha_1}(T)}{\Gamma(\alpha_1+1)}\mathcal{K}^{*}_1(\upsilon) + \dfrac{\Psi_{a}^{\gamma_{1,i}+\alpha_1}(T)}{\Gamma(\gamma_{1,i}+\alpha_1+1)}\mathcal{M}^{*}_1(\upsilon)\bigg)\|x(\cc,\upsilon) - \bb x(\cc,\upsilon)\|_X\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%                                                   
 & + 2\bigg(\dfrac{\Psi_{a}^{\alpha_1}(T)}{\Gamma(\alpha_1+1)}\mathcal{L}^{*}_1(\upsilon) + \dfrac{\Psi_{a}^{\gamma_{1,i}+\alpha_1}(T)}{\Gamma(\gamma_{1,i}+\alpha_1+1)}\mathcal{N}^{*}_1(\upsilon)\bigg)\|y(\cc,\upsilon) - \bb y(\cc,\upsilon)\|_X
\end{align*}

Therefore, for each $t\in J$, and $\upsilon\in \Omega$,
\begin{align*}
 & \big\| \big(T_1(x,y)\big)(\cc,\upsilon) - \big(T_1(\bb x,\bb y)\big)(\cc,\upsilon)\big\|_X\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%                                                   
 & \leq 2\bigg(\dfrac{\Psi_{a}^{\alpha_1}(T)}{\Gamma(\alpha_1+1)}\mathcal{K}^{*}_1(\upsilon) + \dfrac{\Psi_{a}^{\gamma_{1,i}+\alpha_1}(T)}{\Gamma(\gamma_{1,i}+\alpha_1+1)}\mathcal{M}^{*}_1(\upsilon)\bigg)\|x(\cc,\upsilon) - \bb x(\cc,\upsilon)\|_X\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%                                                   
 & + 2\bigg(\dfrac{\Psi_{a}^{\alpha_1}(T)}{\Gamma(\alpha_1+1)}\mathcal{L}^{*}_1(\upsilon) + \dfrac{\Psi_{a}^{\gamma_{1,i}+\alpha_1}(T)}{\Gamma(\gamma_{1,i}+\alpha_1+1)}\mathcal{N}^{*}_1(\upsilon)\bigg)\|y(\cc,\upsilon) - \bb y(\cc,\upsilon)\|_X
\end{align*}


In a similar way, we can find that

\begin{align*}
 & \big\| \big(T_2(x,y)\big)(\cc,\upsilon) - \big(T_2(\bb x,\bb y)\big)(\cc,\upsilon)\big\|_X\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 & \leq 2\bigg(\dfrac{\Psi_{a}^{\alpha_2}(T)}{\Gamma(\alpha_2+1)}\mathcal{K}^{*}_2(\upsilon) + \dfrac{\Psi_{a}^{\gamma_{2,i}+\alpha_2}(T)}{\Gamma(\gamma_{2,i}+\alpha_2+1)}\mathcal{M}^{*}_2(\upsilon)\bigg)\|x(\cc,\upsilon) - \bb x(\cc,\upsilon)\|_X\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%                                                   
 & + 2\bigg(\dfrac{\Psi_{a}^{\alpha_2}(T)}{\Gamma(\alpha_2+1)}\mathcal{L}^{*}_2(\upsilon) + \dfrac{\Psi_{a}^{\gamma_{2,i}+\alpha_2}(T)}{\Gamma(\gamma_{2,i}+\alpha_2+1)}\mathcal{N}^{*}_2(\upsilon)\bigg)\|y(\cc,\upsilon) - \bb y(\cc,\upsilon)\|_X
\end{align*}

Thus,

\[ d\Big(\big(T(x , y)\big)(\cc,\upsilon) , \big(T(\bb x , \bb y)\big)(\cc,\upsilon) \Big)\leq M(\upsilon)
d\Big(\big(x(\cc,\upsilon) , y(\cc,\upsilon)\big) ,\big(\bb x(\cc,\upsilon) , \bb y(\cc,\upsilon)\big)\Big) \]
where
\[ d\Big(\big(x(\cc,\upsilon) , x(\cc,\upsilon)\big) , \big(\bb y(\cc,\upsilon) ,\bb y(\cc,\upsilon)\big)\Big) = \left( \begin{array}{ll}
\|x(\cc,\upsilon) - \bb x(\cc,\upsilon)\|_X\\
\\
\|y(\cc,\upsilon) - \bb y(\cc,\upsilon)\|_X
\end{array} \right) \]

\par  As for every $\upsilon \in \Omega$, the matrix $ M(\upsilon) $ converges to zero, this implies that the operator $ T $ is a $ M(\upsilon)- $contractive operator. Consequently, by theorem \ref{Banach}, we conclude that $ Q $  has a unique fixed point, which is a random solution of systems \eqref{pb}-\eqref{cond}. This completes the proof.


\end{proof}

\subsection{Existence of solutions}

\par In the next result, we prove the existence of solution for the system \eqref{pb}-\eqref{cond} by applying a random version of a Krasnoselskii-type fixed point theorem.
\newline
\begin{Theorem}\label{thm:1}
Assume that (H1)-(H2) and the following hypotheses holds. 

\item[(H4)] there exist measurable functions $ \varphi_j,\chi_j,\omega_j,\lambda_{j,i}, \rho_{j,i}, \mu_{j,i} : J \rightarrow (0,X) ; i=1,2$ and $i = 1,\cdot,m$ such that:
\[\|q_j(t,x,y,\upsilon)\| \leq \varphi_j(t,\upsilon) + \chi_j(t,\upsilon)\| x \|_{[-r,l]} + \omega_j(t,\upsilon)\| y \|_{[-r,l]},\]

\[\|p_{j,i}(t,x,y,\upsilon) \| \leq \lambda_{j,i}(t,\upsilon) + \rho_{j,i}(t,\upsilon)\| x  \|_{[-r,l]} + \mu_{j,i}(t,\upsilon)\| y \|_{[-r,l]},\]
	
for a.e.$t\in J$, and each $x, y \in \mathbb{R}^m $.
\item [(H5)] $ \ww M(\upsilon) \in M_{n \times n}(\mathbb{R}_n)$ is random variable matrix, such that for every $\upsilon \in \Omega$, the matrix
\[ \ww M(\upsilon)=
\left( \begin{array}{lll}
\dfrac{\Psi_{a}^{\alpha_1}(T)}{\Gamma(\alpha_1+1)}\mathcal{K}^{*}_1(\upsilon) \qquad  \dfrac{\Psi_{a}^{\alpha_1}(T)}{\Gamma(\alpha_1+1)}\mathcal{L}^{*}_1(\upsilon)\\
\dfrac{\Psi_{a}^{\alpha_2}(T)}{\Gamma(\alpha_2+1)}\mathcal{K}^{*}_2(\upsilon) \qquad \dfrac{\Psi_{a}^{\alpha_2}(T)}{\Gamma(\alpha_2+1)}\mathcal{L}^{*}_2(\upsilon)
\end{array} \right).\]

converges to zero.\\
Then the coupled system \eqref{pb}-\eqref{cond} has at least a random solution.\\
             %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\end{Theorem}
\begin{proof}

 Let us subdivide the operator $ T $ into two operators $A,B: X \times X \times \Omega \rightarrow X \times X$ as follows:

\[\big((T(x,y)\big)(t,\upsilon)= \big(A(x,y)\big)(t,\upsilon) + \big(B(x,y)\big)(t,\upsilon) \qquad (x,y) \in X \times X,(t,\upsilon)\in J\times\Omega\]
where
\[\big(A(x,y)\big)(t,\upsilon) = \Big(\big(A_1(x,y)\big)(t,\upsilon),\big(A_2(x,y)\big)(t,\upsilon)\Big)\]
and
\[\big(B(x,y)\big)(t,\upsilon) = \Big(\big(B_1(x,y)\big)(t,\upsilon),\big(B_2(x,y)\big)(t,\upsilon)\Big)\]
with


\begin{align}\label{Aj}
\begin{split}
& \big(A_j(x,y)\big)(t,\upsilon)=\\
 & \dfrac{1}{\Gamma(\alpha_j)}\left\lbrace\begin{array}{lll}
 & \eta_j(t,\upsilon),\quad\text{ if } t\in [a-r,a],\\
 \\
 & I^{\alpha_j;\psi}q_j(t,x^t(\upsilon),y^t(\upsilon),\upsilon) - \dfrac{\Psi_a(t)}{\Psi_a(T)} I^{\alpha_j;\psi}q_j(T,x^T(\upsilon),y^T(\upsilon),\upsilon)\bigg),\text{ if } t\in J,\\
 \\
 & \xi_j(t,\upsilon),\quad\text{ if } t\in [T,T+l].
\end{array}\right.
\end{split}
\end{align}
and
\begin{align}\label{Bj}
\begin{split}
& \big(B_j(x,y)\big)(t,\upsilon)=\\
 & \dfrac{1}{\Gamma(\alpha_j)}\left\lbrace\begin{array}{lll}
 & 0,\quad\text{ if } t\in [a-r,a],\\
 \\
 & \displaystyle{\dfrac{\Psi_a(t)}{\Psi_a(T)}}\sum_{i=1}^m I^{\gamma_{j,i}+\alpha_j;\psi} p_{j,i}(T,x^T(\upsilon),y^T(\upsilon),\upsilon) - \sum_{i=1}^m I^{\gamma_{j,i}+\alpha_j;\psi} p_{j,i}(t,x^t(\upsilon),y^t(\upsilon),\upsilon) ,\text{ if } t\in J,\\
 \\
 & 0,\quad\text{ if } t\in [T,T+l].
\end{array}\right.
\end{split}
\end{align}


We need to prove that the operators $ A $ and $ B $ satisfies all conditions of the Theorem \ref{kras}.

The proof is divided into several steps.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{itemize}
\item[\textbf{step 1.}] $A$ is $\ww M(\upsilon)-$contraction operator:
\end{itemize}
As in the previous proof of Theorem \ref{thm:1}, we can obtain\\
For all $\upsilon \in \Omega$,  $(x,y) , (\bb x,\bb y) \in X \times X$, and $t\in J$\\
\begin{align*}
 &\Big\| \big(A_1(u,v)\big)(\cc,\vartheta) - \big(A_1(\bb u,\bb v)\big)(\cc,\vartheta) \Big\|_X \\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 & \leq 2\dfrac{\Psi_{a}^{\alpha_1}(T)}{\Gamma(\alpha_1+1)}\mathcal{K}^{*}_1(\upsilon)\|x(\cc,\upsilon) - \bb x(\cc,\upsilon)\|_X + 2\dfrac{\Psi_{a}^{\alpha_1}(T)}{\Gamma(\alpha_1+1)}\mathcal{L}^{*}_1(\upsilon)\|y(\cc,\upsilon) - \bb y(\cc,\upsilon)\|_X
\end{align*}
And

\begin{align*}
 &\Big\| \big(A_2(u,v)(\cc,\vartheta)\big) - \big(A_2(\bb u,\bb v)\big)(\cc,\vartheta) \Big\|_X \\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 & \leq 2\dfrac{\Psi_{a}^{\alpha_2}(T)}{\Gamma(\alpha_2+1)}\mathcal{K}^{*}_2(\upsilon)\|x(\cc,\upsilon) - \bb x(\cc,\upsilon)\|_X + 2\dfrac{\Psi_{a}^{\alpha_2}(T)}{\Gamma(\alpha_2+1)}\mathcal{L}^{*}_2(\upsilon)\|y(\cc,\upsilon) - \bb y(\cc,\upsilon)\|_X
\end{align*}
Thus

\[d\Big(\big(A(x,y)\big)(\cc,\upsilon) , \big(A(\bb x,\bb y)\big)(\cc,\upsilon)\Big) \leq \ww M(\vartheta)
d\Big(\big(x(\cc,\upsilon),y(\cc,\upsilon)\big) , \big(\bb x(\cc,\upsilon) ,\bb y(\cc,\upsilon)\big)\Big).\]

The Matrix $\ww M(\upsilon)$ converges to zero, then the operator $A$ is a $\ww M(\upsilon)$-contractive  operator.\\

For the following steps, we show that $B$ is completely continuous.
\begin{itemize}
\item[\textbf{step 2.}] $B(\cc,\cc,\upsilon)$ is continuous operator:

\end{itemize}
Let $(x_n , y_n)$ be a sequence such that

\[ (x_n , y_n) \rightarrow (x , y) \in X \times X \qquad \text{ as } n \rightarrow \infty\]

then, for each $\upsilon \in \Omega, t \in J$ and for $j=1,2$, we have
 
\begin{align*}
 & \Big\|  \big(B_j(x_n , y_n)\big)(t,\upsilon) - \big(B_j(x , y)\big)(t,\upsilon)\Big\| \\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 & \leq\bigg( 2\dfrac{\Psi_{a}^{\gamma_{j}+\alpha_j}(T)}{\Gamma(\gamma_{j}+\alpha_j+1)}\mathcal{M}^{*}_j(\upsilon)\|x(\cc,\upsilon) -  x_n(\cc,\upsilon)\|_X\\
 & + 2 \dfrac{\Psi_{a}^{\gamma_{j}+\alpha_j}(T)}{\Gamma(\gamma_{j}+\alpha_j+1)}\mathcal{N}^{*}_j(\upsilon)\|y(\cc,\upsilon) -  y_n(\cc,\upsilon)\|_X\bigg)\rightarrow 0\quad\textbf{ as } n\rightarrow \infty
 \end{align*}

%On the other hand, for any $\vartheta \in \Omega$ and each $t \in J$, we obtain

%\begin{align*}
% & \Big\|  \big(B_2(x_n , y_n)\big)(t,\upsilon) - \big(B_2(x , y)\big)(t,\upsilon)\Big\| \\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% & \leq\bigg( 2\dfrac{\Psi_{a}^{\gamma_{2,i}+\alpha_2}(T)}{\Gamma(\gamma_{2,i}+\alpha_2+1)}\mathcal{M}^{*}_2(\upsilon)\|x(\cc,\upsilon) -  x_n(\cc,\upsilon)\|_{\infty}\\
% & + 2 \dfrac{\Psi_{a}^{\gamma_{2,i}+\alpha_2}(T)}{\Gamma(\gamma_{2,i}+\alpha_2+1)}\mathcal{N}^{*}_2(\upsilon)\|y(\cc,\upsilon) -  y_n(\cc,\upsilon)\|_{\infty}\bigg)\rightarrow 0\quad\textbf{ as } n\rightarrow \infty
% \end{align*}

Hence, $B(\cc,\cc)(t,\upsilon)$ is continuous

\begin{itemize}
\item[\textbf{step 3.}] $ B(\cc,\cc,\upsilon) $ maps bounded sets into bounded sets in $ E \times F $. 
\end{itemize}
First, we set

\[\varphi^{*}_j(\upsilon)=\|\varphi_j(\cc,\upsilon)\|_{[a,T]},\chi^{*}_j(\upsilon)=\|\chi_j(\cc,\upsilon)\|_{[a,T]},\omega^{*}_j(\upsilon)=\|\omega_j(\cc,\upsilon)\|_{[a,T]}\] ,\[\lambda^{*}_j(\upsilon)=\sum^m_{i=1}\|\lambda_{j,i}(\cc,\upsilon)\|_{[a,T]},\rho^{*}_j(\upsilon)=\sum^m_{i=1}\|\rho_{j,i}(\cc,\upsilon)\|_{[a,T]},\mu^{*}_j(\upsilon)=\sum^m_{i=1}\|\mu_{j,i}(\cc,\upsilon)\|_{[a,T]}\]

Indeed, it is enough to show that for any $ r>0 $ there exists a positive constant $ R $ such that

\[ \|(B(x , y)))(\cc,\upsilon)\|_{X \times X} \leq R=(R_1,R_2) \]

for each $(x , y) \in B_r = \lbrace (x , y) \in X \times X : \quad \|x\|_X \leq r, \|y\|_X\leq r\rbrace$, for all $ t \in J$ and for $j=1,2$, we get
 
 \begin{align*} 
\| & \big(B_j(x,y)\big)(t,\upsilon)\|  \\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 & \leq  \sum_{i=1}^m I^{\gamma_{j,i}+\alpha_j;\psi} \big\|p_{j,i}(s,x^s(\upsilon),y^s(\upsilon),\upsilon)\big\|(t) + \sum_{i=1}^m I^{\gamma_{j,i}+\alpha_j;\psi} \big\|p_{j,i}(s,x^s(\upsilon),y^s(\upsilon),\upsilon)\big\|(t)\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 & \leq 2\dfrac{\Psi^{\gamma_{j,i}+\alpha_j}_a(T)}{\Gamma(\gamma_{j,i}+\alpha_j+1)}\bigg(\lambda^*_j(\upsilon)+\rho^*_j(\upsilon)\|x(\cc,\upsilon)\|_X + \mu^*_j(\upsilon)\|y(\cc,\upsilon)\|_X\bigg)\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 & \leq 2\dfrac{\Psi^{\gamma_{j,i}+\alpha_j}_a(T)}{\Gamma(\gamma_{j,i}+\alpha_j+1)}\bigg(\lambda^*_j(\upsilon)+r\Big(\rho^*_j(\upsilon) + \mu^*_j(\upsilon)\Big)\bigg)=R_j
\end{align*}
  
Hence

\[ \Big\| \big(B(x , y)\big)(\cc,\upsilon)\Big\|_{X \times X} = \Big\| \Big(\big(B_1(x , y)\big)(\cc,\upsilon) , \big(B_2(x , y)\big)(\cc,\upsilon)\Big)\Big\|_{X \times X}\leq (R_1 , R_2)=R .\]

\begin{itemize}
\item[\textbf{step 4.}] $B(\cc,\cc,\upsilon)$ maps bounded sets into equicontinuous sets of $ X \times X $.
\end{itemize}
 Let $ B_r$ be a bounded set of $X \times X$ as in \textbf{step.2}, let $t_1,t_2\in J$, where $t_1>t_2$, and any $(x , y)\in B_r$, $ \upsilon \in \Omega $ and for $j=1,2$, we have

\begin{align*}
& \bigg\|  \big(B_j(x , y)\big)(t_1,\upsilon) - \big(B_j(x , y)\big)(t_2,\upsilon)\bigg\|\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%       
 & \leq \bigg\|\dfrac{\Psi_a(t_1)-\Psi_a(t_2)}{\Psi_a(T)}\bigg\|\sum_{i=1}^m I^{\gamma_{j,i}+\alpha_j;\psi} \big\|p_{j,i}(s,x^s(\upsilon),y^s(\upsilon),\upsilon)\big\|(T)\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%       
 & + \sum_{i=1}^m \int^{t_1}_{t_2}  \dfrac{\psi'(s)\big(\psi(t_1) - \psi(s)\big)^{\gamma_{j,i}+\alpha_j-1}}{\Gamma(\gamma_{j,i}+\alpha_j)}  \Big\| p_{j,i}(s,x^s(\upsilon),y^s(\upsilon),\upsilon) \Big\| ds\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%       
  & + \sum_{i=1}^m \int^{t_2}_{a}  \dfrac{\psi'(s)\Big(\big(\psi(t_1) - \psi(s)\big)^{\gamma_{j,i}+\alpha_j-1} - \big(\psi(t_2) - \psi(s)\big)^{\gamma_{j,i}+\alpha_j-1}\Big)}{\Gamma(\gamma_{j,i}+\alpha_j)} \Big\| p_{j,i}(s,x^s(\upsilon),y^s(\upsilon),\upsilon) \Big\| ds \\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 & \leq \Bigg(\bigg\|\dfrac{\Psi_a(t_1)-\Psi_a(t_2)}{\Psi_a(T)}\bigg\|\dfrac{\Psi^{\gamma_{j,i}+\alpha_j}_a(T)}{\Gamma(\gamma_{j,i}+\alpha_j+1)}+\dfrac{\big(\psi(t_1) - \psi(t_2)\big)^{\gamma_{j,i}+\alpha_j}}{\Gamma(\gamma_{j,i}+\alpha_j+1)}\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 & + \bigg\|\dfrac{\big(\psi(t_1) - \psi(a)\big)^{\gamma_{j,i}+\alpha_j}}{\Gamma(\gamma_{j,i}+\alpha_j+1)} - \dfrac{\big(\psi(t_2) - \psi(a)\big)^{\gamma_{j,i}+\alpha_j}}{\Gamma(\gamma_{j,i}+\alpha_j+1)}\bigg\|\Bigg)\bigg(\lambda^*_j(\upsilon)+r\Big(\rho^*_j(\upsilon) + \mu^*_j(\upsilon)\Big)\bigg)\rightarrow 0 \text{ as } t_2 \rightarrow t_1 ,
\end{align*}

Thus the operators $ B_1 $ and $ B_2 $ are equicontinuous, and then $ B $ is also equicontinuous. Hence by the Ascoli-Arzila theorem, we deduce that $ B $ is compact. Therefore we conclude that $ B $ is completely continuous.
Now, it remains to show that the set

\[ \Lambda(\upsilon) = \big\lbrace(x , y) : \Omega \rightarrow X \in X \text{ is measurable }|\zeta(\upsilon)A(x,y) + \zeta(\upsilon) B\Big(\frac{x}{\zeta(\upsilon)} , \frac{y}{\zeta(\upsilon)} , \upsilon\Big)=(x , y)\big\rbrace\]
is bounded for some measurable mapping $\zeta : \Omega \rightarrow \mathbb{R}$ with $0<\zeta(\upsilon)<1$ on $\Omega$, let $(x , y)\in \Lambda $. Then

\begin{align*}
\big\| x(t,\upsilon) \big\|_X & 
%\leq 2\dfrac{\Psi^{\alpha_1}_a(T)}{\Gamma(\alpha_1+1)}\bigg(\varphi^*_1(\upsilon)+\chi^*_1(\upsilon)\|x(\cc,\upsilon)\|_X + \omega^*_1(\upsilon)\|y(\cc,\upsilon)\|_X\bigg)\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% & + 2\dfrac{\Psi^{\gamma_{1,i}+\alpha_1}_a(T)}{\Gamma(\gamma_{1,i}+\alpha_1+1)}\bigg(\lambda^*_1(\upsilon)+\rho^*_1(\upsilon)\|x(\cc,\upsilon)\|_X + \mu^*_1(\upsilon)\|y(\cc,\upsilon)\|_X\bigg)\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
  \leq 2\dfrac{\Psi^{\alpha_1}_a(T)}{\Gamma(\alpha_1+1)}\varphi^*_1(\upsilon)+2\dfrac{\Psi^{\gamma_{1}+\alpha_1}_a(T)}{\Gamma(\gamma_{1}+\alpha_1+1)}\lambda^*_1(\upsilon)\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
 & + 2\bigg( \dfrac{\Psi^{\alpha_1}_a(T)}{\Gamma(\alpha_1+1)}\chi^*_1(\upsilon) + \dfrac{\Psi^{\gamma_{1}+\alpha_1}_a(T)}{\Gamma(\gamma_{1}+\alpha_1+1)}\rho^*_1(\upsilon)\bigg)\|x(\cc,\upsilon)\|_X\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
 & + 2\bigg(\dfrac{\Psi^{\alpha_1}_a(T)}{\Gamma(\alpha_1+1)}\omega^*_1(\upsilon)+\dfrac{\Psi^{\gamma_{1}+\alpha_1}_a(T)}{\Gamma(\gamma_{1}+\alpha_1+1)}\mu^*_1(\upsilon)\bigg)\|y(\cc,\upsilon)\|_X
\end{align*}                        
and
\begin{align*}
\big\| y(t,\upsilon) \big\|_X & \leq 2\dfrac{\Psi^{\alpha_2}_a(T)}{\Gamma(\alpha_2+1)}\varphi^*_2(\upsilon)+2\dfrac{\Psi^{\gamma_{2}+\alpha_2}_a(T)}{\Gamma(\gamma_{2}+\alpha_2+1)}\lambda^*_2(\upsilon)\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
 & + 2\bigg( \dfrac{\Psi^{\alpha_2}_a(T)}{\Gamma(\alpha_2+1)}\chi^*_2(\upsilon) + \dfrac{\Psi^{\gamma_{2}+\alpha_2}_a(T)}{\Gamma(\gamma_{2}+\alpha_2+1)}\rho^*_2(\upsilon)\bigg)\|x(\cc,\upsilon)\|_X\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
 & + 2\bigg(\dfrac{\Psi^{\alpha_2}_a(T)}{\Gamma(\alpha_2+1)}\omega^*_2(\upsilon)+\dfrac{\Psi^{\gamma_{2}+\alpha_2}_a(T)}{\Gamma(\gamma_{2}+\alpha_2+1)}\mu^*_2(\upsilon)\bigg)\|y(\cc,\upsilon)\|_X\\ 
\end{align*}                        
by summing up the above inequalities together, we get

\begin{align*}
 & \big\| x(t,\upsilon) \big\|_X + \big\| y(t,\upsilon) \big\|_X \\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 & \leq 2\sum_{k=1}^2\bigg(\dfrac{\Psi^{\alpha_k}_a(T)}{\Gamma(\alpha_k+1)}\varphi^*_k(\upsilon)+\dfrac{\Psi^{\gamma_{k}+\alpha_k}_a(T)}{\Gamma(\gamma_{k}+\alpha_k+1)}\lambda^*_k(\upsilon)\bigg) \\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 & + 2\sum_{k=1}^2\bigg(\dfrac{\Psi^{\alpha_k}_a(T)}{\Gamma(\alpha_k+1)}\chi^*_k(\upsilon) + \dfrac{\Psi^{\gamma_{k}+\alpha_k}_a(T)}{\Gamma(\gamma_{k}+\alpha_k+1)}\rho^*_k(\upsilon) \bigg)\|x(\cc,\upsilon)\|_X\\
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 & + 2\sum_{k=1}^2\bigg(\dfrac{\Psi^{\alpha_k}_a(T)}{\Gamma(\alpha_k+1)}\omega^*_1(\upsilon)+\dfrac{\Psi^{\gamma_{k}+\alpha_k}_a(T)}{\Gamma(\gamma_{k}+\alpha_k+1)}\mu^*_k(\upsilon)\bigg)\|y(\cc,\upsilon)\|_X   
\end{align*}

Thus,
\[ \| x(\cc,\upsilon) \|_{\infty} + \| y(.,\upsilon) \| \leq \dfrac{C(\upsilon)}{1-S(\upsilon)}. \]  

where
\[C(\upsilon) =2 \sum_{k=1}^2\bigg(\dfrac{\Psi^{\alpha_k}_a(T)}{\Gamma(\alpha_k+1)}\varphi^*_k(\upsilon)+\dfrac{\Psi^{\gamma_{k}+\alpha_k}_a(T)}{\Gamma(\gamma_{k}+\alpha_k+1)}\lambda^*_k(\upsilon)\bigg),\]
                         %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{align*}
S(\upsilon) = 2\max \Bigg \lbrace & \bigg(\sum_{k=1}^2\bigg(\dfrac{\Psi^{\alpha_k}_a(T)}{\Gamma(\alpha_k+1)}\chi^*_k(\upsilon) + \dfrac{\Psi^{\gamma_{k}+\alpha_k}_a(T)}{\Gamma(\gamma_{k}+\alpha_k+1)}\rho^*_k(\upsilon) \bigg) \\
 & , \sum_{k=1}^2\bigg(\dfrac{\Psi^{\alpha_k}_a(T)}{\Gamma(\alpha_k+1)}\omega^*_1(\upsilon)+\dfrac{\Psi^{\gamma_{k}+\alpha_k}_a(T)}{\Gamma(\gamma_{k}+\alpha_k+1)}\mu^*_k(\upsilon)\bigg) \Bigg \rbrace \neq 1.
\end{align*}
This shows that set $\Lambda(\upsilon)$ is bounded, consequently of steps (1-3) and theorem \ref{Banach}, we conclude the operator $ T $ has at least one random fixed point, which is a solution of the system \eqref{pb}-\eqref{cond}.
\end{proof}
\section{An example}

\begin{thebibliography}{9999}
\bibitem{Abbas.A} S. Abbas, N. Al Arifi, M. Benchohra, Y. Zhou;  \emph{Random coupled Hilfer and Hadamard fractional differential systems in generalized Banach spaces.}. Mathematics 2019, 7, 285.
    
\bibitem{Abbas.B} S. Abbas, M. Benchohra, Y. Zhou: \emph{ Coupled Hilfer fractional differential systems with random effects}. Adv. Differ. Equ. 2018, 369(2018).


\bibitem{Almeida} R. Almeida; \emph{ A Caputo fractional derivative of a function with respect to another function }. Commun. Nonlinear Sci. Numer. Simulat. 44 (2017) 460-481.

\bibitem{Bharucha} A. T. Bharucha-Reid; \emph{ Random Integral Equations}. Academic Press, New York, 1972.

\bibitem{Carvalho} A. Carvalho, C. M. A. Pinto; \emph{ A delay fractional order model for the co-infection of malaria and HIV/AIDS}. Int. J. Dyn. Control 5, 168–186 (2017).



\bibitem{Graef} J. R. Graef and A. Petrusel; \emph{ Some Krasnosel'skii type random fixed point theorem}. J. Nonlinear Funct. Anal. 2017, 2017, 46.



\bibitem{Petre} I. R. Petre, A. Pertrusel; \emph{ Kranoselskii's theorem in generalized Banach spaces and applications}. Electron. J. Qual. Theory Differ. Equ. 2012, 85, 1-20.

\bibitem{Sincer} M. L. Sincer, J. J. Nieto, A. Ouahab; \emph{Random fixed point theorems in generalized Banach spaces and applications}. Random Oper. Stoch. Equ. 2016, 24, 93-112.

\bibitem{Varga} R. S. Varga; \emph{Matrix Iterative Analysis }. Second Revised and Expanded, Springer Series in Computational Mathematics 27, Springer: Berlin, Germany, 2000.

 
\end{thebibliography}

\end{document}