An improved fundamental invariant neural network (FI-NN) approach for representing a potential energy surface (PES) involving permutation symmetry is introduced in this work. In this approach, FIs are regarded as symmetric neurons, thus avoiding complex preprocessing of training set data, especially when the training set contains gradient data. In this work, the improved FI-NN method, combined with simultaneous fitting of the energy and gradient strategy, is used for constructing a global accurate PES of a Li2Na system (root-mean-square error of 12.20 cm-1). The potential energies and the corresponding gradients are calculated by a UCCSD(T) method with effective core potentials. Based on the new PES, the vibrational energy levels and the corresponding wave functions of Li2Na molecules are calculated using an accurate quantum mechanics method. To accurately describe the cold or ultracold reaction dynamics of the Li + LiNa(v = 0, j = 0) → Li2(v', j') + Na reaction, the long-range region of the PES in both the reactant and product asymptotes is represented by an asymptotically correct form. A statistical quantum model (SQM) is used to study the dynamics of the ultracold Li + LiNa reaction. The calculated results are in good agreement with the exact quantum dynamics results (B. K. Kendrick, J. Chem. Phys., 2021, 154, 124303), which indicates that the dynamics of the ultracold Li + LiNa reaction can be well described by the SQM approach. The time-dependent wave packet calculations are performed for the Li + LiNa reaction at thermal energies, and the characteristic of differential cross-sections confirms that the reaction follows the complex-forming reaction mechanism.