#  __
Quasi-Gramian Solution of a Noncommutative Extension of the Higher-Order
Nonlinear Schrödinger Equation

__
H. W. A. Riaz1 and J. Lin1
1 department of physics, zhejiang normal university, jinhua 321004, pr china

###### Abstract.

The nonlinear Schrödinger (NLS) equation, which incorporates higher-order
dispersive terms, is widely employed in the theoretical analysis of various
physical phenomena. In this study, we explore the non-commutative extension of
the higher-order NLS equation (HNLS). We treat real or complex-valued
functions, such as $g_{1}=g_{1}(x,\;t)$ and $g_{2}=g_{2}(x,\;t)$, as non-
commutative, and employ the Lax pair associated with the evolution equation as
in the commutation case. We derive the quasi-Gramian solution of the system by
employing a binary Darboux transformation (DT). Moreover, the solution can be
used to study the stability of plane waves and to understand the generation of
periodic patterns in the context of modulational instability.

<EMAIL_ADDRESS><EMAIL_ADDRESS>
##  1\. Introduction

The nonlinear Schrödinger (NLS) equation, which incorporates higher-order
dispersive terms, is widely used in the theoretical analysis of various
physical phenomena, including nonlinear optics, molecular systems, and fluid
dynamics [1, 2, 3]. With the addition of fourth-order terms, known as the
Lakshmanan-Porsezian-Daniel (LPD) equation, it describes higher-order
molecular excitations with quadruple-quadruple coefficients and possesses
integrability [4, 5, 6]. Lakshmanan et al. investigated its application to
study nonlinear spin excitations involving bilinear and biquadratic
interactions [5]. In recent years, Ankiewicz et al. introduced a further
extension of the NLSE by incorporating third-order (odd) and fourth-order
(even) dispersion terms [7]. The integrability of this extended NLSE, with
certain parameter values, was confirmed in Ref. [8], where Lax operators were
introduced. We now write this equation appears in the mentioned references
with some modification as

$\displaystyle\mathrm{i}u_{t}+\alpha_{2}\left(u_{xx}+2u|u|^{2}\right)+\mathrm{i}\alpha_{1}\left(u_{xxx}+6u_{x}|u|^{2}\right)+\gamma(u_{xxxx}+6\bar{u}u_{x}^{2}+4u|u_{x}|^{2}$
(1.1) $\displaystyle+8|u|^{2}u_{xx}+2u^{2}\bar{u}_{xx}+6u|u|^{4})=0.$

where $u=u(x,\;t)$ is a complex-valued scalar function, and $\bar{u}$
represents its complex conjugate. This equation includes several particular
cases, such as the standard nonlinear Schrödinger equation (NLSE) with
$\alpha_{1}=\gamma=0$ [11], the Hirota equation with $\gamma=0$ [12], and the
Lakshmanan-Porsezian-Daniel equation with $\alpha_{1}=0$ [5].

In this study, we explore the non-commutative extension of the higher-order
NLS (HNLS) equation (1.1). Noncommutative integrable systems have received
considerable attention for their relevance in quantum field theories,
$D$-brane dynamics, and string theories [13, 14, 15]. Non-commutativity often
arises from phase-space quantization, introducing non-commutativity among
independent variables through a star product [16, 17]. Our approach to
inducing non-commutativity in a given nonlinear evolution equation parallels
the methods employed by Lechtenfeld et al. [18], Gilson and Nimmo [19], and
Gilson and Macfarlane [20] for the non-commutative generalization of the sine-
Gordon, Kadomtsev-Petviashvili (KP), and Davey-Stewartson (DS) equations,
respectively.

We adopt a systematic method to extend the chosen equation to its non-
commutative form, without explicitly specifying the nature of non-
commutativity. We consider real or complex-valued functions, such as
$g_{1}=g_{1}(x,\;t)$ and $g_{2}=g_{2}(x,\;t)$, as non-commutative and take
advantage of the same Lax pair as in the commutative scenario to describe the
equation of non-linear evolution.

In this paper, we investigate a non-commutative (nc) version of the HNLS
equation. We define the Lax pair for the nc-HNLS equation in this context. To
find solutions to the nc-HNLS equation, we construct the Darboux matrix and
the binary Darboux matrix. We present explicit quasi-Gramian solutions for the
non-commutative fields of the nc-HNLS equation, which, after reducing the non-
commutativity limit, can be reduced to a ratio of Gramian solutions.

##  2\. Modulation instability

For analyzing the modulation instability, we give a plane-wave solution of the
system (1.1)

$u(x,t)=ce^{i(6c^{4}\gamma+2\alpha_{2}c^{2})t},$ (2.1)

The solution provided by equation (2.1) holds significant importance in the
realm of optics, particularly within the context of the nonlinear Schrödinger
equation. This solution represents a wave with constant amplitude that
undergoes a non-linear evolution over time. The dynamics are determined by
parameters such as the amplitude $c$, the constant $\gamma$, and $\alpha_{2}$.
This solution’s application extends to studying the stability of plane waves
and comprehending the generation of periodic patterns through modulational
instability. It serves as a prime example of how the nonlinear Schrödinger
equation can lead to complex behavior in optical systems, thereby making it a
crucial area of research in this field.

An approach to assess the stability of the plane wave solution involves
introducing perturbations to the solution and examining the linearized
evolution of these perturbations. To simplify the analysis, the common phase
can be factored out of the equation. This leads to a first-order ordinary
differential equation (ODE) that couples the complex field with its complex
conjugate as a result of the perturbation. Substituting the perturbed function
$v(x,\;t)$ into equation (2.1), we obtain

$u(x,t)=(c+v(x,\;t))e^{i(6c^{4}\gamma+2\alpha_{2}c^{2})t},$ (2.2)

Substituting (2.2) into (1.1) and after linearization, we have

$iv_{t}+\alpha_{2}v_{xx}+i\alpha_{1}(v_{xxx}+6c^{2}v_{x})+\gamma
v_{xxxx}+2\gamma c^{2}(4v_{xx}+\bar{v}_{xx})+c^{2}(12\gamma
c^{2}+2\alpha_{2})(v+\bar{v})=0,$ (2.3)

In order to analyze the stability of the plane wave solution, the Fourier
transform of the equation is taken. This results in a first-order ODE that
governs the real and imaginary parts of evolution. The stability of the
solution can then be determined by looking at the eigenvalues of this ODE. In
particular, the eigenvalues represent the exponent in the time evolution of
the solution. Thus, the Fourier transform of the evolution equation (2.3) is

$i\hat{v}_{t}-\frac{\alpha_{2}}{2}k^{2}\hat{v}-\alpha_{1}k(-k^{2}+6c^{2})\hat{v}+\gamma
k^{4}\hat{v}-2\gamma c^{2}k^{2}(4\hat{v}+\bar{\hat{v}})+c^{2}(12\gamma
c^{2}+2\alpha_{2})(\hat{v}+\bar{\hat{v}})=0,$ (2.4)

The linear evolution equation for $\hat{v}$ can be evaluated by separating it
into its real and imaginary components. Thus, for $\hat{v}=v_{1}+iv_{2}$, we
have system of differential equation

$y_{t}=\left[\begin{array}[]{cc}0&\beta-6\gamma\,c^{2}k^{2}\\\
-\beta-10\gamma\,c^{2}k^{2}+2c^{2}\left(12c^{2}\gamma+2\alpha_{2}\right)&0\end{array}\right]y,$
(2.5)

where $\beta=\frac{\alpha
2\,k^{2}}{2}+\alpha_{1}k\left(6c^{2}-k^{2}\right)-\gamma\,k^{4}$, and
$y=[\begin{array}[]{cc}v_{1}&v_{2}\end{array}]^{T}$. One can evaluate the
stability of a system by analyzing exponential solutions in the form of $y=\nu
e^{\omega t}$, which leads to an eigenvalue problem where the eigenvalues are
then given by

$\displaystyle\omega(k)={\frac{|k|}{2}\sqrt{\left(\alpha_{2}\left(-8c^{2}+k^{2}\right)-2\alpha_{1}k\left(-6c^{2}+k^{2}\right)-2\gamma\left(24c^{4}-10c^{2}k^{2}+k^{4}\right)\right)\beta_{1}}},$
(2.6)
$\displaystyle\beta_{1}=k\alpha_{2}-2\alpha_{1}\left(-6c^{2}+k^{2}\right)-2\gamma
k\left(6c^{2}+k^{2}\right).$

The plot of the equation (2.6) is shown in Fig. 1

Figure 1.

The stability of the solution becomes evident when examining a graphical
representation of the real parts of eigenvalues plotted against different
frequencies. If the real part is positive, the solution will exhibit growth;
conversely, negative values indicate decay. The overall stability is
established by observing whether the real eigenvalue remains positive or
negative across various frequency ranges. This phenomenon is referred to as
modulational instability.

##  3\. Non-commutative HNLS equation

The spectral problem associated with the nc-HNLS equation is given by

$\displaystyle\Gamma=\partial_{x}-\lambda\mathcal{J}-\mathcal{U},$ (3.1)
$\displaystyle\Delta=\partial_{t}-\mathcal{B}-\mathcal{V}_{p},$ (3.2)

where

$\displaystyle\mathcal{J}=\left[\begin{array}[]{cc}\mathrm{-I}&0\\\
0&\mathrm{I}\end{array}\right],\quad\mathcal{U}=\left[\begin{array}[]{cc}0&{u}\\\
-{\mathit{u^{{\dagger}}}}&0\end{array}\right],\quad\mathcal{V}_{p}=\gamma\left[\begin{array}[]{cc}\mathcal{V}_{1}&\mathcal{V}_{2}\\\
\mathcal{V}_{3}&\mathcal{V}_{4}\end{array}\right],$ (3.9)
$\displaystyle\mathcal{B}=\small{\begin{pmatrix}\rho_{1}+\rho_{2}{u}{\mathit{u^{{\dagger}}}}-\alpha_{1}\left(\left(\frac{\partial}{\partial
x}{u}\right){\mathit{u^{{\dagger}}}}-{u}\left(\frac{\partial}{\partial
x}{\mathit{u^{{\dagger}}}}\right)\right)&\mathcal{A}_{2}\\\
\mathcal{A}_{1}&-\rho_{1}-\rho_{2}{\mathit{u^{{\dagger}}}}{u}-\alpha_{1}\left(-{\mathit{u^{{\dagger}}}}\left(\frac{\partial}{\partial
x}{u}\right)+\left(\frac{\partial}{\partial
x}{\mathit{u^{{\dagger}}}}\right){u}\right)\end{pmatrix}},$ (3.10)
$\displaystyle\mathcal{A}_{1}=\small{-4\lambda^{2}\alpha_{1}{\mathit{u^{{\dagger}}}}+2\lambda\left(-\alpha_{2}{\mathit{u^{{\dagger}}}}+\mathrm{I}\alpha_{1}\left(\frac{\partial}{\partial
x}{\mathit{u^{{\dagger}}}}\right)\right)+\mathrm{I}\alpha_{2}\left(\frac{\partial}{\partial
x}{\mathit{u^{{\dagger}}}}\right)-\alpha_{1}\left(-2{\mathit{u^{{\dagger}}}}{u}{\mathit{u^{{\dagger}}}}-\frac{\partial^{2}}{\partial
x^{2}}{\mathit{u^{{\dagger}}}}\right)},$
$\displaystyle\mathcal{A}_{2}=4\lambda^{2}\alpha_{1}{u}+2\lambda\left(\alpha_{2}{u}+\mathrm{I}\alpha_{1}\left(\frac{\partial}{\partial
x}{u}\right)\right)+\mathrm{I}\alpha_{2}\left(\frac{\partial}{\partial
x}{u}\right)-\alpha_{1}\left(2{u}{\mathit{u^{{\dagger}}}}{u}+\frac{\partial^{2}}{\partial
x^{2}}{u}\right),$
$\displaystyle\mathcal{V}_{1}=\mathrm{I}\left(\left(\frac{\partial^{2}}{\partial
x^{2}}{u}\right){\mathit{u^{{\dagger}}}}+{u}\left(\frac{\partial^{2}}{\partial
x^{2}}{\mathit{u^{{\dagger}}}}\right)-\left(\frac{\partial}{\partial
x}{u}\right)\left(\frac{\partial}{\partial
x}{\mathit{u^{{\dagger}}}}\right)\right)-2\lambda\left({u}\left(\frac{\partial}{\partial
x}{\mathit{u^{{\dagger}}}}\right)-\left(\frac{\partial}{\partial
x}{u}\right){\mathit{u^{{\dagger}}}}\right),$
$\displaystyle\mathcal{V}_{2}=-4\,\mathrm{I}\lambda^{2}\left(\frac{\partial}{\partial
x}{u}\right)+3\,\mathrm{I}\left({u}{\mathit{u^{{\dagger}}}}\left(\frac{\partial}{\partial
x}{u}\right)+\left(\frac{\partial}{\partial
x}{u}\right){\mathit{u^{{\dagger}}}}{u}\right)+\mathrm{I}\left(\frac{\partial^{3}}{\partial
x^{3}}{u}\right)+2\lambda\left(\frac{\partial^{2}}{\partial x^{2}}{u}\right),$
$\displaystyle\mathcal{V}_{3}=-4\,\mathrm{I}\lambda^{2}\left(\frac{\partial}{\partial
x}{\mathit{u^{{\dagger}}}}\right)+3\,\mathrm{I}\left({\mathit{u^{{\dagger}}}}{u}\left(\frac{\partial}{\partial
x}{\mathit{u^{{\dagger}}}}\right)+\left(\frac{\partial}{\partial
x}{\mathit{u^{{\dagger}}}}\right){u}{\mathit{u^{{\dagger}}}}\right)+\mathrm{I}\left(\frac{\partial^{3}}{\partial
x^{3}}{\mathit{u^{{\dagger}}}}\right)-2\lambda\left(\frac{\partial^{2}}{\partial
x^{2}}{\mathit{u^{{\dagger}}}}\right),$
$\displaystyle\mathcal{V}_{4}=-\mathrm{I}\left(\left(\frac{\partial^{2}}{\partial
x^{2}}{\mathit{u^{{\dagger}}}}\right){u}+{\mathit{u^{{\dagger}}}}\left(\frac{\partial^{2}}{\partial
x^{2}}{u}\right)-\left(\frac{\partial}{\partial
x}{\mathit{u^{{\dagger}}}}\right)\left(\frac{\partial}{\partial
x}{u}\right)\right)+2\lambda\left(-{\mathit{u^{{\dagger}}}}\left(\frac{\partial}{\partial
x}{u}\right)+\left(\frac{\partial}{\partial
x}{\mathit{u^{{\dagger}}}}\right){u}\right),$

where
$\rho_{1}=-4\,\mathrm{I}\alpha_{1}\lambda^{3}-2\,\mathrm{I}\lambda^{2}\alpha_{2},\;\rho_{2}=\,\mathrm{I}(2\lambda\alpha_{1}+\alpha_{2})$,
$\rho_{3}=3\,\mathrm{I}{u}{\mathit{u^{{\dagger}}}}{u}{\mathit{u^{{\dagger}}}}-4\,\mathrm{I}\lambda^{2}{u}{\mathit{u^{{\dagger}}}}+8\,\mathrm{I}\lambda^{4},\;\rho_{4}=-3\,\mathrm{I}{\mathit{u^{{\dagger}}}}{u}{\mathit{u^{{\dagger}}}}{u}+4\,\mathrm{I}\lambda^{2}{\mathit{u^{{\dagger}}}}{u}-8\,\mathrm{I}\lambda^{4}$,
$\rho_{5}=-8\lambda^{3}{u}+4\lambda{u}{\mathit{u^{{\dagger}}}}{u},\;\rho_{6}=8\lambda^{3}{\mathit{u^{{\dagger}}}}-4\lambda{\mathit{u^{{\dagger}}}}{u}{\mathit{u^{{\dagger}}}}$.
The equation of motion for the system can be derived by setting the commutator
of $\Gamma$ and $\Delta$ equal to zero and equating the coefficients at
$\lambda$

$\displaystyle
i{u}_{t}+i\left({u}_{xxx}+3({u}_{x}{\mathit{u^{\dagger}}}{u}+{u}{\mathit{u^{\dagger}}}{u}_{x})\right)\alpha_{1}+\left(2{u}{\mathit{u^{\dagger}}}{u}+{u}_{xx}\right)\alpha_{2}+({u}_{xxxx}$
(3.11)
$\displaystyle+2({u}_{x}{\mathit{u^{\dagger}}}_{x}{u}+{u}{\mathit{u^{\dagger}}}_{x}{u}_{x}+{u}{\mathit{u^{\dagger}}}_{xx}{u})+4({u}_{xx}{\mathit{u^{\dagger}}}{u}+{u}{\mathit{u^{\dagger}}}{u}_{xx})+6({u}_{x}{\mathit{u^{\dagger}}}{u}_{x}+{u}{\mathit{u^{\dagger}}}{u}{\mathit{u^{\dagger}}}{u}))\gamma=0,$

where $u=u(x,\;t)$ is an nc object, † denotes the adjoint (Hermitian
conjugate), $\alpha_{1},\;\alpha_{2}$ and $\gamma$ are real parameters, and
$\lambda$ is a spectral parameter (real or complex). The equation presented in
(3.11) is a non-commutative generalization of the HNLS equation, as given in
(1.1). This equation exhibits several interesting properties. For instance,
when both $\alpha_{2}$ and $\gamma$ are zero, it reduces to an nc
generalization of the complex modified Korteweg–de Vries (KdV) equation and to
the standard modified KdV equation for real-valued $u$. Moreover, setting
$\gamma$ to zero results in an nc generalization of the Hirota equation, while
setting $\alpha_{1}$ and $\alpha_{2}$ to zero simultaneously yields the well-
known Lakshmanan–Porsezian–Daniel (LPD) equation. Finally, when $\alpha_{1}$
and $\gamma$ are set to zero, the equation reduces to an nc generalization of
the nonlinear Schrödinger (NLS) equation. After relaxing the noncommutativity
condition, Eq. (3.11) corresponds to the commutative counterpart. The spectral
problem linked with (1.1) remains the same as that of (3.1) and (3.2), with
the exception that $u$ and $\bar{u}$ are now perceived as commutative
functions.

###  3.1. Quasideterminants

In non-commutative algebra, quasideterminants serve as a replacement for
ordinary determinants of matrices. They hold a similar significance in non-
commutative algebra as ordinary determinants do in commutative algebra and
have found vast applications in the domain of non-commutative integrable
systems [19, 21, 22].

The quasideterminant $|M|_{ij}$ for $i,\;j=1,\;...,\;n$ of an $n\times n$
matrix over a non-commutative ring R, expanded about the matrix $m_{ij}$, is
defined as

$\left|M\right|_{ij}=\left|\begin{array}[]{cc}M^{ij}&c^{i}_{j}\\\
r^{j}_{i}&\framebox{$m_{ij}$}\end{array}\right|=m_{ij}-r^{j}_{i}\left(M^{ij}\right)^{-1}c^{i}_{j},$
(3.12)

where $m_{ij}$ is referred to as the expansion point and represents the
$ij$-th entry of $M$, $r^{j}_{i}$ denotes the $i$-th row of $M$ without the
$j$-th entry, $c^{i}_{j}$ represents the $j$-th column of $M$ without the
$i$-th row, and $M^{ij}$ is the submatrix of $M$ obtained by removing the
$i$-th row and the $j$-th column from $M$.

Quasideterminants are not merely a generalization of usual commutative
determinants but are also related to inverse matrices. The inverse of a matrix
$M=\left(\begin{smallmatrix}m_{11}&m_{12}\\\
m_{21}&m_{22}\end{smallmatrix}\right)$ is defined as

$\displaystyle
M^{-1}=\left(|M|_{ji}^{-1}\right)=\left(\begin{array}[]{cc}\left|\begin{array}[]{cc}\framebox{$m_{11}$}&m_{12}\\\
m_{21}&m_{22}\end{array}\right|^{-1}&\left|\begin{array}[]{cc}{m_{11}}&\framebox{$m_{12}$}\\\
{m_{21}}&m_{22}\end{array}\right|^{-1}\\\
\left|\begin{array}[]{cc}{m_{11}}&m_{12}\\\
\framebox{$m_{21}$}&m_{22}\end{array}\right|^{-1}&\left|\begin{array}[]{cc}{m_{11}}&m_{12}\\\
m_{21}&\framebox{$m_{22}$}\end{array}\right|^{-1}\end{array}\right).$ (3.23)

##  4\. Darboux transformation

In this section, a Darboux transformation is introduced for the system of nc-
HNLS equation (3.11) through the definition of the Darboux matrix

$\displaystyle D(\mathcal{Y})=\lambda I-\mathcal{Y}\Lambda\mathcal{Y}^{-1},$
(4.1)

and the Lax operators $\Gamma$ and $\Delta$

$\Gamma=\partial_{x}-\lambda\mathcal{J}-\mathcal{U},\quad\Delta=\partial_{t}-\mathcal{B}-\mathcal{V}_{p}.$
(4.2)

The spectral parameter $\lambda$, which can be real or complex, is
incorporated along with the constant $q\times q$ matrix $\Lambda$, and the nc
objects $\mathcal{U},\;\mathcal{V}_{p}$, and $\mathcal{B}$ from (3.9) and
(3.10) respectively, are utilized as entries in the Lax operators. Now,
consider a function $\varphi=\varphi(x,\;t)$ that is an eigenfunction of the
Lax operators $\Gamma$ and $\Delta$ such that $\Gamma(\varphi)=0$ and
$\Delta(\varphi)=0$. We can define a new function $\tilde{\varphi}$ using the
Darboux matrix $D(\mathcal{Y})$:

$\displaystyle\tilde{\varphi}$ $\displaystyle:=$ $\displaystyle
D_{\mathcal{Y}}(\varphi),$ (4.5) $\displaystyle=$
$\displaystyle\lambda\varphi-\mathcal{Y}\Lambda\mathcal{Y}^{-1}\phi=\left|\begin{array}[]{cc}\mathcal{Y}&\varphi\\\
\mathcal{Y}\Lambda&\framebox{$\lambda\varphi$}\end{array}\right|$

The function $\tilde{\varphi}$ is a generic function of the new operators
$\tilde{\Gamma}=D_{\mathcal{Y}}\Gamma D_{\mathcal{Y}}^{-1}$ and
$\tilde{\Delta}=D_{\mathcal{Y}}\Delta D_{\mathcal{Y}}^{-1}$. Next, we define
$\mathcal{Y}_{[1]}=\mathcal{Y}_{1}$ and $\varphi_{[1]}=\varphi$ as generic
eigenfunctions of $\Gamma_{[1]}=\Gamma$ and $\Delta_{[1]}=\Delta$, where
$\mathcal{Y}_{1}$ and $\varphi_{1}$ are $2\times 2$ matrices. We then define
$\varphi_{[2]}=D_{\mathcal{Y}_{[1]}}\left(\varphi_{[1]}\right)$ and
$\mathcal{Y}_{[2]}=\varphi_{[2]}|_{\varphi\rightarrow\mathcal{Y}_{2}}$ to be
eigenfunctions of the new operators
${\Gamma}_{[2]}=D_{\mathcal{Y}_{[1]}}{\Gamma}_{[1]}D_{\mathcal{Y}_{[1]}}^{-1}$
and
${\Delta}_{[2]}=D_{\mathcal{Y}_{[1]}}{\Delta}_{[1]}D_{\mathcal{Y}_{[1]}}^{-1}$.

We define $\mathcal{Y}_{i},\;i=1,\;...,\;n$ (where each $\mathcal{Y}_{i}$ is a
matrix of size $2\times 2$) be the set of particular eigenfunctions of
$\Gamma$ and $\Delta$. In the non-commutative case, we consider each entry of
$\mathcal{Y}_{i}$ as a matrix. We define the $(n+1)$th eigenfunction as
$\varphi_{[n+1]}=D_{\mathcal{Y}_{[n]}}\left(\varphi_{[n]}\right)$, which is a
generic eigenfunction of the new operators
${\Gamma}_{[n+1]}=D_{\mathcal{Y}_{[n]}}{\Gamma}_{[n]}D_{\mathcal{Y}_{[n]}}^{-1}$,
and
${\Delta}_{[n+1]}=D_{\mathcal{Y}_{[n]}}{\Delta}_{[n]}D_{\mathcal{Y}_{[n]}}^{-1}$.
Here, $\mathcal{Y}_{[n]}$ is a DT from $\Gamma_{[n]}$ to $\Gamma_{[n+1]}$ and
$\Delta_{[n]}$ to $\Delta_{[n+1]}$. Thus, the $n$th-order DT is given by
$\varphi_{[n+1]}=D_{\mathcal{Y}_{[n]}}\left(\varphi_{[n]}\right)=\lambda\varphi_{[n]}-\mathcal{Y}_{[n]}\Lambda_{n}\mathcal{Y}_{[n]}^{-1}\varphi_{[n]}$,
where $\mathcal{Y}_{[j]}=\varphi_{[j]}|_{\varphi\rightarrow\mathcal{Y}_{j}}$.
Let us define a matrix $\Xi=(\mathcal{Y}_{1},\;...,\;\mathcal{Y}_{n})$
comprising of the eigenfunctions $\mathcal{Y}_{i},\;i=1,\;...,\;n$, and set
$\varphi_{[1]}=\varphi$. Then, we can represent the $n$th iteration of the DT
in quasideterminant form as follows:

$\varphi_{[n+1]}=\left|\begin{array}[]{cc}\Xi&\varphi\\\ \vdots&\vdots\\\
\Xi^{(n-1)}&\varphi^{(n-1)}\\\
\Xi^{(n)}&\framebox{$\varphi^{(n)}$}\end{array}\right|,$ (4.6)

Here, $\varphi^{(n)}=\lambda^{n}\varphi$ and $\Xi^{(n)}=\Xi{\Lambda}^{n}$.,
where each $\Lambda^{i},\;i=1,\;...,\;n$, is a constant matrix. Hence, we
expressed a quasideterminant formula for $\varphi_{n+1}$ in terms of the known
eigenfunctions $\mathcal{Y}_{i},\;i=1,\;...,\;n$ and the eigenfunction
$\varphi$ of the “seed” Lax pair $\Gamma=\Gamma_{1}$, $\Delta=\Delta_{1}$.

##  5\. Quasi-Wronskian solutions

In the upcoming analysis, we will examine how the DT $D_{\mathcal{Y}}=\lambda
I-\mathcal{Y}\Lambda\mathcal{Y}^{-1}$ affects the Lax operator
$\Gamma=\Gamma_{1}$, where $\mathcal{Y}$ is an eigenfunction of $\Gamma$
(since $\Gamma(\mathcal{Y})=0$ by definition) and $\Lambda$ is an eigenvalue
matrix. It is important to note that the same results apply to the operator
$\Delta=\Delta_{1}$. As a result of this transformation, the operator $\Gamma$
is converted to a new operator $\tilde{\Gamma}=\Gamma_{[2]}$, which can be
expressed as $\tilde{\Gamma}=D_{\mathcal{Y}}{\Gamma}D_{\mathcal{Y}}^{-1}$. By
substituting (3.1) and (4.1) into the latter equation and equating the
coefficients at $\lambda^{j}$, we obtain two equations,

$\displaystyle\mathcal{U}_{[2]}-\mathcal{U}-[\mathcal{J},\;\mathcal{Y}\Lambda\mathcal{Y}^{-1}]=0,$
(5.1)

and

$-\partial_{x}(\mathcal{Y}\Lambda\mathcal{Y}^{-1})+[\mathcal{U},\;\mathcal{Y}\Lambda\mathcal{Y}^{-1}]+[\mathcal{J},\;\mathcal{Y}\Lambda\mathcal{Y}^{-1}]\mathcal{Y}\Lambda\mathcal{Y}^{-1}=0.$
(5.2)

To confirm the validity of equation (5.2), we express equation (3.1) using a
particular eigenfunction $\mathcal{Y}$ as
$\mathcal{Y}_{x}=\mathcal{J}\mathcal{Y}\Lambda+\mathcal{U}\mathcal{Y}$. By
utilizing this equation, we can easily check that the condition expressed in
(5.2) is satisfied. To simplify the notation, a matrix $\mathcal{F}$ is
introduced such that $\mathcal{U}=[\mathcal{F},\;\mathcal{J}]$. This equation
is satisfied if
$\mathcal{F}=\frac{1}{2\dot{\imath}}\left(\begin{smallmatrix}0&{u}\\\
{u}^{{\dagger}}&0\end{smallmatrix}\right)$. Then, equation (5.1) with
$\mathcal{U}=[\mathcal{F},\;\mathcal{J}]$ can be used to obtain
$\mathcal{F}_{[2]}=\mathcal{F}-\mathcal{Y}^{(1)}\mathcal{Y}^{-1}$, where
$\mathcal{Y}^{(1)}$ is defined as $\mathcal{Y}\Lambda$. After $n$ repeated
applications of the DT $D_{\mathcal{Y}}$, we have

$\displaystyle\mathcal{F}_{[n+1]}=\mathcal{F}_{[n]}-\mathcal{Y}^{(1)}_{[n]}\mathcal{Y}^{-1}_{[n]}=\mathcal{F}-\sum_{k=1}^{n}\mathcal{Y}^{(1)}_{[k]}\mathcal{Y}^{-1}_{[k]},$
(5.3)

where
$\mathcal{F}_{[1]}=\mathcal{F},\;\mathcal{Y}_{[1]}=\mathcal{Y}_{1}=\mathcal{Y}$
and $\Lambda_{1}=\Lambda$.

Because our nc-HNLS equation (3.11) is expressed in terms of $u$ and
$u^{{\dagger}}$, it is more appropriate to express the quasi-Wronskian
solution in terms of these objects. For this, we express each
$\mathcal{Y}_{i},\;i=1,\;...,\;n$ as a $2\times 2$ matrix as
$\mathcal{Y}_{i}=\left(\begin{smallmatrix}\varphi_{2i-1}&\varphi_{2i}\\\
\chi_{2i-1}&\chi_{2i}\end{smallmatrix}\right)$. For $\varphi=\varphi(x,\;t)$
and $\chi=\chi(x,\;t)$, we can express $\mathcal{F}_{[n+1]}$ as

$\mathcal{F}_{[n+1]}=\mathcal{F}+\left(\begin{array}[]{cc}\left|\begin{array}[]{cc}\widehat{\Xi}&f_{2n-1}\\\
\varphi^{(n)}&\framebox{$0$}\end{array}\right|&\left|\begin{array}[]{cc}\widehat{\Xi}&f_{2n}\\\
\varphi^{(n)}&\framebox{$0$}\end{array}\right|\\\ \\\
\left|\begin{array}[]{cc}\widehat{\Xi}&f_{2n-1}\\\
\chi^{(n)}&\framebox{$0$}\end{array}\right|&\left|\begin{array}[]{cc}\widehat{\Xi}&f_{2n}\\\
\chi^{(n)}&\framebox{$0$}\end{array}\right|\end{array}\right),$ (5.4)

where $\widehat{\Xi}=(\mathcal{Y}_{j}^{(i-1)})_{i,\;j=1,\;...,\;n}$ is the
$2n\times 2n$ matrix of $\mathcal{Y}_{1},\;...,\;\mathcal{Y}_{n}$, and
$f_{2n-1}$ and $f_{2n}$ are the column vectors $2n\times 1$ with a 1 in the
$(2n-1)$th and $(2n)$th row, respectively, and zeros elsewhere, while
$\varphi^{(n)}=\left(\varphi_{1}^{(n)},\;...,\;\varphi_{2n}^{(n)}\right),\;\chi^{(n)}=\left(\chi_{1}^{(n)},\;...,\;\chi_{2n}^{(n)}\right)$
denote the $1\times 2n$ row vectors. We can express the quasi-Wronskian
solutions for $u$ and $u^{{\dagger}}$ by utilizing
$\mathcal{F}=\frac{1}{2\dot{\imath}}\left(\begin{smallmatrix}0&u\\\
u^{\dagger}&0\\\ \end{smallmatrix}\right)$ which leads to:

$u_{[n+1]}=u+2\dot{\imath}\left|\begin{array}[]{cc}\widehat{\Xi}&f_{2n}\\\
\varphi^{(n)}&\framebox{$0$}\end{array}\right|,\qquad
u_{[n+1]}^{\dagger}=u^{\dagger}+2\dot{\imath}\left|\begin{array}[]{cc}\widehat{\Xi}&f_{2n-1}\\\
\chi^{(n)}&\framebox{$0$}\end{array}\right|.$ (5.5)

We proceed in the next section to construct the binary DT for the nc-HNLS
equation, using a strategy similar to that employed in [9].

##  6\. Binary Darboux transformation

We introduce $\mathcal{Y}_{1},\;...,\;\mathcal{Y}_{n}$ as eigenfunctions of
the Lax operators $\Gamma$ and $\Delta$, and
$\mathcal{Z}_{1},\;...,\;\mathcal{Z}_{n}$ as eigenfunctions of the adjoint Lax
operators $\Gamma^{{\dagger}}$ and $\Delta^{\dagger}$. Assuming
$\varphi_{[1]}=\varphi$ to be a generic eigenfunction of the Lax operators
$\Gamma$ and $\Delta$, and $\psi_{[1]}=\psi$ to be a generic eigenfunction of
the adjoint Lax operators $\Gamma^{{\dagger}}$ and $\Delta^{{\dagger}}$, we
define the binary DT and its adjoint as

$D_{\mathcal{Y}_{[1]},\;\mathcal{Z}_{[1]}}=I-\mathcal{Y}_{[1]}\Upsilon(\mathcal{Y}_{[1]},\;\mathcal{Z}_{[1]})^{-1}\Omega^{-{\dagger}}\mathcal{Z}_{[1]}^{{\dagger}},$
(6.1)
$D_{\mathcal{Y}_{[1]},\;\mathcal{Z}_{[1]}}^{{\dagger}}=I-\mathcal{Z}_{[1]}\Upsilon(\mathcal{Y}_{[1]},\;\mathcal{Z}_{[1]})^{-{\dagger}}\Lambda^{-{\dagger}}\mathcal{Y}_{[1]}^{{\dagger}}.$
(6.2)

In the context of the binary DT, we use $\mathcal{Y}_{[1]}=\mathcal{Y}_{1}$ as
the initial eigenfunction that characterizes the transformation from the Lax
operators $\Gamma$ and $\Delta$ to the new operators $\tilde{\Gamma}$ and
$\tilde{\Delta}$. Similarly, we define $\mathcal{Z}_{[1]}=\mathcal{Z}_{1}$ to
represent the adjoint transformation, where the potential $\Upsilon$

$\displaystyle\Omega^{{\dagger}}\Upsilon(\mathcal{Y}_{[1]},\;\mathcal{Z}_{[1]})+\Upsilon(\mathcal{Y}_{[1]},\;\mathcal{Z}_{[1]})\Lambda=\mathcal{Z}_{[1]}^{{\dagger}}\mathcal{Y}_{[1]},$
(6.3) $\displaystyle\left(\lambda
I+\Omega^{{\dagger}}\right)\Upsilon(\varphi_{[1]},\;\mathcal{Z}_{[1]})=\mathcal{Z}_{[1]}^{{\dagger}}\varphi_{[1]},$
(6.4)
$\displaystyle\Upsilon(\mathcal{Y}_{[1]},\;\psi_{[1]})\left(\mu^{{\dagger}}I+\Lambda\right)=\psi_{[1]}^{{\dagger}}\mathcal{Y}_{[1]}.$
(6.5)

The transformed operators $\widetilde{\Gamma}=\Gamma_{[2]}$ and
$\widetilde{\Delta}={\Delta}_{[2]}$ is defined as

$\displaystyle{\Gamma}_{[2]}$ $\displaystyle=$ $\displaystyle
D_{\mathcal{Y}_{[1]},\;\mathcal{Z}_{[1]}}{\Gamma}_{[1]}D_{\mathcal{Y}_{[1]},\;\mathcal{Z}_{[1]}}^{-1},$
(6.6) $\displaystyle{\Delta}_{[2]}$ $\displaystyle=$ $\displaystyle
D_{\mathcal{Y}_{[1]},\;\mathcal{Z}_{[1]}}{\Delta}_{[1]}D_{\mathcal{Y}_{[1]},\;\mathcal{Z}_{[1]}}^{-1},$
(6.7)

with generic eigenfunction

$\displaystyle\varphi_{[2]}$ $\displaystyle:=$ $\displaystyle
D_{\mathcal{Y}_{[1]},\;\mathcal{Z}_{[1]}}(\varphi_{[1]}),$ (6.8)
$\displaystyle=$
$\displaystyle\varphi_{[1]}-\mathcal{Y}_{[1]}\Upsilon(\mathcal{Y}_{[1]},\;\mathcal{Z}_{[1]})^{-1}(I+\lambda
I\Omega^{-\dagger})\Upsilon(\varphi_{[1]},\;\mathcal{Z}_{[1]}).$

and with generic adjoint eigenfunction,

$\displaystyle\psi_{[2]}$ $\displaystyle:=$ $\displaystyle
D_{\mathcal{Y}_{[1]},\;\mathcal{Z}_{[1]}}^{-{\dagger}}(\psi_{[1]}),$ (6.9)
$\displaystyle=$
$\displaystyle\psi_{[1]}-\mathcal{Z}_{[1]}\Upsilon(\mathcal{Y}_{[1]},\;\mathcal{Z}_{[1]})^{-{\dagger}}(I+\mu
I\Lambda^{-\dagger})\Upsilon(\mathcal{Y}_{[1]},\;\psi_{[1]})^{{\dagger}}.$

For the nth iteration of the binary DT, we choose the eigenfunction
$\mathcal{Y}_{[n]}$ that defines the transformation from $\Gamma_{[n]}$,
$\Delta_{[n]}$ to $\Gamma_{[n+1]}$, $\Delta_{[n+1]}$. Similarly, we choose the
eigenfunction $\mathcal{Z}_{[n]}$ for the adjoint transformation from
$\Gamma_{[n]}^{\dagger}$, $\Delta_{[n]}^{\dagger}$ to
$\Gamma_{[n+1]}^{\dagger}$, $\Delta_{[n+1]}^{\dagger}$. The Lax operators
$\Gamma_{[n]}$ and $\Delta_{[n]}$ exhibit covariance under the binary Darboux
transformation.

$D_{\mathcal{Y}_{[n]},\;\mathcal{Z}_{[n]}}=I-\mathcal{Y}_{[n]}\Upsilon(\mathcal{Y}_{[n]},\;\mathcal{Z}_{[n]})^{-1}\Omega^{-{\dagger}}\mathcal{Z}_{[n]}^{\dagger},$
(6.10)

and adjoint binary DT

$D_{\mathcal{Y}_{[n]},\;\mathcal{Z}_{[n]}}^{-{\dagger}}=I-\mathcal{Z}_{[n]}\Upsilon(\mathcal{Y}_{[n]},\;\mathcal{Z}_{[n]})^{-1}\Lambda^{-{\dagger}}\mathcal{Y}_{[n]}^{\dagger}.$
(6.11)

The transformed operators

$\displaystyle{\Gamma}_{[n+1]}$ $\displaystyle=$ $\displaystyle
D_{\mathcal{Y}_{[n]},\;\mathcal{Z}_{[n]}}{\Gamma}_{[n]}D_{\mathcal{Y}_{[n]},\;\mathcal{Z}_{[n]}}^{-1},$
(6.12) $\displaystyle{\Delta}_{[n+1]}$ $\displaystyle=$ $\displaystyle
D_{\mathcal{Y}_{[n]},\;\mathcal{Z}_{[n]}}{\Delta}_{[n]}D_{\mathcal{Y}_{[n]},\;\mathcal{Z}_{[n]}}^{-1},$
(6.13)

have generic eigenfunction

$\varphi_{[n+1]}=\varphi_{[n]}-\mathcal{Y}_{[n]}\Upsilon(\mathcal{Y}_{[n]},\;\mathcal{Z}_{[n]})^{-1}(I+\lambda
I\Omega^{-{\dagger}})\Upsilon(\varphi_{[n]},\;\mathcal{Z}_{[n]}),$ (6.14)

and generic adjoint eigenfunction;

$\psi_{[n+1]}=\psi_{[n]}-\mathcal{Z}_{[n]}\Upsilon(\mathcal{Y}_{[n]},\;\mathcal{Z}_{[n]})^{-1}(I+\mu
I\Lambda^{-{\dagger}})\Upsilon(\mathcal{Y}_{[n]},\;\psi_{[n]})^{\dagger}.$
(6.15)

By introducing the matrices $\Xi=(\mathcal{Y}_{1},\;...,\;\mathcal{Y}_{n})$
and $Z=(\mathcal{Z}_{1},\;...,\;\mathcal{Z}_{n})$, we can represent these
findings in the framework of quasi-Gramians, yielding the following
expression:

$\varphi_{[n+1]}=\left|\begin{array}[]{cc}\Upsilon(\Xi,\;Z)&(I+\lambda
I\hat{\Omega}^{-\dagger})\Upsilon(\varphi,\;Z)\\\
\Xi&\framebox{$\varphi$}\end{array}\right|,\quad\psi_{[n+1]}=\left|\begin{array}[]{cc}\Upsilon(\Xi,\;Z)^{{\dagger}}&(I+\mu
I\hat{\Lambda}^{-\dagger})\Upsilon(\Xi,\;\psi)^{{\dagger}}\\\
Z&\framebox{$\psi$}\end{array}\right|,$ (6.16)

where
$\hat{\Omega}=\text{diag}(\Omega,\;...,\;\Omega),\;\hat{\Lambda}=\text{diag}(\Lambda,\;...,\;\Lambda)$,
with both $\Omega$ and $\Lambda$ being 2 × 2 matrices.

##  7\. Quasi-Gramian solutions

In this section, we now determine the effect of binary DT
$D_{\mathcal{Y},\;\mathcal{Z}}=I-\xi\Upsilon(\mathcal{Y},\;\mathcal{Z})^{-1}\Omega^{-{\dagger}}\mathcal{Z}^{{\dagger}}$
on the Lax operator $\Gamma$, with $\mathcal{Y}_{1},\;...,\;\mathcal{Y}_{n}$
being eigenfunctions of $\Gamma$. Similarly, let
$\mathcal{Z}_{1},\;...,\;\mathcal{Z}_{n}$ denote the eigenfunctions of the
adjoint Lax operator ${\Gamma}^{{\dagger}}$. The same results apply to the
operators ${\Delta}$ and ${\Delta}^{{\dagger}}$.

As the binary DT
$D_{\mathcal{Y}_{[1]},\;\mathcal{Z}_{[1]}}=D_{\mathcal{Y},\;\mathcal{Z}}$
combines two ordinary DTs, $D_{\mathcal{Y}_{[1]}}=D_{\mathcal{Y}}$ and
$D_{\hat{\mathcal{Y}}_{[1]}}=D_{\hat{\mathcal{Y}}}$, it follows that the Lax
operator $\Gamma$ is transformed into a new Lax operator $\widehat{\Gamma}$
under the binary DT, given by:

$\widehat{\Gamma}={D_{\widehat{\mathcal{Y}}}}{\Gamma}{D_{\widehat{\mathcal{Y}}}^{-1}},$
(7.1)

we get

$\widehat{\mathcal{U}}=\mathcal{U}+[\mathcal{J},\;\mathcal{Y}\Upsilon(\mathcal{Y},\;\mathcal{Z})^{-1}\mathcal{Z}^{{\dagger}}].$
(7.2)

Since $\mathcal{U}=[\mathcal{F},\;\mathcal{J}]$, so that

$\widehat{\mathcal{F}}=\mathcal{F}-\mathcal{Y}\Upsilon(\mathcal{Y},\;\mathcal{Z})^{-1}\mathcal{Z}^{{\dagger}}.$
(7.3)

After $n$ iterations of applying the binary DT
$D_{\mathcal{Y},\;\mathcal{Z}}$, the resulting expression is given by:

$\displaystyle\mathcal{F}_{[n+1]}=\mathcal{F}_{[n]}-\mathcal{Y}_{[n]}\Upsilon(\mathcal{Y}_{[n]},\;\mathcal{Z}_{[n]})^{-1}\mathcal{Z}_{[n]}^{{\dagger}}=\mathcal{F}-\sum_{i=1}^{n}\mathcal{Y}_{[i]}\Upsilon(\mathcal{Y}_{[i]},\;\mathcal{Z}_{[i]})^{-1}\mathcal{Z}_{[i]}^{{\dagger}}.$
(7.4)

Using the notation
$\mathcal{F}_{[1]}=\mathcal{F},\;\mathcal{F}_{[2]}=\widehat{\mathcal{F}},\;\mathcal{Y}_{[1]}=\mathcal{Y}_{1}$
and $\mathcal{Z}_{[1]}=\mathcal{Z}_{1}$, and introducing the matrices
$\Xi=(\mathcal{Y}_{1},\;...,\;\mathcal{Y}_{n})$ and
$Z=(\mathcal{Z}_{1},\;...,\;\mathcal{Z}_{n})$, the result (7.4) can be
expressed in the form of quasi-Gramian as follows:

$\mathcal{F}_{[n+1]}=\mathcal{F}+\left|\begin{array}[]{cc}\Upsilon(\Xi,\;Z)&Z^{{\dagger}}\\\
\Xi&\framebox{$0_{2}$}\end{array}\right|,$ (7.5)

It is worth noting that each
$\mathcal{Y}_{i},\;\mathcal{Z}_{i},\;i=1,\;...,\;n$ is a $2\times 2$ matrix.
Given that our system of the nc-HNLS equation is formulated in terms of non-
commutative objects $u,\;u^{{\dagger}}$, we find it more appropriate to
express the quasi-Gramian solution (7.5) in terms of these objects. Thus, we
introduce the matrices $\mathcal{Y}_{i}$ following a similar approach as in
the quasi-Wronskian case. We also define $Z=\Xi Q^{{\dagger}}$, where $Q$
represents a constant matrix of size $2n\times 2n$ and † denotes the Hermitian
conjugate. It is noted that $\Xi$ and $Z$ adhere to the same dispersion
relation and remain unchanged when multiplied by a constant matrix.
Consequently, the quasi-Gramian solution (7.5) can also be represented as
follows:

$\mathcal{F}_{[n+1]}=\mathcal{F}+\left(\begin{array}[]{cc}\left|\begin{array}[]{cc}\Upsilon(\Xi,\;Z)&Q\varphi^{{\dagger}}\\\
\varphi&\framebox{$0$}\end{array}\right|&\left|\begin{array}[]{cc}\Upsilon(\Xi,\;Z)&Q\chi^{{\dagger}}\\\
\varphi&\framebox{$0$}\end{array}\right|\\\ \\\
\left|\begin{array}[]{cc}\Upsilon(\Xi,\;Z)&Q\varphi^{{\dagger}}\\\
\chi&\framebox{$0$}\end{array}\right|&\left|\begin{array}[]{cc}\Upsilon(\Xi,\;Z)&Q\chi^{{\dagger}}\\\
\chi&\framebox{$0$}\end{array}\right|\end{array}\right),$ (7.6)

where $\varphi=(\varphi_{1},\;...,\;\varphi_{n})$ and
$\chi=(\chi_{1},\;...,\;\chi_{n})$ are row vectors. Thus, quasi-Gramian
expression for $u$ and $u^{{\dagger}}$ are given by

$\displaystyle
u_{[n+1]}=u+2\dot{\imath}\left|\begin{array}[]{cc}\Upsilon(\Xi,\;Z)&Q\chi^{{\dagger}}\\\
\varphi&\framebox{$0$}\end{array}\right|,\quad
u_{[n+1]}^{{\dagger}}=u^{{\dagger}}+2\dot{\imath}\left|\begin{array}[]{cc}\Upsilon(\Xi,\;Z)&Q\varphi^{{\dagger}}\\\
\chi&\framebox{$0$}\end{array}\right|.$ (7.11)

Equation (7.11) represents the quasi-Gramian solutions for the nc-HNLS
equation (3.11). If we relax the non-commutativity condition, the equation can
be simplified and expressed as a ratio of simple Gramians. In the limit of
commutativity, we obtain the following expression

$\displaystyle
u_{[n+1]}=u-2\dot{\imath}\frac{\left|\begin{array}[]{cc}\Upsilon(\Xi,\;Z)&Q\chi^{{\dagger}}\\\
\phi&0\end{array}\right|}{\left|\Delta(\Xi,\;Z)\right|},\quad\bar{u}_{[n+1]}=\bar{u}-2\dot{\imath}\frac{\left|\begin{array}[]{cc}\Upsilon(\Xi,\;Z)&Q\varphi^{{\dagger}}\\\
\chi&0\end{array}\right|}{\left|\Delta(\Xi,\;Z)\right|}.$ (7.16)

These expressions define the Gramian solutions of the higher-order NLS
equation.

##  8\. Explicit solutions

When $u=0$, spectral problem (3.1) and (3.2) has the solution

$\varphi_{j}=e^{\dot{\imath}\zeta_{j}},\;\quad\chi_{j}=e^{-\dot{\imath}\zeta_{j}},\quad\zeta_{j}=-\lambda_{j}x+2\lambda_{j}^{2}(4\gamma\,\lambda_{j}^{2}-2\lambda_{j}\alpha_{1}-\alpha_{2})t,$
(8.1)

To simplify the notation and work with only $\varphi_{1},\ldots,\varphi_{n}$
and $\chi_{1},\ldots,\chi_{n}$, we introduce the following relabeling. We
redefine $\varphi_{i}$ as $\varphi_{\frac{i+1}{2}}$ for odd values of $i$
(i.e., $i=1,\ldots,2n-1$), and set $\varphi_{i}=0$ for even values of $i$
(i.e., $i=2,4,\ldots,2n$). Similarly, we relabel $\chi_{i}$ as
$\chi_{\frac{i}{2}}$ for even values of $i$, and $\chi_{i}=0$ for odd values
of $i$. We then have

$\varphi=(\varphi_{1},\;0,\;\varphi_{2},\;0,\;...,\;\varphi_{n},\;0),\quad\chi=(0,\;\chi_{1},\;0,\;\chi_{2},\;...,\;0,\;\chi_{n}),$
(8.2)

Using the notation $\mathcal{Y}_{i}=\text{diag}(\varphi_{i},\chi_{i})$ for
$i=1,\ldots,n$, where $\varphi_{i}$ and $\chi_{i}$ are given in (8.1), we now
focus on the commutative case. In this case, we can express the Gramian
solution (7.16) as follows.

$u_{[n+1]}=-2\dot{\imath}\frac{\left|\begin{array}[]{cc}\Upsilon(\Xi,\;Z)&Q\chi^{{\dagger}}\\\
\varphi&0\end{array}\right|}{\left|\Upsilon(\Xi,\;Z)\right|}=-2\dot{\imath}\frac{\mathcal{G}}{\mathcal{R}},\;\text{say}$
(8.3)

where $\mathcal{R}=\text{det}(\mathcal{W})$ and

$\mathcal{W}=\Upsilon(\Xi,\;Z)=\int Z^{{\dagger}}\mathcal{J}\Xi dx+I_{2n},$
(8.4)

Here, $I_{2n}$ is the identity matrix with dimensions $2n\times 2n$.
Constructing the matrix $\Xi$ involves arranging the eigenfunctions
$\mathcal{Y}_{1},\;\mathcal{Y}_{2},\;...,\;\mathcal{Y}_{n}$, where each
$\mathcal{Y}_{i}$ represents an eigenfunction of the Lax operators $\Gamma$
and $\Delta$, presented as a $2\times 2$ matrix. Similarly, assembling the
matrix $Z$ involves the eigenfunctions
$\mathcal{Z}_{1},\;\mathcal{Z}_{2},\;...,\;\mathcal{Z}_{n}$, with
$\mathcal{Z}_{i}$ serving as eigenfunctions of the adjoint Lax operators
$\Gamma^{\dagger}$ and $\Delta^{\dagger}$. The matrix $\Upsilon(\Xi,\;Z)$ is a
matrix $2n\times 2n$, its entries being scalar components $(1\times 1)$. As we
proceed to discuss the non-commutative case, we will consider every component
of $\mathcal{Y}_{i}$ and $\mathcal{Z}_{i}$ as a matrix. Presenting the matrix
$Q$, a constant matrix $2n\times 2n$, we define $Z$ as the outcome of
multiplying $\Xi$ by the Hermitian adjoint of $Q$, denoted as $Z=\Xi
Q^{\dagger}$, allowing us to express the equation (8.4) as

$\mathcal{W}=Q\int\Xi^{{\dagger}}\mathcal{J}\Xi dx+I_{2n}=Q\Theta+I_{2n},$
(8.5)

where

$\Theta=\left(\begin{array}[]{ccccc}-\dot{\imath}\int^{x}\varphi_{1}^{\ast}\varphi_{1}dx&0_{2}&\ldots&-\dot{\imath}\int^{x}\varphi_{1}^{\ast}\varphi_{n}dx&0_{2}\\\
0_{2}&\dot{\imath}\int^{x}\chi_{1}^{\ast}\chi_{1}dx&\ldots&0_{2}&\dot{\imath}\int\chi_{1}^{\ast}\chi_{n}dx\\\
\vdots&\vdots&\ddots&\vdots&\vdots\\\
-\dot{\imath}\int^{x}\varphi_{n}^{\ast}\varphi_{1}dx&0_{2}&\ldots&-\dot{\imath}\int^{x}\varphi_{n}^{\ast}\varphi_{n}dx&0_{2}\\\
0_{2}&\dot{\imath}\int^{x}\chi_{n}^{\ast}\chi_{1}dx&\ldots&0_{2}&\dot{\imath}\int^{x}\chi_{n}^{\ast}\chi_{n}dx\\\
\end{array}\right).$ (8.6)

To derive one soliton ($n=1,\;u_{[2]}\equiv u_{1}$) solution for the
commutative Hirota equation. We choose
$Q=\left(\begin{smallmatrix}q_{1}&q_{2}\\\
q_{2}&q_{1}\end{smallmatrix}\right)$, one soliton solution is given by

$u_{1}=\frac{8\,\mathrm{I}c_{1}q_{\mathit{2}}\lambda_{I}^{2}{\mathrm{e}}^{2\,\mathrm{I}\xi_{1}}}{4\,\mathrm{I}\lambda_{I}c_{1}q_{\mathit{1}}\cosh\\!\left(\xi_{2}\right)-4\lambda_{I}^{2}c_{1}^{2}+q_{\mathit{1}}^{2}-q_{\mathit{2}}^{2}}$
(8.7)

where

$\displaystyle\xi_{1}=[8\left(\lambda_{R}^{4}-6\lambda_{R}^{2}\lambda_{I}^{2}+\lambda_{I}^{4}\right)\gamma+2\left(-\lambda_{R}^{2}+\lambda_{I}^{2}\right)\alpha_{2}+4\left(-\lambda_{R}^{3}+3\lambda_{R}\lambda_{I}^{2}\right)\alpha_{1}]t-\lambda_{R}x,$
(8.8)
$\displaystyle\xi_{2}=8[8\left(-\lambda_{R}^{3}+\lambda_{R}\lambda_{I}^{2}\right)\gamma+\lambda_{R}\alpha_{2}+\left(-\lambda_{I}^{2}+3\lambda_{R}^{2}\right)\alpha_{1}]\lambda_{I}t+2x\lambda_{I},$

To visualize this soliton solution, we plot the propagation of the $u_{1}$
soliton in the commutative case with a velocity of
$8\left(-\lambda_{R}^{3}+\lambda_{R}\lambda_{I}^{2}\right)\gamma+\lambda_{R}\alpha_{2}+\left(-\lambda_{I}^{2}+3\lambda_{R}^{2}\right)\alpha_{1}$,
where $\lambda=\lambda_{R}+i\lambda_{I}$. Fig. 2 illustrates the behavior of
the soliton over time.

(a)

(b)

Figure 2. Evolution of the solution (8.7) with parameters
$\alpha_{1}=1.5,\;\alpha_{2}=\gamma=1,\;c_{1}=0.5,\;q_{1}=2,\;q_{2}=-1,\;(a)\;\lambda=0.1+0.5\dot{\imath},\;(b)\;\lambda=0.5\dot{\imath}$.

###  8.1. Noncommutative case

We now discuss the noncommutative case. It has been shown in [10] that the
behavior of matrix solitons differs from their scalar counterparts. Unlike
scalar solitons, which maintain their amplitudes unchanged during
interactions, matrix solitons undergo transformations that depend on certain
rules. These transformations affect the amplitudes, which are determined by
vectors rather than individual values in the noncommutative setting. When
considering the case of $n=1$, we choose the solutions $\varphi$ and $\chi$ of
the Lax pair to be $2\times 2$ matrices, given by

$\varphi_{j}=e^{\dot{\imath}\zeta_{j}}I_{2},\;\quad\chi_{j}=e^{-\dot{\imath}\zeta_{j}}I_{2},\quad\zeta_{j}=-\lambda_{j}x+2\lambda_{j}^{2}(4\gamma\,\lambda_{j}^{2}-2\lambda_{j}\alpha_{1}-\alpha_{2})t,$
(8.9)

Each entry in $\varphi$ and $\chi$ and constant matrix $Q$ is a $2\times 2$
matrix, so that these matrices are given by

$\displaystyle\varphi=\left(\begin{array}[]{cccc}\varphi_{1}&0&0&0\\\
0&\varphi_{1}&0&0\end{array}\right),\;\chi=\left(\begin{array}[]{cccc}0&0&\chi_{1}&0\\\
0&0&0&\chi_{1}\end{array}\right),\;Q=\left(\begin{array}[]{cccc}q_{11}&q_{12}&q_{13}&q_{14}\\\
q_{12}&q_{11}&q_{14}&q_{13}\\\ q_{13}&q_{14}&q_{33}&q_{34}\\\
q_{14}&q_{13}&q_{34}&q_{33}\end{array}\right),$ (8.18)

Therefore, the quasi-Gramian expression for $u_{2}$ (which we will now denote
as $u^{1}$ for the non-commutative case) can be expressed as follows

$\displaystyle u^{1}$ $\displaystyle=$
$\displaystyle{2\dot{\imath}}\left(\begin{array}[]{cc}\left|\begin{array}[]{ccccc}&&&&\\\
&&\Upsilon(\Xi,\;Z)&&\chi^{11}\\\ &&&&\\\
\varphi_{1}&0&0&0&\framebox{$0$}\end{array}\right|&\left|\begin{array}[]{ccccc}&&&&\\\
&&\Upsilon(\Xi,\;Z)&&\chi^{12}\\\ &&&&\\\
\varphi_{1}&0&0&0&\framebox{$0$}\end{array}\right|\\\
\left|\begin{array}[]{ccccc}&&&&\\\ &&\Upsilon(\Xi,\;Z)&&\chi^{11}\\\ &&&&\\\
0&\varphi_{1}&0&0&\framebox{$0$}\end{array}\right|&\left|\begin{array}[]{ccccc}&&&&\\\
&&\Upsilon(\Xi,\;Z)&&\chi^{12}\\\ &&&&\\\
0&\varphi_{1}&0&0&\framebox{$0$}\end{array}\right|\end{array}\right),$ (8.37)
$\displaystyle=$
$\displaystyle{2\dot{\imath}}\left(\begin{array}[]{cc}u_{11}&u_{12}\\\
u_{21}&u_{22}\end{array}\right),\;\text{say},$ (8.40)

where
$\chi^{11}=(q_{13}\chi_{1}^{\ast}\;q_{14}\chi_{1}^{\ast}\;q_{33}\chi_{1}^{\ast}\;q_{34}\chi_{1}^{\ast})^{{\dagger}},\;\chi^{12}=(q_{14}\chi_{1}^{\ast}\;q_{13}\chi_{1}^{\ast}\;q_{34}\chi_{1}^{\ast}\;q_{33}\chi_{1}^{\ast})^{{\dagger}}$
and $\Upsilon$ is the potential defined in (8.4) with each entry being a
$2\times 2$ matrix.

In the context of non-commutative systems, the soliton solution (8.37) is
intricately influenced by both the spectral parameter $\lambda$ and the
elements composing the matrix $Q$. When specific entries, such as $q_{13}$ and
$q_{14}$, are deliberately set to zero, the resulting solutions for $u_{11}$,
$u_{12}$, $u_{21}$, and $u_{22}$ appear trivial. And, where
$q_{13}=q_{14}=-2$, the graphical representations of solutions $u_{11}$,
$u_{12}$, $u_{21}$, and $u_{22}$ manifesting as a single consolidated plot
instead of the originally intended four. Noteworthy is the fact that under
this condition, all solitons propagate with a consistent amplitude of $0.2163$
units. Unlike this symmetry, when $q_{13}\neq q_{14}$, the resulting graphs
exhibit a variety of double- and single-peaked patterns for each component of
the matrix $u^{1}$ (see Fig. 4). It is worth highlighting that solitons
$u_{12}$ and $u_{21}$ advance with an amplitude of 0.3860 units of large peak
and 0.2317 units of small peak, while the solitons $u_{11}$ and $u_{22}$
display amplitude of 1.9797 units. Additionally, an intriguing situation
unfolds when we choose $q_{12}=-q_{14}=-2$ and
$q_{11}=q_{13}=q_{33}=q_{34}=0$. Under these conditions, we notice a
consistent single-peaked soliton for $u_{11}$ and $u_{22}$ with an amplitude
of 0.2214 units, while simultaneously observing a kink pattern in $u_{12}$ and
$u_{21}$ of maximum height 0.8858 units (as seen in Fig. 5).

Figure 3. Evolution of the solution (8.37) with parameters
$\alpha_{1}=1.5,\;\alpha_{2}=\gamma=1,\;c_{1}=0.5,\;q_{11}=0.5,\;q_{12}=0,\;q_{13}=q_{14}=-1,\;q_{33}=-0.2,\;q_{34}=-0.1,\;\lambda=-0.1+0.5\dot{\imath}$.

(a)

(b)

(c)

(d)

Figure 4. Evolution of the solution (8.37) with parameters
$\alpha_{1}=1.5,\;\alpha_{2}=\gamma=1,\;c_{1}=0.5,\;q_{11}=0.5,\;q_{12}=-0.2,\;q_{13}=-1,\;q_{14}=-0.6,\;q_{33}=0,\;q_{34}=0.1,\;\lambda=-0.1+0.5\dot{\imath}$.

(a)

(b)

(c)

(d)

Figure 5. Evolution of the solution (8.37) with parameters
$\alpha_{1}=1.5,\;\alpha_{2}=\gamma=1,\;c_{1}=0.5,\;q_{11}=0,\;q_{12}=-2,\;q_{13}=0,\;q_{14}=2,\;q_{33}=0,\;q_{34}=0,\;\lambda=-0.1+0.5\dot{\imath}$.

The solution (8.37) includes several particular cases. In the case where both
$\alpha_{2}$ and $\gamma$ are zero, the solution (8.37) takes a different
form. It becomes a solution of a non-commutative generalization of the complex
modified Korteweg-de Vries (KdV) equation and is further reduced to the
standard modified KdV equation when the variable $u$ is real-valued.
Additionally, if we set $\gamma$ to zero, we obtain the solution of the non-
commutative extension of the Hirota equation. On the other hand, when we
simultaneously set $\alpha_{1}$ and $\alpha_{2}$ to zero, we have the solution
of the non-commutative Lakshmanan–Porsezian–Daniel (LPD) equation. Lastly, if
we set $\alpha_{1}$ and $\gamma$ to zero, we get the solution of non-
commutative generalization of the nonlinear Schrödinger (NLS) equation.

(a)

(b)

(c)

(d)

Figure 6. The profiles of $u_{11}$ and $u_{12}$ with $\nu=\alpha_{2}=0$ in (a)
and (b) and, with $\nu=0$ in (c) and (d). All other parameters are the same as
in Fig. 5.

(a)

(b)

(c)

(d)

Figure 7. The profiles of $u_{11}$ and $u_{12}$ with $\alpha_{1}=\alpha_{2}=0$
in (a) and (b) and, with $\alpha_{1}=\nu=0$ in (c) and (d). All other
parameters are the same as in Fig. 5.

In summary, studying the non-commutative version is important because it gives
us different choices for arranging solitons. These arrangements depend not
only on the spectral parameter $\lambda$ but also on values in a matrix.
Similarly, two soliton solutions for the non-commutative case are depicted in
Figs. 8-10

(a)

(b)

(c)

(d)

Figure 8. The profiles of $u_{ij},\;i,\;j=1,\;2$ with
$\alpha_{1}=0.5,\;\alpha_{2}=\nu=1,\;\lambda_{1}=0.5i,\;\lambda_{2}=-0.1-0.1i$.
All other parameters are the same as in Figure 5.

(a)

(b)

(c)

(d)

Figure 9. The profiles of $u_{ij},\;i,\;j=1,\;2$ with
$\alpha_{1}=0.5,\;\alpha_{2}=\nu=1,\;\lambda_{1}=0.6i,\;\lambda_{2}=-1.1-1.1i$.

(a)

(b)

(c)

(d)

Figure 10. The profiles of $u_{ij},\;i,\;j=1,\;2$ with
$\alpha_{1}=0.5,\;\alpha_{2}=\nu=1,\;\lambda_{1}=0.1+0.5i,\;\lambda_{2}=-0.5i$.

From Fig. 8, it can be seen that solitons with breather-like structures can be
seen in $u_{ij},i,\;j=1,\;2$ components. This is because when the force
between the two solitons is strong enough, bound-state solitons can merge into
breather-like solitons. This process is known as soliton fusion, and it occurs
when two solitons combine to form a single soliton.

Figs. 9-10 show two-soliton solutions of nc-HNLS equation, which correspond to
the interaction of two individual solitons, that is, two distinct lumps of
energy moving at different velocities that interact and scatter without
changing their shape. Other multisoliton expressions can be obtained by
repeatedly applying the Darboux transformation to the seed solution. The
three-soliton configuration represents three distinct amplitudes of soliton
scattering. Note that, we have omitted the explicit expression of soliton
solutions for non-commutative cases as it is long and cumbersome.

Acknowledgment: We acknowledge the support of the National Natural Science
Foundation of China, Nos. 11835011 and 12375006.

Data Availability: Not Applicable.

##  9\. Concluding remarks

This study explored the non-commutative extension of the higher-order
nonlinear Schrödinger equation. Through Darboux and binary Darboux
transformations, we expressed in quasi-Wronskian and quasi-Gramian forms.
These solutions were intricately linked to the non-commutative higher-order
nonlinear schrödidnger equation and its associated Lax pair. Importantly, we
demonstrated single-, double-peaked, kink, and bright solitons in non-
commutative settings. Our approach offers a powerful tool not only for
understanding the non-commutative higher order nonlinear Schrödinger equation
but also for constructing multisolitons in related non-commutative integrable
systems.

## References

  * [1] Ablowitz, M. J., Prinari, B., and Trubatch, A. D. (2004). Discrete and continuous nonlinear Schrödinger systems (Vol. 302). Cambridge University Press.
  * [2] Kivshar, Y. S., and Agrawal, G. P. (2003). Optical solitons: from fibers to photonic crystals. Academic press.
  * [3] Agrawal, G. P. (2000). Nonlinear fiber optics. In Nonlinear Science at the Dawn of the 21st Century (pp. 195-211). Berlin, Heidelberg: Springer Berlin Heidelberg.
  * [4] Daniel, M., and Deepamala, K. (1995). Davydov soliton in alpha helical proteins: higher order and discreteness effects. Physica A: Statistical Mechanics and its Applications, 221(1-3), 241-255.
  * [5] Lakshmanan, M., Porsezian, K., and Daniel, M. (1988). Effect of discreteness on the continuum limit of the Heisenberg spin chain. Physics Letters A, 133(9), 483-488.
  * [6] Porsezian, K., Daniel, M., andLakshmanan, M. (1992). On the integrability aspects of the one‐dimensional classical continuum isotropic biquadratic Heisenberg spin chain. Journal of mathematical physics, 33(5), 1807-1816.
  * [7] Ankiewicz, A., Wang, Y., Wabnitz, S., and Akhmediev, N. (2014). Extended nonlinear Schrödinger equation with higher-order odd and even terms and its rogue wave solutions. Physical Review E, 89(1), 012907.
  * [8] Ankiewicz, A., and Akhmediev, N. (2014). Higher-order integrable evolution equation and its soliton solutions. Physics Letters A, 378(4), 358-361.
  * [9] Nimmo, J. J. C., Gilson, C. R., and Ohta, Y. (2000). Applications of Darboux transformations to the self-dual Yang-Mills equations. Theoretical and Mathematical Physics, 122(2), 239-246.
  * [10] Goncharenko, V. M. (2001). Multisoliton solutions of the matrix KdV equation. Theoretical and Mathematical Physics, 126(1), 81-91.
  * [11] Lakshmanan, M. (1977). Continuum spin system as an exactly solvable dynamical system. Physics Letters A, 61(1), 53-54.
  * [12] Hirota, R. (1973). Exact envelope‐soliton solutions of a nonlinear wave equation. Journal of Mathematical Physics, 14(7), 805-809.
  * [13] Minwalla, S., Van Raamsdonk, M., and Seiberg, N. (2000). Noncommutative perturbative dynamics. Journal of High Energy Physics, 2000(02), 020.
  * [14] Furuta, K., and Inami, T. (2000). Ultraviolet property of noncommutative Wess–Zumino–Witten model. Modern Physics Letters A, 15(15), 997-1002.
  * [15] Seiberg, N., and Witten, E. (1999). String theory and noncommutative geometry. Journal of High Energy Physics, 1999(09), 032.
  * [16] Dimakis, A., and Muller-Hoissen, F. (2000). Bicomplexes, integrable models, and noncommutative geometry. arXiv preprint hep-th/0006005.
  * [17] Lechtenfeld, O., and Popov, A. D. (2001). Noncommutative multi-solitons in 2+ 1 dimensions. Journal of High Energy Physics, 2001(11), 040.
  * [18] Lechtenfeld, O., Mazzanti, L., Penati, S., Popov, A. D., and Tamassia, L. (2005). Integrable noncommutative sine-Gordon model. Nuclear Physics B, 705(3), 477-503.
  * [19] Gilson, C. R., and Nimmo, J. J. C. (2007). On a direct approach to quasideterminant solutions of a noncommutative KP equation. Journal of Physics A: Mathematical and Theoretical, 40(14), 3839.
  * [20] Gilson, C. R., and Macfarlane, S. R. (2009). Dromion solutions of noncommutative Davey–Stewartson equations. Journal of Physics A: Mathematical and Theoretical, 42(23), 235202.
  * [21] Etingof, P., Gelfand, I., and Retakh, V. (1997). Factorization of differential operators, quasideterminants, and nonabelian Toda field equations. arXiv preprint q-alg/9701008.
  * [22] Riaz, H. W. A., and ul Hassan, M. (2018). Multisoliton solutions of integrable discrete and semi-discrete principal chiral equations. Communications in Nonlinear Science and Numerical Simulation, 54, 416-427.