Data-driven 2D stationary quantum droplets and wave propagations in the amended GP equation with two potentials via deep neural networks learning
Jin Song Zhenya Yan




































language=python


    numbers=left,
    numberstyle= ,
    basicstyle=,
    keywordstyle= ,
    commentstyle= ,
    frame=shadowbox, 
    rulesepcolor=  ,
    escapeinside=“, 
    breaklines,
    showspaces=false,
    xleftmargin=2em, aboveskip=1em,
    framexleftmargin=2em



    corollaryCorollarylemmaLemmapropositionPropositionremarkRemark=-0.6in =-0.80in
        =-0.3in =0.00in
        =225mm =175mm
        =0.1in
        [4]
        =13pt plain
        
        
        
        
        
        
        
        
        
        
        
        =14pt plain
        =16pt 
        Data-driven 2D stationary quantum droplets and wave propagations in the amended GP equation with two potentials via deep neural networks learning
        
        
        Jin Song^1,2 and Zhenya Yan^1,2,*
        [^*Corresponding author. Email address:
        zyyan@mmrc.iss.ac.cn]
        ^1KLMM, Academy of Mathematics and Systems Science,
        Chinese Academy of Sciences, Beijing 100190, China 
        ^2School of Mathematical Sciences, University of Chinese Academy of
        Sciences, Beijing 100049, China
        
        Abstract: In this paper, we develop a systematic deep learning approach to solve two-dimensional (2D) stationary quantum droplets (QDs) and investigate their wave propagation in the 2D amended Gross–Pitaevskii equation with Lee-Huang-Yang correction and two kinds of potentials. Firstly, we use the initial-value iterative neural network (IINN) algorithm for 2D stationary quantum droplets of stationary equations. Then the learned stationary QDs are used as the initial value conditions for physics-informed neural networks (PINNs) to explore their evolutions in the some space-time region. Especially, we consider two types of potentials, one is the 2D quadruple-well Gaussian potential and the other is the -symmetric HO-Gaussian potential, which lead to spontaneous symmetry breaking and the generation of multi-component QDs. The used deep learning method can also be applied to study wave propagations of other nonlinear physical models.
        
        
        
        
        
        =13pt
        
        
        
        
        § INTRODUCTION
        
        
        Recent intensive researches have focused on quantum droplets (QDs), a new state of liquid matter <cit.>. QDs are characterized by a delicate balance between mutual attraction and repulsion, leading to their unique properties. QDs have potential applications in ultracold atoms and superfluids and have been studied widely <cit.>. As the ultra-dilute liquid matter, QDs are nearly incompressible, self-sustained liquid droplets, with distinctive properties such as extremely low densities and temperatures <cit.>. The Lee-Huang-Yang (LHY) effect <cit.>, driven by quantum fluctuations, has been introduced to prevent QDs from collapsing due to mean-field approximation, enabling the prediction of stable QDs in weakly interacting Bose-Einstein condensates (BECs) <cit.>.
        
        Experimental realizations of QDs have been achieved in various systems, including single-component dipolar bosonic gases, binary Bose-Bose mixtures of different atomic states in ^39K, and in the heteronuclear mixture of ^41K and ^87Rb atoms <cit.>. The accurate description of QDs has been made possible by the amended Gross-Pitaevskii (GP) equation with Lee-Huang-Yang (LHY) correction, which has been shown to agree with experimental observations <cit.>.
        The reduction of dimensionality from 3D to 2D has a significant impact on the form of the LHY term. In this case, the repulsive quartic nonlinearity is replaced by a cubic nonlinearity with an additional logarithmic factor <cit.>
        
         
        
        such that the 2D amended GP equation in the binary BECs with two mutually symmetric components trapped in a potential can be written as the following dimensionless form after scaling
        
            iψ_t=-12∇_^2ψ +2ln(2|ψ|^2)|ψ|^2ψ+U()ψ,
        
        where the complex wave function ψ=ψ(, t),
        =(x,y) stands for the 2D rescaled coordinates and t∈ℝ, ∇_^2=∂^2/∂ x^2+∂^2/∂ y^2, U() is an external potential, which can be real or complex.
        A variety of trapping configurations in BECs have allowed for the direct observation of fundamental manifestations of QDs. For instance, stable 2D anisotropic vortex QDs have been predicted in effectively 2D dipolar BECs <cit.>.  More importantly, vortical QDs have been found to be stable without the help of any potential by a systematic numerical investigation and analytical estimates <cit.>. Additionally, vortex-carrying QDs can be experimentally generated in systems with attractive inter-species and repulsive intra-species interactions, confined in a shallow harmonic trap with an additional repulsive Gaussian potential at the center <cit.>.
        Furthermore, the exploration of QDs trapped in -symmetric potentials has also been pursued <cit.>.
        
        
        Recently, there has been a surge in the development of deep neural networks for studying partial differential equations (PDEs). Various approaches, such as physics-informed neural networks (PINNs) <cit.>, deep Ritz method <cit.>, and PDE-net <cit.>, have been proposed to effectively handle PDE problems. Among them, the PINNs method incorporates the physical constraints into the loss functions, allowing the models to learn and represent the underlying physics more accurately <cit.>.
        
        Moreover, these deep learning methods have been extended to solve a wide range of PDEs in various fields <cit.>.
        
        For the general 2D stationary QDs in the form Ψ(x,y,t)=ϕ(x,y)e^-iμ t, solving for ϕ(x,y) is an important problem because ϕ serves as an initial-value condition of PINNs.
        In general, the traditional numerical methods were developed thus far to compute solitary waves, including the Petviashvili method, accelerated imaginary-time evolution (AITEM)
        method, squared-operator iteration (SOM) method, and Newton-conjugate-gradient (NCG) method <cit.>. More recently,
        we proposed a new deep learning approach called the initial-value iterative neural network (IINN) for solitary wave computations of many types of nonlinear wave euqations <cit.>, which offers a mesh-free approach by taking advantage of automatic differentiation, and could overcomes the curse of dimensionality.
        
        Motivated by the aforementioned discussions, the main objective of this paper is to develop a systematic deep learning approach to solve 2D stationary QDs, and investigate their evolutions in the amended Gross–Pitaevskii equation with potentials. Especially, we consider two types of potentials, one is 2D quadruple-well Gaussian potential and the other is -symmetric HO-Gaussian potential, which lead to spontaneous symmetry breaking and the generation of multi-component QDs.
        The remainder of this paper is arranged as follows. Firstly, we introduce the PINNs deep learning
        framework for the evolution of QDs. And then the IINN framework for stationary QDs is presented in Sec 2.
        In Sec. 3, data-driven 2D QDs in the amended GP equation with two types of potential are exhibited, respectively. Finally, we give some conclusions and discussions in Sec. 4.
        
        
        
        
        
        § THE FRAMEWORK OF  DEEP LEARNING METHOD
        
        
        In the following, we focus on the trapped stationary QDs to Eq. (<ref>) in the form ψ(, t) = ϕ()e^-iμ t, where μ stands for the chemical potential <cit.>, and lim_||→∞ϕ()=0 for ϕ()∈ℝ[]. Substituting the stationary solution into Eq. (<ref>) yields the following nonlinear stationary equation obeyed by the nonlinear localized eigenmode ϕ():
        
            μϕ=-1/2∇^2_ϕ+2ln(2|ϕ|^2)|ϕ|^2ϕ+U()ϕ.
        
        In general, it is difficult to get the explicit, exact solutions of Eq. (<ref>) with the potentials. For general parametric conditions, one can usually use the numerical iterative methods to solve Eq. (<ref>) with zero-boundary conditions by choosing the proper initial value, such as Newton-conjugate-gradient (NCG) method <cit.>, the spectral renormalization method <cit.>, and the squared-operator iteration method <cit.>.
        In this paper, we extend the deep learning IINN method for the computations of stationary QDs, and then we use the stationary QDs as the initial data to analyze the evolutions of QDs with the aid of PINNs.
        
        
        
        
         §.§ The IINN framework for stationary QDs
        
        
        Based on traditional numerical iterative methods and physics-informed neural networks (PINNs), recently we proposed the initial value iterative neural network (IINN) algorithm for solitary wave computations <cit.>.
        In the following, we will introduce the main idea of IINN method.
        Two identical fully connected neural networks are employed to learn the desired solution ϕ^* of Eq. (<ref>).
        
        For the first network, we choose an appropriate initial value ϕ_0 such that it is sufficiently close to ϕ^*. Then we randomly select N training points
        {_i}_i=1^N within the region and train the network parameters θ by minimizing the mean squared error loss ℒ_1, aiming to make the output of network ϕ̅ sufficiently close to initial value ϕ_0, where loss function ℒ_1 is defined as follows
        
            ℒ_1:=1/N∑_i=1^N|ϕ̅(_i)-ϕ_0(_i)|^2.
        
        
        For the second network, we initialize the network parameters θ={W,B} with the learned weights and biases from the first network,
        that is
        
            θ_0=argmin ℒ_1(θ).
        
        For the network output ϕ̂,
        we define the loss function ℒ_2 as follows and utilize SGD or Adam optimizer to minimize it.
        
            ℒ_2:=1/N∑_i=1^N|Lϕ̂(_i)|^2/max_i(|ϕ̂(_i)|).
        
        It should be noted that ℒ_2 is different from the loss function ℒ_0 defined in PINNs. Here we are not taking boundaries into consideration, instead we incorporate max(|ϕ̂|) to ensure that ϕ̂ does not converge to trivial solution.
        
        
        
        
        
        
        
        
        
         §.§ The PINNs framework for the evolution of QDs
        
        
        Base on the obtained the stationary QDs in Sec. 3.1, we utilize the PINNs deep learning framework <cit.> to address the data-driven solutions of Eq. (<ref>). The core concept of PINNs involves training a deep neural network to satisfy the physical laws and accurately represent the solutions for various nonlinear partial differential equations. In the case of the 2D amended GP equation (<ref>), we incorporate initial-boundary value conditions.
        
            {[ iψ_t+12∇_^2ψ -2ln(2|ψ|^2)|ψ|^2ψ-U()ψ=0,  
            (, t)∈Ω× (0, T),ψ(,0)=ϕ(),  ∈Ω, ψ(,t)|_∈∂Ω=ϕ_b(t),   t∈ [0,T], ].
        
        where ϕ() is the solution of stationary equation (<ref>), and we take ϕ_b(t)≡ 0, which is solved by the IINN method  in Sec. 3.1.
        
        We rewrite the wave-function as ψ(,t)=p(,t)+iq(,t) with p(,t) and q(,t) being its real and imaginary parts, respectively. Then the complex-valued PINNs ℱ(, t)=ℱ_p(, t)+iℱ_q(, t) with ℱ_p(, t) and ℱ_q(, t) being its real and imaginary parts, respectively can be defined as
        
            [ ℱ(, t):= iψ_t+12∇_^2ψ -2ln(2|ψ|^2)|ψ|^2ψ-U()ψ,ℱ_p(, t):= -q_t+12∇^2_p-2ln[2(p^2+q^2)](p^2+q^2)p-real(U)p+imag(U)q,ℱ_q(, t):= p_t+12∇^2_q-2ln[2(p^2+q^2)](p^2+q^2)q-real(U)q-imag(U)p, ]
        
        where real(U) and imag(U) represent the real and imaginary parts of the external potential U(), respectively.
        Therefore, a fully-connected neural network NN(, t; W, B) with i hidden Layers and n neurons in every hidden layer can be constructed, where initialized parameters W = {w_j}_1^i+1 and B = {b_j}_1^i+1 being the weights and bias. Then, by the given activation function σ, one can obtain the expression in the form
        
            A_j=σ(w_j· A_j-1+b_j),
        
        where the w_j is a dim(A_j)×dim(A_j-1) matrix, A_0, A_i+1, b_i+1∈ℝ^2 and A_j, b_j∈ℝ^n.
        
        
        
        
        Furthermore, a Python library for PINNs, DeepXDE, was designed to serve a research tool for solving problems in computational science and engineering <cit.>. Using DeepXDE, we can conveniently define the physics-informed neural network ℱ(, z) as
        
        import deepxde as dde
        def pde(x, psi):
            p = psi[:, 0:1]
            q = psi[:, 1:2]
            p_xx = dde.grad.hessian(psi, x, component=0, i=0, j=0)
            q_xx = dde.grad.hessian(psi, x, component=1, i=0, j=0)
            p_yy = dde.grad.hessian(psi, x, component=0, i=1, j=1)
            q_yy = dde.grad.hessian(psi, x, component=1, i=1, j=1)
            p_t = dde.grad.jacobian(psi, x, i=0, j=2)
            q_t = dde.grad.jacobian(psi, x, i=1, j=2)
            F_p = -q_t + 0.5*(p_xx+p_yy) - 2*tf.log(2*(p**2+q**2))*(p**2+q**2)*p - (V*p - W*q)
            F_q = p_t + 0.5*(q_xx+q_yy) - 2*tf.log(2*(p**2+q**2))*(p**2+q**2)*q - (V*q + W*p)
            return [F_p, F_q]
        
        
        In order to train the neural network to fit the solutions of Eq. (<ref>), the total mean squared error (MSE) is defined as the following loss function containing three parts
        
            ℒ_0=MSE_F+MSE_I+MSE_B,
        
        with
        
            [ MSE_F=1/N_f∑_ℓ=1^N_f(|ℱ_p(_f^ℓ,t_f^ℓ)|^2
                 +|ℱ_q(_f^ℓ,t_f^ℓ)|^2),   MSE_I=1/N_I∑_ℓ=1^N_I(|p(_I^ℓ,0)-p_0^ℓ|^2
                 +|q(_I^ℓ,0)-q_0^ℓ|^2),   MSE_B=1/N_B∑_ℓ=1^N_B(|p(_B^ℓ,t_B^ℓ)|^2
                 +|q(_B^ℓ,t_B^ℓ)|^2), ]
        
        where {_f^ℓ,t_f^ℓ}_ℓ^N_f are connected with the marked points in Ω×[0,T] for the PINNs ℱ(, t)=ℱ_p(, t)+iℱ_q(, t),
        {_I^ℓ,p_0^ℓ,q_0^ℓ}_ℓ^N_I represent the initial data with ϕ(_I^ℓ)=p_0^ℓ+iq_0^ℓ, and {_B^ℓ,t_B^ℓ}_ℓ^N_B are linked with the randomly selected boundary training points in domain ∂Ω×[0,T].
        
        And then, we choose a hyperbolic tangent function tanh(·) as the activation function (of course one can also choose other nonlinear functions as the activation functions), and use Glorot normal to initialize variate. Therefore, the fully connected neural network can be written in Python as follows
        
        data = dde.data.TimePDE(
            geomtime, pde,
            initial-boundary value conditions,
            N_f, N_B, N_I,
        )
        net = dde.maps.FNN(layer_size, "tanh", Glorot normal)
        model = dde.Model(data, net)
        
        And then, with the aid of some optimization approaches (e.g., Adam & L-BFGS) <cit.>, we minimize the whole MSE ℒ_0 to make the approximated solution satisfy Eq. (<ref>) and initial-boundary value conditions.
        
        model.compile("adam", lr=1.0e-3)
        model.train(epochs)
        model.compile("L-BFGS")
        model.train()
        
        
        Therefore, for the given initial condition ϕ() solved by the IINNs, we can use PINNs to obtain solutions for the whole space-time region.
        
        Therefore, the main steps by the combination of the IINN and PINNs deep learning methods to solve the amended GP equation (<ref>) with potentials and the initial-boundary value conditions are presented as follows:
        
          1) Given an initial value that is sufficiently close to the stationary QDs we want to obtain. And a fully connected network NN_1 is trained to fit it.
        Then the IINN method is used to solve Eq. (<ref>).
        
        
          2) We initialize the network parameters of the second network NN_2 with the learned weights and biases from
        the first network NN_1, that is θ_0=argmin ℒ_1(θ). And train the NN_2 by minimizing the loss function ℒ_2 in terms of the optimization algorithm.
        
        
          3) Construct a fully-connected neural network NN(, t; θ) with randomly initialized parameters, and the PINNs ℱ(, t) is given by Eq. (<ref>).
        
        
          4) Generate the training data sets for the initial value condition given by IINN method, and considered model respectively from the initial boundary and within the region.
        
        
          5) Construct a training loss function ℒ_0 given by Eq. (<ref>) by summing the MSE of both the ℱ(, t) and initial-boundary value residuals.
        And train the NN to optimize the parameters θ={W, B} by minimizing the loss function in terms of
        the Adam & L-BFGS optimization algorithm.
        
        
        
        
        
        In what follows, the deep learning scheme is used to investigate the data-driven QDs of the 2D amended GP equation (<ref>) with two types of potential (quadruple-well Gaussian potential and -symmetric HOG potential).
        
        
        
        
        
        
        
        § DATA-DRIVEN 2D QDS IN AMENDED GP EQUATION WITH POTENTIALS
        
        
         §.§ Data-driven QDs in amended GP equation with quadruple-well Gaussian potential
        
        
        Firstly, we consider the 2D quadruple-well Gaussian potential <cit.>
            U()= V_0∑_j=1^4exp[-k(-_j)^2 ],   V_0<0,   k>0,
        
        where _j=(± x_0,± y_0), j=1,2,3,4 control the locations of these four potential wells, and |V_0| and k regulate the depths and widths of potential wells, respectively. Recently, based on the usual numerical methods, the spontaneous symmetry breaking (SSB) of 2D QDs was considered for the amended GP equation with the potential (<ref>) <cit.>, in which the complete pitchfork symmetry breaking bifurcation diagrams were presented for the possible stationary states with four modes, which involve twelve different real solution branches and one complex solution branches (for the complex one, the norm N=∫_ℝ^2|ϕ|^2d^2 r is the same as for one real branch),
        see Fig. <ref> for diagrams about the norm as a function of μ, and stable/unstable modes.
        
        In the following, we use the deep learning method to consider the four branches, that is, branches A0, A1, A3 and A4 (see Fig. <ref>).
        It should be noted that for the same potential parameters and chemical potential μ, Eq. (<ref>) can admit different solutions, which cannot be solved by general deep learning methods.
        Here we take potential parameters as V_0 =-0.5 and k=0.1, and consider Ω = [-12, 12]×[-12,12], T = 5 and μ=-0.5.
        If not otherwise specified, we choose a 4-hidden-layer deep neural network with 100 neurons per layer, and set learning rate α = 0.001.
        
        
        
        
        
        
        Čǎšě ̌1̌.̌—In branch A0, we firstly obtain the stationary QDs by the IINN method. We set N = 20000, and take the initial value as
        
            ϕ_0=∑_j=1^4a_jexp[-k(-_j)^2 ],
        
        where a_j=0.46  (j=1,2,3,4) and k=0.1.
        Through the IINN method, the learned QDs can be obtained at μ=-0.5, whose intensity diagram |ϕ()| and 3D profile
        are shown in Figs. <ref>(a1, a2), after 20000 steps of iterations with NN_1 and 3000 steps of iterations with NN_2. The relative L_2 error is 8.255472e-03 compared to the exact solution (numerically obtained). And the module of absolute error is exhibited in Fig. <ref>(a3). The loss-iteration plot of NN_1 is displayed in Fig. <ref>(a1).
        
        Then according to PINNs method, we take random sample points N_f=20000, N_B=150 and N_I=1000, respectively. Then, by using 40000 steps Adam and 10000 steps L-BFGS optimizations, we obtain the learned QDs solution ψ̂(,t) in the whole space-time region.
        Figs. <ref>(b1, b2, b3) exhibit the magnitude of the predicted solution at different time t = 0,  2.5, and 5.0, respectively. And the initial state (ϕ()=ψ(,t=0)) of the learned solution by IINN method and PINNs method is shown in Figs. <ref>(c1, c2), respectively. Furthermore, nonlinear propagation simulation of the learned 2D QDs is displayed by the isosurface of learned soliton at values 0.1, 0.5 and 0.9 hereinafter (see Fig. <ref>(c3)). The relative L_2 norm errors of ψ(, t), p(, t) and q(, t), respectively, are 1.952e-02, 1.396e-02 and 1.061e-02. And the loss-iteration plot is displayed in Fig. <ref>(a2).
        
        We should mention that the training stops in each step of the L-BFGS optimization when
        
            L_k-L_k+1/max{|L_k|,|L_k+1|,1}≤ np.finfo(float).eps,
        
        where L_k denotes loss in the n-th step L-BFGS optimization, and np.finfo(float).eps represent Machine Epsilon. Here we always set the default float type to `float64'. When the relative error between L_k and L_k+1 is less than Machine Epsilon, the iteration stops.
        
        
        
        
        
        
        Čǎšě ̌2̌.̌—In branch A1, similarly we get the stationary QDs by IINN method. We set N = 20000, and take the initial value as
        
            ϕ_0=∑_j=1^4a_jexp[-k(-_j)^2 ],
        
        where a_1=0.46,  a_2=a_3=a_4=0 and k=0.1.
        According to the IINN method, the learned QDs can be obtained at μ=-0.5, whose intensity diagram |ϕ()| and 3D profile
        are shown in Figs. <ref>(a1, a2), after 10000 steps of iterations with NN_1 and 5000 steps of iterations with NN_2. The relative L_2 error is 8.821019e-03 compared to the exact solution. And the module of absolute error is exhibited in Fig. <ref>(a3). The loss-iteration plot of NN_1 is displayed in Fig. <ref>(b1).
        
        Then through the PINNs method, we take random sample points N_f=20000, N_B=150 and N_I=1000, respectively. Then, by using 40000 steps Adam and 10000 steps L-BFGS optimizations, we obtain the learned QDs solution ψ̂(,t) in the whole space-time region.
        Figs. <ref>(b1, b2, b3) exhibit the magnitude of the predicted solution at different time t = 0,  2.5, and 5.0, respectively. And the initial state (ϕ()=ψ(,t=0)) of the learned solution by IINN method and PINNs method is shown in Figs. <ref>(c1, c2), respectively. Besides, nonlinear propagation simulation of the learned 2D QDs is displayed by the isosurface of learned soliton (see Fig. <ref>(c3)). The relative L_2 norm errors of ψ(, t), p(, t) and q(, t), respectively, are 3.364e-02, 3.767e-02 and 4.309e-02. And the loss-iteration plot is displayed in Fig. <ref>(b2).
        
        
        
        
        
        Čǎšě ̌3̌.̌—In branch A3, we firstly obtain the stationary QDs by IINN method. We set N = 20000, and take the initial value as
        
            ϕ_0=∑_j=1^4a_jexp[-k(-_j)^2 ],
        
        where a_1=a_3=0.3,  a_2=a_4=0 and k=0.1.
        Through the IINN method, the learned QDs can be obtained at μ=-0.5, whose intensity diagram |ϕ()| and 3D profile
        are shown in Figs. <ref>(a1, a2), after 15000 steps of iterations with NN_1 and 5000 steps of iterations with NN_2. The relative L_2 error is 7.201800e-03 compared to the exact solution. And the module of absolute error is exhibited in Fig. <ref>(a3).
        
        Then according to PINNs method, we take random sample points N_f=20000, N_B=150 and N_I=1000, respectively. Then, by using 40000 steps Adam and 10000 steps L-BFGS optimizations, we obtain the learned QDs solution ψ̂(,t) in the whole space-time region.
        Figs. <ref>(b1, b2, b3) exhibit the magnitude of the predicted solution at different time t = 0,  2.5, and 5.0, respectively. And the initial state (ϕ()=ψ(,t=0)) of the learned solution by IINN method and PINNs method is shown in Figs. <ref>(c1, c2), respectively. Furthermore, nonlinear propagation simulation of the learned 2D QDs is displayed by the isosurface of learned soliton (see Fig. <ref>(c3)). The relative L_2 norm errors of ψ(, t), p(, t) and q(, t), respectively, are 2.002e-02, 2.458e-02 and 2.356e-02.
        
        
        Čǎšě ̌4̌.̌—In branch A4, we get the stationary QDs by IINN method. We set N = 20000, and take the initial value as
        
            ϕ_0=∑_j=1^4a_jexp[-k(-_j)^2 ],
        
        where a_1=a_2=a_3=0.46, a_4=0 and k=0.1.
        Through the IINN method, the learned QDs can be obtained at μ=-0.5, whose intensity diagram |ϕ()| and 3D profile
        are shown in Figs. <ref>(a1, a2), after 20000 steps of iterations with NN_1 and 3000 steps of iterations with NN_2. The relative L_2 error is 3.380430e-03 compared to the exact solution (numerically obtained). And the module of absolute error is exhibited in Fig. <ref>(a3).
        
        Then according to PINNs method, we take random sample points N_f=20000, N_B=150 and N_I=1000, respectively. Then, by using 40000 steps Adam and 10000 steps L-BFGS optimizations, we obtain the learned QDs solution ψ̂(,t) in the whole space-time region.
        Figs. <ref>(b1, b2, b3) exhibit the magnitude of the predicted solution at different time t = 0,  2.5, and 5.0, respectively. And the initial state (ϕ()=ψ(,t=0)) of the learned solution by IINN method and PINNs method is shown in Figs. <ref>(c1, c2), respectively. Furthermore, nonlinear propagation simulation of the learned 2D QDs is displayed by the isosurface of learned soliton at values 0.1, 0.5 and 0.9 (see Fig. <ref>(c3)). The relative L_2 norm errors of ψ(, t), p(, t) and q(, t), respectively, are 1.256e-02, 1.899e-02 and 1.499e-02.
        
        
        
        
        
         §.§ Data-driven QDs in amended GP equation with -symmetric HOG potential
        
        
        
        
        
        In this subsection, we consider the following -symmetric HO-Gaussian (HOG) potential with the real and imaginary parts being <cit.>
            V()=r^2( 1+e^-r^2) +V_0 (e^-2x^2+e^-2y^2),    
            W()=W_0( xe^-x^2+ye^-y^2) ,
        
        where r^2=x^2+y^2, the coefficient in front of the HO potential is set to be 1,
        the real parameter V_0 modulates the profile of the external potential V(), and real W_0 is the strength of gain-loss distribution W(). The vortex solitons were produced for a variety of 2D spinning QDs in the -symmetric potential, modeled by the amended GP equation with Lee–Huang–Yang corrections <cit.>, where the dependence of norm N on chemical potential μ was illustrated for different families of droplet modes in -symmetric HOG potential (see Fig. <ref>).
        
        In the following, we use the deep learning method to consider the multi-component QDs under different chemical potential.
        Here we take potential parameters as V_0 =-1/16 and W_0=1, and consider Ω = [-8, 8]×[-8,8], T = 3.
        Considering that the solution ϕ() of Eq. (<ref>) is a complex-valued function, similarly we set the network's output ϕ()=p()+iq() and then separate Eq. (<ref>) into its real and imaginary parts.
        
            [ ℱ_p():= -12∇^2_p+2ln[2(p^2+q^2)](p^2+q^2)p+real(U)p-imag(U)q-μ p,ℱ_q():= -12∇^2_q+2ln[2(p^2+q^2)](p^2+q^2)q+real(U)q+imag(U)p-μ q. ]
        
        Then the loss function ℒ_2 becomes
        
            ℒ_2:=1/N∑_i=1^N(|ℱ_p(_i)|^2+|ℱ_q(_i)|^2)/max_i(√((p(_i)^2+q(_i)^2)).
        Čǎšě ̌1̌.̌ ̌ ̌Q̌Ďš ̌w̌ǐťȟ ̌ťȟě ̌ǒňě-̌čǒm̌p̌ǒňěňť ̌šťřǔčťǔřě.̌—Firstly, we consider the -symmetric droplets with the simplest structure.
        We can obtain the initial conditions by computing the spectra and eigenmodes in the linear regime, which can be given as follows
        
            ℋΦ (𝐫)=λΦ (𝐫),  ℋ
            =-∇ _𝐫^2+U(),
        
        where λ and Φ (𝐫) are the eigenvalue and localized
        eigenfunction, respectively. The linear spectral problem (<ref>) can be solved
        numerically by dint of the Fourier spectral method <cit.>.
        
        We take the initial value as the linear mode Φ at ground state and N=10000.
        Through the IINN method, the learned QDs can be obtained at μ=2, after 10000 steps of iterations with NN_1 and 10000 steps of iterations with NN_2.
        Figs. <ref>(a1, a2, a3) exhibit the intensity diagram of real part, imaginary part and |ϕ()|.  The module of absolute error is shown in Fig. <ref>(a4).
        The relative L_2 errors of ϕ(), p() and q(), respectively, are 1.992564e-02, 5.547692e-02 and 2.075972e-02.
        
        Then according to PINNs method, we take random sample points N_f=20000, N_B=150 and N_I=1000, respectively. Then, by using 30000 steps Adam and 10000 steps L-BFGS optimizations, we obtain the learned  QDs solution ψ̂(,t) in the whole space-time region.
        Figs. <ref>(b1, b2, b3) exhibit the magnitude of the predicted solution at different time t = 0,  1.5, and 3.0, respectively. And, nonlinear propagation simulation of the learned 2D QDs is displayed by the isosurface of learned soliton at values 0.1, 0.5 and 0.9 hereinafter (see Fig. <ref>(b4)). The relative L_2 norm errors of ψ(, t), p(, t) and q(, t), respectively, are 2.312e-02, 1.995e-02 and 2.331e-02.
        
        
        
        
        
        
        Čǎšě ̌2̌.̌ ̌ ̌Q̌Ďš ̌w̌ǐťȟ ̌ťȟě ̌ťw̌ǒ-̌čǒm̌p̌ǒňěňť ̌šťřǔčťǔřě—Second, we consider the -symmetric droplets with the two-component structure.
        We take the initial value as the linear mode Φ at the first excited state and N=10000.
        Through the IINN method, the learned QDs can be obtained at μ=2.8, after 10000 steps of iterations with NN_1 and 10000 steps of iterations with NN_2.
        Figs. <ref>(a1, a2, a3) exhibit the intensity diagram of real part, imaginary part and |ϕ()|.  The module of absolute error is shown in Fig. <ref>(a4).
        The relative L_2 errors of ϕ(), p() and q(), respectively, are 5.231775e-02, 1.546516e-02 and 6.117475e-02.
        
        Then according to PINNs method, we take random sample points N_f=20000, N_B=150 and N_I=1000, respectively. Then, by using 30000 steps Adam and 15000 steps L-BFGS optimizations, we obtain the learned  QDs solution ψ̂(,t) in the whole space-time region.
        Figs. <ref>(b1, b2, b3) exhibit the magnitude of the predicted solution at different time t = 0,  1.5, and 3.0, respectively. And, nonlinear propagation simulation of the learned 2D QDs is displayed by the isosurface of learned soliton at values 0.1, 0.5 and 0.9 hereinafter (see Fig. <ref>(b4)). The relative L_2 norm errors of ψ(, t), p(, t) and q(, t), respectively, are 4.326e-02, 5.158e-02 and 5.182e-02.
        
        
        
        
        
        
        
        
        Čǎšě ̌3̌.̌ ̌ ̌Q̌Ďš ̌w̌ǐťȟ ̌ťȟě ̌ťȟřěě-̌čǒm̌p̌ǒňěňť ̌šťřǔčťǔřě—Then, we consider the -symmetric droplets with the three-component structure.
        We take the initial value as the linear mode Φ at the second excited state and N=10000.
        Through the IINN method, the learned QDs can be obtained at μ=4.3, after 10000 steps of iterations with NN_1 and 10000 steps of iterations with NN_2.
        Figs. <ref>(a1, a2, a3) exhibit the intensity diagram of real part, imaginary part and |ϕ()|.  The module of absolute error is shown in Fig. <ref>(a4).
        The relative L_2 errors of ϕ(), p() and q(), respectively, are 3.606203e-02, 5.148407e-02 and 4.234037e-02.
        
        Then according to PINNs method, we take random sample points N_f=20000, N_B=150 and N_I=1000, respectively. Then, by using 30000 steps Adam and 15000 steps L-BFGS optimizations, we obtain the learned  QDs solution ψ̂(,t) in the whole space-time region.
        Figs. <ref>(b1, b2, b3) exhibit the magnitude of the predicted solution at different time t = 0,  1.5, and 3.0, respectively. And, nonlinear propagation simulation of the learned 2D QDs is displayed by the isosurface of learned soliton at values 0.1, 0.5 and 0.9 hereinafter (see Fig. <ref>(b4)). The relative L_2 norm errors of ψ(, t), p(, t) and q(, t), respectively, are 5.410e-02, 6.771e-02 and 6.189e-02.
        
        
        
        
        
        Čǎšě ̌4̌.̌ ̌ ̌Q̌Ďš ̌ ̌w̌ǐťȟ ̌ťȟě ̌f̌ǒǔř-̌čǒm̌p̌ǒňěňť ̌šťřǔčťǔřě—Finally, we consider the -symmetric droplets with the four-component structure.
        We take the initial value as the linear mode Φ at the three excited state and N=10000.
        Through the IINN method, the learned QDs can be obtained at μ=4.2, after 10000 steps of iterations with NN_1 and 20000 steps of iterations with NN_2.
        Figs. <ref>(a1, a2, a3) exhibit the intensity diagram of real part, imaginary part and |ϕ()|.  The module of absolute error is shown in Fig. <ref>(a4).
        The relative L_2 errors of ϕ(), p() and q(), respectively, are 1.186940e-02, 4.332541e-02 and 1.185871e-02.
        
        Then according to PINNs method, we take random sample points N_f=20000, N_B=150 and N_I=1000, respectively. Then, by using 30000 steps Adam and 15000 steps L-BFGS optimizations, we obtain the learned  QDs solution ψ̂(,t) in the whole space-time region.
        Figs. <ref>(b1, b2, b3) exhibit the magnitude of the predicted solution at different time t = 0,  1.5, and 3.0, respectively. And, nonlinear propagation simulation of the learned 2D QDs is displayed by the isosurface of learned soliton at values 0.1, 0.5 and 0.9 hereinafter (see Fig. <ref>(b4)). The relative L_2 norm errors of ψ(, t), p(, t) and q(, t), respectively, are 5.128e-02, 5.841e-02 and 5.669e-02.
        
        
        
         Řěm̌ǎřǩ.̌ It should be noted that the systematic results shown in Figs. <ref>and <ref> have been investigated by numerical methods in Refs. <cit.>. In this paper, we mainly consider partial solutions and their short evolutions by using the machine learning method. For their stability, in general it can be solved by linear eigenvalue problems and long time evolution.   However, solving the eigenvalue problem via machine learning methods may be more difficult because it also involves multi-solution problems.  This will be our future work to consider.
          Furthermore, many methods already exist to return stable and accurate predictions across long temporal horizons by training multiple individual networks in different temporal sub-domains <cit.>. These approaches inevitably lead to a large computational cost and more complex network structure. We use the parallel PINNs to investigate the longer-time evolutions for both solutions via the domain decomposition (see Fig. <ref>). We can see that the solutions are still stable after the longer-time evolutions,  compared with Fig. <ref>(c3) and Fig. <ref>(b4), respectively. 
        
        
        § CONCLUSIONS AND DISCUSSIONS
        
        
        In conclusion, we have investigated the 2D stationary QDs and their evolutions in amended Gross–Pitaevskii equation with potentials via deep learning neural networks. Firstly, we use the IINN method for learning 2D stationary QDs. Then the learned 2D stationary QDs are used as the initial-value conditions for PINNs to display  their evolutions in the some space-time regions. Especially, we consider two types of potentials, one is 2D quadruple-well Gaussian potential and the other is -symmetric HO-Gaussian potential, which lead to spontaneous symmetry breaking and the generation of multi-component QDs.
        
        On the other hand, in order to study the stability of the QDs, we can use deep learning methods to study the interactions between droplets. Furthermore, we can investigate the spinning QDs in terms of spinning coordinates, x^'=xcos (ω t)+ysin (ω t), y^'=ycos (ω t)-xsin (ω t) with angular velocity ω. We will investigate these issues in future.
        
        
         ̌̌ AcknowledgementŤhe work  was supported by the National Natural Science Foundation of China  under Grant No. 11925108.
        
        
        991 D. S. Petrov, Quantum mechanical stabilization of a collapsing Bose-Bose mixture, Phys. Rev. Lett. 115 (2015) 155302.
        
        
        4 G. E. Astrakharchik, B. A. Malomed, Dynamics of one-dimensional quantum droplets, Phys. Rev. A 98 (2018) 013631.
        
        3 L. Chomaz, S. Baier, D. Petter, M. J. Mark, F. Wachtler, L. Santos, F. Ferlaino, Quantum-fluctuation-driven crossover from a dilute Bose–Einstein
        condensate to a macrodroplet in a dipolar quantum fluid, Phys. Rev. X 6 (2016) 041039.
        
        6 E. Shamriz, Z. Chen, B. A. Malomed, Suppression of the quasitwo-dimensional quantum collapse in the attraction field by the
        Lee-Huang-Yang effect, Phys. Rev. A 101 (2020) 063628.
        
        5 Y. V. Kartashov, B. A. Malomed, L. Torner, Metastability of quantum droplet clusters, Phys. Rev. Lett. 122 (2019) 193902.
        
        8 Z. Luo, W. Pang, B. Liu, Y. Li, B. A. Malomed, A new form of liquid matter: Quantum droplets, Front. Phys. 16 (2021) 32201.
        
        
        7 M. Tylutki, G. E. Astrakharchik, B. A. Malomed, D. S. Petrov, Collective excitations of a one-dimensional quantum droplet, Phys. Rev. A 101 (2020)
        051601(R).
        
        10 Y. V. Kartashov, B. A. Malomed, L. Torner, Structured heterosymmetric quantum droplets, Phys. Rev. Res. 2 (2020) 033522.
        
        
        9 Y. V. Kartashov, V. V. Konotop, D. A. Zezyulin, L. Torner, Bloch oscillations in optical and Zeeman lattices in the presence of spin–orbit coupling, Phys.
        Rev. Lett. 117 (2016) 215301.
        
        
        11 C. Cabrera, L. Tanzi, J. Sanz, B. Naylor, P. Thomas, P. Cheiney, L. Tarruell,
        Quantum liquid droplets in a mixture of Bose–Einstein condensates, Science 359 (2018) 301.
        
        12 P. Cheiney, C. R. Cabrera, J. Sanz, B. Naylor, L. Tanzi, L. Tarruell, Bright
        soliton to quantum droplet transition in a mixture of Bose–Einstein condensates, Phys. Rev. Lett. 120 (2018) 135301.
        
        13 I. Ferrier-Barbut, H. Kadau, M. Schmitt, M. Wenzel, T. Pfau, Observation
        of quantum droplets in a strongly dipolar Bose gas, Phys. Rev. Lett. 116 (2016) 215301.
        
        
        14 D. Edler, C. Mishra, F. Wächtler, R. Nath, S. Sinha, L. Santos, Quantum
        fuctuations in quasi-one-dimensional dipolar Bose–Einstein condensates, Phys. Rev. Lett. 119 (2017) 050403.
        
        15 T. D. Lee, K. Huang, C. N. Yang, Eigenvalues and eigenfunctions of a Bose
        system of hard spheres and its lowtemperature properties, Phys. Rev. 106 (1957) 1135.
        
        ob1 H. Kadau, M. Schmitt, M. Wentzel, C. Wink, T. Maier, I. Ferrier-Barbut, and T. Pfau, Observing the Rosenzweig instability of a quantum ferrofluid, Nature 530,
        194 (2016).
        
        ob2 M. Schmitt, M. Wenzel, F. Büttcher, I. Ferrier-Barbut, and T. Pfau, Self-bound droplets of a dilute magnetic quantum liquid, Nature 539, 259 (2016).
        
        
        22 V. Cikojević, L. V. Markić, G. E. Astrakharchik, J. Boronat, Universality in ultradilute liquid Bose-Bose mixtures, Phys. Rev. A 99 (2019) 023618.
        
        23 V. Cikojević, K. Dzelalija, P. Stipanovic, L. V. Markić, J. Boronat, Ultradilute quantum liquid drops, Phys. Rev. B 97 (2018) 140502(R).
        
        
        Li G. Li, X. Jiang, B. Liu, Z. Chen, B. A. Malomed, and Y. Li, Two-dimensional anisotropic vortex quantum droplets in dipolar Bose-Einstein condensates. Front. Phys., 19 (2024) 22202.
        
        
        liyy  Y. Li, Z. Chen, Z. Luo, C. Huang, H. Tan, W. Pang, B.A. Malomed, Two-dimensional vortex quantum droplets, Phys. Rev. A 98 (2018) 063602.
        
        25  M.N. Tengstrand, P. Stürmer, E.Ö. Karabulut, S.M. Reimann, Rotating binary
        Bose–Einstein condensates and vortex clusters in quantum droplets, Phys. Rev. Lett. 123 (2019) 160405.
        
        26 Z. Zhou, X. Yu, Y. Zou, H. Zhong, Dynamics of quantum droplets in a
        one-dimensional optical lattice, Commun. Nonlinear Sci. Numer. Simul. 78 (2019) 104881.
        
        27 B. Liu, H. Zhang, R. Zhong, X. Zhang, X. Qin, C. Huang, Y. Li, B.A. Malomed,
        Symmetry breaking of quantum droplets in a dual-core trap, Phys. Rev. A 99 (2019) 053602.
        
        28 Z. Zhou, B. Zhu, H. Wang, H. Zhong, Stability and collisions of quantum
        droplets in -symmetric dual-core couplers, Commun. Nonlinear Sci. Numer. Simul. 91 (2020) 105424.
        
        29 J. Song, Z. Yan, Dynamics of 1D and 3D quantum droplets in
        parity-time-symmetric harmonic-Gaussian potentials with two competing nonlinearities, Physica D 442 (2022) 133527.
        
        
        
        
        pinn M. Raissi, P. Perdikaris, G. E. Karniadakis. Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations, J. Comput. Phys. 378 (2019) 686.
        
        pinn1 S. Goswami, C. Anitescu, S. Chakraborty, T. Rabczuk, Transfer learning enhanced
        physics informed neural network for phase-field modeling of fracture, Theor. Appl. Fract. Mech. 106 (2020) 102447.
        
        pinn2 A. D. Jagtap, K. Kawaguchi, G. E. Karniadakis, Adaptive activation functions accelerate convergence in deep and physics-informed neural networks, J. Comput. Phys.
        404 (2020) 109136.
        
        pinn3 X. Meng, Z. Li, D. Zhang, G. E. Karniadakis, PPINN: parareal physics-informed neural network for time-dependent PDEs, J. Comput. Phys. 370 (2020) 1132.
        
        
        deepxde L. Lu, X. Meng, Z. Mao, G.E. Karniadakis, DeepXDE: a deep learning library for solving differential equations, SIAM Rev. 63 (2021) 208–228.
        
        
        deepritz W. E, B. Yu, The deep Ritz method: a deep learning-based numerical algorithm for
        solving variational problems, Commun. Math. Stat. 6 (2018) 1–12.
        
        pnet  Z. Long, Y. Lu, X. Ma, B. Dong, PDE-net: learning PDEs from data, in: Proceedings
        of the 35th International Conference on Machine Learning, in: PMLR, vol. 80, 2018,
        pp. 3208–3216.
        pnet1 Z. Long, Y. Lu, B. Dong, PDE-net 2.0: learning PDEs from data with a numericsymbolic hybrid deep network, J. Comput. Phys. 399 (2019) 108925.
        
        twostage S. Lin, Y. Chen, A two-stage physics-informed neural network method based on
        conserved quantities and applications in localized wave solutions, J. Comput. Phys. 457 (2022) 111053.
        
        pinndq J. Pu, Y. Chen, Complex dynamics on the one-dimensional quantum droplets via
        time piecewise PINNs, Physica D 454 (2023) 133851.
        
        
        yan1 L. Wang, Z. Yan, Data-driven peakon and periodic peakon solutions and parameter discovery of some nonlinear dispersive equations via deep learning, Physica D 428 (2021) 133037.
        
        yan2  J. Song, Z. Yan, Deep learning soliton dynamics and complex potentials recognition for 1D and 2D -symmetric saturable nonlinear Schrödinger equations, Physica D 448 (2023) 133729.
        
        yan3 Z. Zhou, Z. Yan, Solving forward and inverse problems of the logarithmic nonlinear Schrödinger equation with -symmetric harmonic potential via deep learning,
        Phys. Lett. A 387 (2021) 127010.
        
        yan4 L. Wang, Z. Yan, Data-driven rogue waves and parameter discovery in the defocusing
        NLS equation with a potential using the PINN deep learning, Phys. Lett. A 404 (2021) 127408.
        
        yan5 M. Zhong, S. Gong, S.-F. Tian,  Z. Yan, Data-driven rogue waves and parameters discovery in nearly
        integrable PT-symmetric Gross–Pitaevskii equations via PINNs deep learning,  Physica D 439 (2022) 133430.
        
        yan6 Z. Zhou, Z. Yan, Is the neural tangent kernel of PINNs deep learning general partial
        differential equations always convergent? Physica D 457 (2024) 133987.
        
        chen1  H. Zhou, J. Pu, Y. Chen, Data-driven forward–inverse problems for the variable
        coefficients Hirota equation using deep learning method, Nonlinear Dyn. 111 (2023) 14667–14693.
        
        li1 J.H. Li, B. Li, Mix-training physics-informed neural networks for the rogue waves of nonlinear Schrödinger equation, Chaos, Solitons and Fractals 164 (2022) 112712.
        
        pet T. I. Lakoba and J. Yang, A generalized Petviashvili iteration method forscalar and vector Hamiltonian equations with arbitrary form of
        nonlinearity, J. Comp. Phys. 226 (2007) 1668-1692.
        
        it  J. Yang and T. I. Lakoba, Accelerated imaginary-time evolution methods for the computation of solitary waves, Stud. Appl. Math. 120 (2008)
        265-292.
        
        49 J. Yang, Newton-conjugate-gradient methods for solitary wave computations,
        J. Comput. Phys. 228 (2009) 7007–7024.
        
        51 J. Yang and T. I. Lakoba, Universally-convergent squared-operator iteration
        methods for solitary waves in general nonlinear wave equations, Stud. Appl. Math.
        118 (2007) 153–197.
        
        52 J. Yang, Nonlinear Waves in Integrable and
        Nonintegrable Systems (SIAM, 2010).
        
        
        IINN J. Song, M. Zhong, G. E. Karniadakis, and Z. Yan, Two-stage initial-value iterative physics-informed neural networks
        for simulating solitary waves of nonlinear wave equations, J. Comput. Phys.  505 (2024) 112917.
        
        
        
        fourwell J. Song, H. Dong, D. Mihalache,  Z. Yan, Spontaneous symmetry breaking, stability and adiabatic changes of 2D
        quantum droplets in amended Gross–Pitaevskii equation with multi-well potential, Physica D 448 (2023) 133732.
        
        spin J. Song, Z. Yan,  B. A. Malomed, Formations and dynamics of two-dimensional spinning asymmetric quantum
        droplets controlled by a -symmetric potential,  Chaos 33 (2023) 033141.
        
        50 M. J. Ablowitz, Z. H. Musslimani, Spectral renormalization method for computing self-localized solutions to nonlinear systems, Opt. Lett. 30 (2005) 2140–2142.
        
        
        
        adam D. Kingma, J. Ba, Adam: a method for stochastic optimization, 2014, arXiv:1412.6980.
        
        bfgs D.C. Liu, J. Nocedal, On the limited memory BFGS method for large scale optimization, Math. Program. 45 (1989) 503–528.
        
        
        xPPINN K. Shukla, A. Jagtap, G. E. Karniadakis, Parallel physics-informed neural networks via domain decomposition. J Comput Phys. 447 (2021) 110683.
        PPINN X. Meng, Z. Li, D. Zhang,  G. E. Karniadakis, PPINN: Parareal physics-informed neural
        network for time-dependent PDEs, Comput. Methods Appl. Mech. Eng. 370 (2020) 113250.
        ednn Y. Du, T. Zaki, Evolutional deep neural network, Phys. Rev. E, 104 (2021) 045303.
        longdeeponet S. Wang, P. Perdikaris, Long-time integration of parametric evolution equations with
        physics-informed DeepONets, J. Comput. Phys. 475 (2023) 111855.
