Heterogeneous expert partial differential equation solving method based on time sequence causality and teacher-student constraints

Abstract

This invention belongs to the field of physics field solving, specifically disclosing a heterogeneous expert partial differential equation (PDE) solution method based on temporal causality and teacher-student constraints. By integrating three core technologies—hybrid expert (MoE), Kolmogorov-Arnold network (KAN), and curriculum learning (CL)—and introducing teacher-student consistency constraints, residual adaptive resampling (RAR), and a two-stage hybrid optimization strategy, a hybrid expert physics information neural network (MoE-KAN-PINN) is constructed to achieve high-precision and high-stability solutions for long-term, strongly nonlinear PDEs. This invention achieves several key technological breakthroughs in neural network architecture design, physical constraint optimization strategies, and long-term temporal progression mechanisms, effectively overcoming the core technical problems faced by existing PINN methods in solving long-term, strongly nonlinear PDEs, such as causal failure, catastrophic forgetting, insufficient local feature capture ability, low solution accuracy, and poor optimization efficiency.

Description

Technical Field

[0001] This invention relates to the field of physical field solution technology, specifically to a method for solving heterogeneous expert partial differential equations based on temporal causality and teacher-student constraints. Background Technology

[0002] In the fields of deep learning-driven scientific computing and physics modeling, solving partial differential equations (PDEs) places stringent demands on the function approximation accuracy, high-order derivative computation capabilities, and parameter utilization efficiency of neural networks. Traditional multi-layer perceptrons (MLPs) rely on fixed linear weight matrices and globally shared nonlinear activation functions to complete feature mapping. When dealing with complex continuous function fitting, high-dimensional physical field simulation, and high-order partial derivative solving, they are prone to numerical rigidity and parameter redundancy, making it difficult to achieve a balance between accuracy and computational efficiency.

[0003] Physics-Informed Neural Networks (PINNs) are PDE solving methods that integrate physical constraints and deep learning capabilities. By embedding physical information such as governing equations and boundary conditions into the loss function, they achieve end-to-end solutions, becoming an important method in the field of scientific computing. Time-marching Physics-Informed Neural Networks (Time-marching PINNs) are a sequential training strategy based on PINNs. The core of this approach is to decompose long-term PDE solving tasks into multiple continuous time sub-intervals for sequential solving, in order to follow the temporal causality of the physical system. Although time-marching methods achieve forward information transmission at a macroscopic scale by discretizing the time domain into independent sub-intervals, they are essentially a fragmented piecewise modeling strategy with the following inherent drawbacks: 1) Single-point dependence and error accumulation in information transmission: Information transmission between adjacent sub-intervals depends only on the single-point predicted value at the end time of the previous interval. This fragile single-point propagation mechanism carries significant risks: if a small error exists at a single point (e.g., due to local oscillations or numerical noise), subsequent intervals will use this error as an initial condition for training, causing the error to amplify exponentially along the time direction. Especially for chaotic systems or long-term problems, this error accumulation effect can rapidly destroy the entire solution process, causing the model training to diverge or converge to a completely incorrect physical solution. 2) Catastrophic forgetting and solution inconsistency: The training of the subsequent sub-interval is completely independent of the previous interval. Once the previous interval is trained, it is frozen and no longer participates in subsequent optimization. This means that while the model is learning a new time domain, it cannot make any corrections or maintain the old time domains it has already learned. Since the evolution of the physical field is continuous, this "learning new and forgetting old" mechanism is very likely to cause discontinuities or abrupt changes in the global solution along the time direction, destroying the smoothness and consistency that the physical field should have. 3) Artificial discontinuities introduced by boundary stitching: Sub-intervals are only stitched together by soft constraints through the solution continuity (C^0) condition, lacking mandatory requirements for the continuity of higher-order derivatives (such as flux, gradient). For convection-dominated problems, shock propagation, or abrupt interface changes, the first and higher derivatives of the physical fields often need to be continuous. Simply satisfying C^0 in a stitching method introduces spurious numerical oscillations and non-physical discontinuities at the time interface, further deteriorating the solution accuracy.

[0004] Causal Weights Physics-Informed Neural Networks (CausalWeights PINNs) is a PINN training strategy that respects causality on continuous time scales by modifying the loss function without performing temporal domain decomposition. The core idea is to introduce a weight factor that decays exponentially with the cumulative residual for the PDE residuals at different times, forcing the model to learn the evolution of the physical field sequentially from the optimization objective level. However, this method still requires calculating the residuals and automatic differentiation of all collocation points in the entire spatiotemporal domain in each iteration, without any reduction or compression of the computational domain. Its memory usage and computational cost are comparable to standard PINNs. More importantly, this method is a single-network global strategy, lacking the ability to combine with domain decomposition or multi-expert architectures, and is ill-suited to the spatial non-uniformity and multi-scale characteristics of the physical field.

[0005] Adaptive time-sampling (PINNs) is a PINN training strategy that indirectly guides the model's learning order by dynamically adjusting the distribution of training data. Its core principle is to mimic the "from easy to difficult" learning logic of a course, gradually expanding the time sampling domain to avoid exposing the model to complex samples in the early stages of training. Essentially, this method is a data sampling initiation strategy and does not fundamentally change the mathematical form of the loss function or the optimization objective. It suffers from insurmountable technical flaws: once all time slices' collocations are finally included in the training set, its loss function is completely consistent with standard PINN, meaning that the risk of causal failure still exists in the later stages of training. Although the introduction order of time slices has a temporal relationship, the collocations within each time slice are still optimized simultaneously during the optimization process. The residuals of later time slices are not suppressed by the insufficient optimization of the residuals of earlier time slices, thus failing to guarantee a strict information propagation order. Furthermore, the time slice division and advancement step size of this method are usually preset fixed values, lacking the ability to adapt to the speed of physical field evolution. In the phase where the solution changes slowly, advancing too quickly may introduce redundant computations; in the phase where the solution changes drastically, advancing too slowly may waste computational resources. The fixed nature of the sampling strategy also makes it difficult to achieve adaptive adjustment of collocation density in regions with drastic local changes, such as shock waves.

[0006] Therefore, there is an urgent need for a solution that can address the technical challenges of existing methods in solving long-term, highly nonlinear PDEs, and provide a more efficient and accurate deep learning solution for simulating complex physical fields in scientific computing. Summary of the Invention

[0007] To address the problems existing in the prior art, this invention provides a heterogeneous expert partial differential equation (PDE) solution method based on temporal causality and teacher-student constraints. It integrates three core technologies: hybrid expert (MoE), Kolmogorov-Arnold network (KAN), and curriculum learning (CL). It introduces teacher-student consistency constraints, residual adaptive resampling (RAR), and a two-stage hybrid optimization strategy to construct a hybrid expert physical information neural network (MoE-KAN-PINN). This achieves high-precision and high-stability solutions for long-term, highly nonlinear PDEs, thus solving the problems mentioned in the background art.

[0008] To achieve the above objectives, the present invention provides the following technical solution: a method for solving heterogeneous expert partial differential equations based on temporal causality and teacher-student constraints, comprising the following steps: S1. Spatiotemporal domain construction and full dataset preloading; S2. Construct a hybrid expert physical information neural network MoE-KAN-PINN and a dual network architecture for teachers and students; S3. Construction of dynamic staged total loss function; S4. Execute a two-stage hybrid optimization strategy to train the MoE-KAN-PINN model.

[0009] S5. Generation of physical field prediction solutions.

[0010] Preferably, step S1 specifically includes the following: S11. Definition of General Spatiotemporal Domain: Let the spatial solution domain of the target PDE be... d is the spatial dimension, and the total time domain is Where T is the maximum solution time as defined by the physical problem; S12. Full Dataset Generation and Preloading: Latin hypercube sampling (LHS) is used to generate PDE internal collocations. The internal PDE collocation set is as follows: In all time and space domains Internal LHS sampling is performed, and the samples are immediately sorted globally in ascending order according to the time coordinate t. S13. Dynamic Adaptive Data Slicing Rules: Based on the time-domain expansion strategy of course learning, a predefined generalized time segmentation plan is used to expand the total time domain. Divided into Each progressive training phase, with a single phase time step. Equal step size settings can be used; for the first Each training phase automatically completes two levels of data slicing; specifically, it includes training set selection and replay set selection. Training set selection: extraction Boundary condition points and PDE collocation points within the range are used as training data for the current stage; Indicates the current maximum training time; Playback set filtering: Extraction PDE collocations within the range, used as the replay dataset for teacher-student consistency constraints, are only used within... Effective immediately; This indicates the maximum time that the preceding sequence has been learned.

[0011] Preferably, step S2 specifically includes the following: S21. Construct a heterogeneous expert library consisting of four parallel heterogeneous KAN expert sub-networks and one FCN expert. Each expert network uses normalized spatiotemporal coordinates (x,t) as input and maps them to intermediate feature maps with consistent dimensions through unique mathematical bases or activation functions to achieve adaptive expression of multi-scale and non-uniform features of physical fields. The four parallel heterogeneous KAN expert sub-networks include: Wav-KAN, RBF-KAN, Chebyshev-KAN, and Fourier-KAN. S22. The gated network adopts a multi-layer fully connected structure with two hidden layers. Each layer has 64 neurons. The activation function is Tanh. The input is the normalized spatiotemporal coordinates (x,t). The output layer is a Softmax layer. Finally, a probability vector with a dimension of 5 is obtained. This vector is the activation weights of the five expert networks and satisfies the normalization condition that the sum of the weights is 1. S23. Definition of Teacher-Student Dual Network Architecture: Student Network: The trainable MoE-KAN-PINN model in the current training phase. During the k-th training phase, the model weights converged in the (k-1)-th phase are loaded as initialization parameters and participate in gradient updates and parameter optimization throughout the process. Teacher Network: Enabled only when k≥2, it completely freezes the model parameters that have been trained and converged in the (k-1)th stage, does not participate in gradient updates, and is only used to output the true values ​​of the physics predictions for the playback set, providing a constraint benchmark of the learned knowledge for the student network.

[0012] Preferably, in step S3, based on physical constraints and teacher-student consistency constraints, a total loss function that dynamically adjusts with the training phase is constructed. The loss function design is divided into two training phases, specifically including the following: S31. Dynamic loss function for serialized time-domain course learning: In the k-th training phase, the current time domain is... The total loss function consists of four parts: PDE residual loss, initial condition loss, boundary condition loss, and teacher-student consistency constraint loss. Its mathematical expression is: ; in, Loss due to initial conditions For boundary condition loss, For the PDE residual loss in the current extended time domain, For the loss due to teacher-student consistency constraints, the weighting coefficients are... , , All are configurable positive real numbers; S32, RAR enhancement and second-order fine-tuning loss function: After completing the full-time domain course, we enter the second-order fine-tuning stage. We remove the teacher-student consistency constraint loss and focus on improving accuracy in locally high-error regions. The loss function is simplified to: ; in, To enhance the residual loss of PDE, the mean square value of the residuals for all PDE collocation points, including high error points automatically extracted by the RAR mechanism, is calculated.

[0013] Preferably, in step S31, a dual dynamic termination mechanism of "stopping upon reaching the target + budget circuit breaker" is introduced during the training process. The dual dynamic termination criteria are as follows: Condition A, Physical Residue-Driven Adaptive Achievement Criterion: When the relative convergence slope, i.e., the absolute value of the change in the loss function over multiple consecutive iterations, is less than... When the current stage converges, it is considered that the current stage has converged. Condition B, computational budget-driven mandatory circuit breaker criterion: setting an iteration upper limit, i.e., the number of training epochs in a single stage reaches a preset upper limit of 20,000 epochs.

[0014] Preferably, in step S4, a two-stage hybrid optimization strategy of "first-stage serialization temporal exploration + second-stage RAR fine-tuning" is adopted to achieve efficient model optimization and accuracy breakthrough. The specific execution includes the following: S41. Serialization Time-Domain Exploration Based on Adam Optimizer: A discrete-time segmented advancement mechanism is adopted to complete N stages of progressive training, resulting in a basic model that converges across the entire time domain. S42. RAR Adaptive Enhancement and Fine-Tuning Based on L-BFGS Optimizer: On the basic model obtained in step S41, the training set is enhanced through the RAR mechanism, and the accuracy is improved by combining the second-order optimizer L-BFGS. By defining the closure function required by the L-BFGS optimizer, the entire process of gradient clearing, forward propagation, loss calculation, and backpropagation is completed within the closure function, and the total loss value is returned. The step method of the optimizer is called to execute the closure function and complete the convergence optimization. S43. After fine-tuning, save the final converged model weights to obtain the trained MoE-KAN-PINN model.

[0015] Preferably, in step S41, the training process includes the following: S411, Training Initialization: Load the full dataset to the computing device and initialize the index during the initialization phase. The preceding model path is empty; S412, Single-stage training loop: For the k-th training stage, perform the following operations: (1) Data slicing: Based on the current T_curr and T_prev, complete the automatic selection of training set and replay set; (2) Model initialization: Initialize the student network. If k≥2, load the model weights saved in stage k-1 as initialization parameters; at the same time, freeze the converged model in stage k-1 as the teacher network and pre-calculate the true values ​​of the replay set. (3) Optimizer and learning rate scheduler configuration: Adam first-order optimizer is adopted, the initial learning rate is a configurable positive real number, and the weight decay coefficient is a configurable positive real number; StepLR step learning rate scheduler is adopted to achieve rapid exploration in the early stage of training and steady convergence in the later stage. (4) Iterative training: The maximum number of iterations in a single stage is set to a configurable positive integer. Each iteration executes the standard process of gradient clearing, forward propagation, residual calculation, back propagation, and parameter update, and regularly prints monitoring indicators to track the training status. (5) Model saving: After a single stage of training is completed, the converged model weights of the current stage are saved as the initial weights and teacher network baselines for the next stage. S413, Phase Advancement: After completing the k-th phase of training, update the phase index. Repeat step S412 of the single-stage training process until... After completing the course learning in the full time domain [0, T], we obtain the basic model that converges in the full time domain.

[0016] Preferably, the enhancement of the training set through the RAR mechanism includes RAR residual adaptive resampling, specifically executed as follows: (1) Construction of high-density reconnaissance grid: Generate a high-density reconnaissance grid covering the entire time and space domain. The grid density is configured according to the computing resources and solution accuracy requirements. The typical value is tens of thousands to hundreds of thousands of sampling points to achieve full coverage without dead angles. (2) Residual calculation: Input the reconnaissance grid into the initialized model, and calculate the absolute value of the PDE residual at each sampling point through automatic differentiation. The larger the absolute value of the residual, the higher the model prediction error and the more drastic the change in the physical field. (3) High error point extraction: Extract the top K points with the largest absolute value of residuals in the reconnaissance grid as "difficult samples" for areas of drastic local change such as shock waves and interface abrupt changes, where K is a configurable positive integer; (4) Training set enhancement: The extracted K high error points are merged into the original PDE collocation set to form an enhanced PDE training set.

[0017] Preferably, in step S5, the collocation points obtained in step S1 are input into the trained MoE-KAN-PINN model to complete the end-to-end generation of the physics prediction solution. The specific process includes the following: S51. Input Normalization: Perform interval normalization on the input spatiotemporal coordinates (x, t); S52. Parallel Feature Extraction: Normalized spatiotemporal coordinates are input in parallel to a gated network and M heterogeneous KAN expert sub-networks. Each expert sub-network performs independent nonlinear feature mapping on the input coordinates, generating intermediate feature maps with consistent dimensions. The gated network synchronously calculates and outputs the normalized activation weights of each expert network. ,satisfy ; S53. Weighted Feature-Level Fusion: The intermediate feature maps output by each expert network are element-wise weighted with the corresponding activation weights output by the gating network. The weighted feature maps are then summed to obtain the fused comprehensive feature vector. The final prediction solution is obtained by mapping the comprehensive feature vector through the output layer. The mathematical expression is:

[0018] in, The first output of the gating network i The activation weights of an expert network For the first i Intermediate feature maps output by an expert network.

[0019] The beneficial effects of this invention are: 1) Ensuring strict adherence to physical causality and improving the stability of long-term solutions: This invention employs a dynamic time-domain expansion strategy based on curriculum learning, decomposing long-term time domains into progressive training phases. This allows the model to learn only the physical laws within the current time domain at each phase, expanding to the next phase only after sufficient convergence. This fundamentally eliminates the problem of missing intra-interval causality in existing Time-marching PINNs, strictly adhering to the temporal causal evolution logic of the physical system. Simultaneously, through a knowledge transfer mechanism for parameter initialization between phases, the model learns new content based on already mastered physical laws, significantly reducing the learning difficulty of long-term solutions and effectively preventing the model from getting trapped in local optima or training divergence due to complex spatiotemporal coupling, thus significantly improving the stability of long-term PDE solutions.

[0020] 2) Effectively suppressing catastrophic forgetting and ensuring the smoothness and consistency of the global solution: This invention innovatively introduces a teacher-student consistency constraint mechanism, freezing the converged model of the previous stage into the teacher network. By constraining the prediction differences between the student network and the teacher network in the learned time domain, it effectively retains the learned knowledge, fundamentally solving the catastrophic forgetting problem of existing sequential training methods. This mechanism allows the model to learn new time-domain physical laws without compromising the prediction accuracy of the learned time domain, ensuring the smoothness and consistency of the physical field solution in the entire spatiotemporal domain in the time direction. It avoids the problems of discontinuous solutions and numerical oscillations caused by multi-network splicing or segmented modeling in traditional methods, thus improving the rationality of the physics solution.

[0021] 3) Enhanced Complex Feature Capture Capability for Extremely High Solution Accuracy: This invention integrates MoE and KAN to construct a MoE-KAN heterogeneous expert network. Five heterogeneous expert sub-networks specifically capture different types and scales of features from the physical field, such as local shock waves, boundary features, and global smoothing features. The gated network enables adaptive expert selection based on the physical field characteristics, completely overcoming the feature representation limitations of traditional single-network structures. Simultaneously, a RAR data-driven shock wave capture mechanism is introduced to automatically locate and densify sampling points in high-error regions such as shock waves, providing high-density local physical constraints for model training. In tests using the one-dimensional Burgers equation, the MSE of the model's prediction results compared to the reference solution is as low as 1.020256 × 10⁻⁶. -8 The global relative L2 error is only 0.016%, achieving ultra-high precision solution for strongly nonlinear PDEs and completely eliminating the Gibbs oscillation phenomenon of traditional high-order interpolation methods.

[0022] 4) Overcoming the bottleneck of optimization accuracy and improving model optimization efficiency: This invention designs a two-stage hybrid optimization strategy of "global exploration + second-order fine-tuning," achieving complementary advantages between first-order and second-order optimizers: In the first stage, the Adam first-order optimizer, combined with a step-rate learning, quickly overcomes the non-convex optimization barrier and establishes a basic understanding of the global physical evolution law; in the second stage, the L-BFGS second-order optimizer utilizes the curvature information of the approximate Hessian matrix, combined with the RAR-enhanced training set, to perform global fine-tuning, completely overcoming the accuracy bottleneck of traditional first-order optimizers in areas of drastic local changes. This strategy reduces the model loss from 10 to 10 with extremely high computational efficiency. -4 The magnitude was directly reduced to 10 -8 This significantly reduces the blind spots and trial-and-error costs of hyperparameter tuning, thereby improving the overall optimization efficiency of the model.

[0023] 5) Enhancing the physical interpretability of the model and achieving region-adaptive modeling: The MoE-KAN heterogeneous expert architecture adopted in this invention effectively breaks the "black box" attribute of traditional deep learning models. The expert activation weights output by the gating network can intuitively reflect the contribution of each expert network in different spatiotemporal regions: in the shock region, the weight of Wav-KAN is significantly increased; in the boundary region, the weight of Chebyshev-KAN has the highest proportion; in the smooth region, Fourier-KAN or FCN becomes dominant. By extracting and analyzing the spatiotemporal heat map of the gating network, the physical space and time period responsible for each mathematical basis function network can be clearly observed, providing a powerful mathematical evidence tool for in-depth research on complex physical laws (such as mode transitions and energy cascades), and realizing the unity of data-driven region-adaptive modeling and physical interpretability.

[0024] 6) Excellent versatility and scalability, adaptable to solving various complex PDEs: The core technology modules of this invention (MoE-KAN, dynamic temporal advancement, teacher-student consistency constraints, RAR, and two-stage optimization) all adopt a modular design, independent of the specific form of PDEs, and can flexibly adapt to solving various PDEs containing shock waves, multiple scales, long time histories, and strong nonlinearities. Simultaneously, the modular structure of the model is naturally adapted to parallel computing environments. By increasing the number of expert networks and adjusting network structure parameters, it is possible to solve higher-dimensional and more complex PDEs, demonstrating excellent versatility and engineering scalability. This provides a novel and efficient deep learning solution for complex physics field simulation in scientific computing. Attached Figure Description

[0025] Figure 1 This is a schematic diagram of the process for solving heterogeneous expert partial differential equations based on temporal causality and teacher-student constraints in an embodiment of the present invention. Figure 2 A schematic diagram of the overall architecture of the hybrid expert physical information neural network MoE-KAN-PINN. Detailed Implementation

[0026] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0027] Kolmogorov-Arnold Networks (KANs) are built upon the classic Kolmogorov-Arnold representation theorem, revolutionizing the traditional network structure design logic of Multi-Level Processing (MLPs). They directly embed learnable single-variable nonlinear functions (such as B-spline curves, polynomials, and wavelet bases) into the connection edges of the neural network, realizing the nonlinear transformation from input to features as edge functions. Compared to traditional fully connected networks, KANs, with their unique structural design, exhibit significant advantages in the compactness of function approximation and parameter efficiency. On the one hand, the single-variable learnable nature of the edge functions allows the network to adaptively capture both local features and global trends of the physical field, effectively alleviating common problems in PDE solving in MLPs, such as gradient vanishing and spectral bias. On the other hand, KANs' mathematical expression is more closely aligned with the analytic characteristics of continuous functions, resulting in stronger numerical stability during the automatic calculation of higher-order derivatives, providing a new deep learning mathematical foundation for high-precision physical field simulation and PDE numerical solutions.

[0028] Mixture of Experts (MoE) models are a type of neural network architecture based on a divide-and-conquer conditional computation strategy. Their core concept is to achieve adaptive partitioning of the input space and dynamic allocation of computational tasks through the collaborative work of multiple complementary, specialized "expert subnetworks" and one or more "gated networks." In complex physical system modeling and PDE solving, physical fields (such as flow fields, stress fields, and neutron flux fields) often exhibit strong spatial non-uniformity and multi-scale characteristics—some regions show smooth, low-frequency changes, while shock waves, interfaces, and gradient abrupt changes exhibit drastic high-frequency jumps. A single network structure cannot simultaneously address the feature representation needs of different regions. The MoE architecture is well-suited to such multi-scale, non-uniform modeling scenarios: each expert subnetwork can be specifically designed to capture a certain type and scale of physical features, while the gated network, based on the spatiotemporal coordinates and feature attributes of the input data, uses learned weight distributions to allocate different inputs to the most suitable expert subnetwork for computation. Finally, the outputs of each expert are weighted and fused to obtain the final prediction. This architectural design allows the model to significantly improve its depth of expression and adaptability to multi-scale physical phenomena while maintaining overall computational efficiency. At the same time, MoE's modular structure is naturally adapted to parallel computing environments, which can increase model capacity without significantly increasing inference costs. Its data-driven region partitioning capability does not require manual pre-setting of boundaries or labels, providing a more flexible and efficient solution for modeling complex physical systems.

[0029] Curriculum Learning (CL) is a machine learning training strategy that mimics the progressive learning patterns of humans and animals, moving from simple to complex and from basic to advanced. Its core is to organize and present training samples in a meaningful order, first allowing the model to learn simple concepts and samples, and then gradually introducing complex concepts and samples, guiding the model to complete the learning process from basic to advanced. This strategy is also known as shaping in the fields of cognitive science and animal training, and is an important means of solving non-convex optimization problems and improving the model's generalization ability and convergence speed.

[0030] In deep learning-driven PDE solving, especially long-term time-varying PDE solving tasks, traditional training strategies employ global spatiotemporal random sampling and point allocation for joint optimization, which suffers from significant learning efficiency and convergence problems: neural networks struggle to directly capture the evolutionary laws of physical fields across vast time spans, and are prone to getting trapped in local optima in the early stages of training due to complex spatiotemporal coupling relationships. Furthermore, they are susceptible to violating the temporal causality of physical systems, leading to model divergence and predictions that contradict actual physical laws. The core adaptability of CL (Constant Learning) lies precisely in its dual solution to the temporal evolution characteristics of PDEs and the non-convex optimization problem in machine learning.

[0031] From the perspective of optimization principles, course learning can be viewed as a special type of extension method—a global optimization strategy for non-convex objective functions. Its core is to first optimize a smooth objective function, then gradually reduce the smoothness until it approximates the original objective function, allowing the model parameters to gradually enter a more favorable region of attraction of the original objective function. Specifically, in course learning, it defines a set of parameters... Changing training distribution The initial distribution focuses on simple samples, then the entropy of the distribution is gradually increased to include more complex samples, eventually converging to the target training distribution. The weighting of samples during this process. Follow Monotonically increasing ensures that the model always learns new content based on the knowledge it has already mastered, greatly reducing the learning difficulty.

[0032] From the perspective of practical applications of solving partial differential equations, the physical field distribution of time-varying partial differential equations evolves gradually over time from initial conditions, and the physical state at subsequent moments is highly dependent on the accumulation of previous states. This characteristic aligns perfectly with the progressive learning logic of the course. Based on this physical law, the course can decompose the long-term spatiotemporal solution domain into progressive learning stages along the time dimension. The teaching syllabus guides the model by first establishing stable physical constraints within a small time range near the initial moment, accurately fitting the initial conditions and the physical evolution laws within a short timeframe; then, as the training rounds increase, the solution domain is gradually pushed towards later time periods until it covers the entire spatiotemporal domain.

[0033] This training strategy strictly adheres to the causal evolution logic of physical systems while addressing the core pain points of traditional global optimization: on the one hand, it avoids the model getting trapped in local optima due to premature exposure to complex spatiotemporal coupling samples, effectively enhancing the convergence stability of solving complex nonlinear dynamic systems; on the other hand, it allows the model to gradually accumulate knowledge of physical laws during training, improving the accuracy and generalization ability of numerical predictions. Furthermore, course learning can significantly accelerate the convergence speed of models to the global optimum in convex optimization problems, making it one of the core optimization methods for solving long-term, multi-scale partial differential equations using deep learning. Its ideas have also been extended to other machine learning fields such as active learning, transfer learning, and lifelong learning.

[0034] This invention addresses the core technical problems faced by existing Physical Information Neural Networks (PINNs) in solving long-term, strongly nonlinear partial differential equations (PDEs), such as causal failure, catastrophic forgetting, and accuracy bottlenecks. It introduces a curriculum learning (CL) framework based on dynamic time-domain extension, fundamentally reshaping the training paradigm of PINNs and achieving strict adherence to the causal evolution sequence of physical systems.

[0035] To address the technical shortcomings of existing time-progression methods that still employ global joint optimization within sub-intervals, leading to a lack of causality within the intervals, this invention employs a progressive time-domain expansion strategy, starting with shorter intervals and gradually increasing the time domain. This ensures that the model learns only the physical laws within the current time domain at each stage, and only expands to the next stage after sufficient convergence, fundamentally eliminating the problem of missing causality within intervals.

[0036] To address the catastrophic forgetting problem prevalent in existing sequential training methods—where the model severely compromises the prediction accuracy of previously learned time domains while learning new ones—this invention innovatively introduces a teacher-student consistency constraint mechanism to effectively retain learned knowledge and ensure the consistency and smoothness of the global solution in the time direction.

[0037] To address the shortcomings of existing PINN methods in capturing features and achieving low accuracy in regions with drastic local changes, such as shock waves, this invention integrates Hybrid Expert (MoE) and Kolmogorov-Arnold Network (KAN) to construct a heterogeneous expert pool. It combines a Residual Adaptive Resampling (RAR) mechanism to accurately capture local high-frequency features. At the same time, it designs a two-stage hybrid optimization strategy of "global exploration + second-order fine-tuning" to break through the accuracy bottleneck of traditional first-order optimizers and improve the overall solution accuracy and optimization efficiency of the model.

[0038] This invention provides a technical solution: a method for solving heterogeneous expert partial differential equations based on temporal causality and teacher-student constraints. The technical solution process is as follows: Figure 1 As shown. This invention integrates three core technologies: Hybrid Expert (MoE), Kolmogorov-Arnold Network (KAN), and Curriculum Learning (CL). It introduces teacher-student consistency constraints, Residual Adaptive Resampling (RAR), and a two-stage hybrid optimization strategy to construct a Hybrid Expert Physical Information Neural Network (MoE-KAN-PINN), achieving high-precision and high-stability solutions for long-term, strongly nonlinear PDEs. The Hybrid Expert Physical Information Neural Network (MoE-KAN-PINN) adopts an end-to-end physical guidance structure, suitable for solving complex PDEs containing shock waves, multi-scale features, and long-term evolution. The overall framework is as follows: Figure 2 As shown. The model input is spatiotemporal coordinates. The output is the predicted physical field value at the corresponding location. This approach integrates a hybrid expert mechanism, curriculum learning strategies, teacher-student consistency constraints, and a two-stage hybrid optimization technique. Its core is divided into three main modules: a feature extraction layer, a feature fusion layer, and an output layer. It also incorporates a teacher-student network interaction mechanism to ensure knowledge retention. The overall solution method includes the following steps: Step 1: Spatiotemporal Domain Construction and Full Dataset Preloading To address the requirements of solving the target PDE, a generalized spatiotemporal domain definition and the construction and preloading of the full training dataset were completed, providing a standardized data foundation for subsequent phased training. The specific execution is as follows: Step 1.1: Definition of the general spatiotemporal domain. Let the spatial solution domain of the target PDE be... (d is the spatial dimension, typically 1, 2, or 3), the total time domain is Where T is the maximum solution time as defined by the physical problem; Step 1.2: Generation and Preloading of the Full Dataset. Latin Hypercube Sampling (LHS) is used to generate collocations within the PDE. The internal PDE collocation set is then generated. In the entire space-time domain Internal LHS sampling. Core optimization: Immediately after sampling, the samples are sorted globally in ascending order according to the time coordinate t to meet the slicing requirements of causal training.

[0039] Step 1.2.1: Initial condition point set :fixed , Perform 1D LHS sampling on the spatial domain x.

[0040] Step 1.2.2: Boundary condition point set At the spatial boundary and At each point, LHS sampling is performed on the time axis t to enforce the boundary values ​​of the physical system.

[0041] Step 1.3: Dynamic Adaptive Data Slicing Rules. Based on the time-domain expansion strategy learned in the course, a generalized time segmentation plan is predefined to divide the total time domain... Divided into Each progressive training phase, with a single phase time step. Equal step size can be used ( ) settings; for the first Each training phase (corresponding to the current maximum training time) The maximum time already learned in the previous steps Automatically complete two-level data slicing: Step 1.3.1: Training set selection: Extraction Boundary condition points and PDE collocation points within the range are used as training data for the current stage; Step 1.3.2: Replay Set Filtering: Extraction PDE collocations within the range, used as the replay dataset for teacher-student consistency constraints, are only used within... It will take effect at that time.

[0042] Step 2: Construction of Hybrid Expert Physical Information Neural Network (MoE-KAN-PINN) and Teacher-Student Dual Network Architecture A generalized heterogeneous expert physical information neural network integrating MoE and KAN is constructed. A teacher-student dual-network architecture is designed to achieve knowledge retention and progressive learning, as detailed below: Step 2.1: Construct a heterogeneous expert library consisting of four parallel heterogeneous KAN expert sub-networks and one FCN expert. Each expert network uses normalized spatiotemporal coordinates (x,t) as input and maps them to intermediate feature maps of consistent dimension through distinctive mathematical bases or activation functions, thereby achieving adaptive representation of multi-scale, non-uniform features of the physical field. Wav-KAN: Based on the Mexican hat wavelet function with learnable translation and scaling parameters, it is specifically designed to capture the characteristics of drastic gradient changes such as shock waves and interface abrupt changes, thanks to its excellent local time-frequency localization capabilities. RBF-KAN: It uses Gaussian kernel radial basis functions as the basis and fits spatial clustering features through local receptive fields, making it suitable for capturing local non-uniform distributions in physical fields; Chebyshev-KAN: Based on Chebyshev polynomial expansion, it utilizes the property of minimizing the maximum error in the domain to improve the model's approximation accuracy of the physical field in the boundary region. Fourier-KAN: Utilizing sine and cosine Fourier basis functions, it leverages the advantage of global orthogonal basis functions to capture the global evolution trend and smoothing components of the physical field; FCN: As a traditional MLP expert, it uses the Tanh activation function to provide a basic nonlinear smooth approximation, complementing the KAN expert.

[0043] Each expert network obtains differentiated feature preferences through random initialization, and achieves natural functional differentiation under the drive of physical loss, jointly constructing a comprehensive feature space representation of complex physical fields.

[0044] Step 2.2: The gated network adopts a multi-layer fully connected structure with two hidden layers, each layer containing 64 neurons. The activation function is Tanh, and the input is normalized spatiotemporal coordinates. The output layer is a Softmax layer, which ultimately yields a probability vector of dimension 5. This vector represents the activation weights of the five expert networks and satisfies the normalization condition that the sum of the weights is 1.

[0045] The core function of the gating network is to dynamically allocate expert weights based on the local physical properties of the spatiotemporal coordinates: for shock regions, the weights of Wav-KAN are automatically amplified; for boundary regions, the weights of Chebyshev-KAN are automatically amplified; and for smooth regions, the weights of Fourier-KAN or FCN are automatically amplified, thereby achieving adaptive expert selection based on physical field characteristics, allowing the model to call the most suitable feature extraction network in different regions. Step 2.3: Definition of Teacher-Student Dual Network Architecture Step 2.3.1: Student Network: The trainable MoE-KAN-PINN model for the current training phase. During the k-th training phase, the model is loaded with the model from the previous phase. The model weights that converge in each stage are used as initialization parameters (random initialization in stage 1) and participate in gradient updates and parameter optimization throughout the process; Step 2.3.2: Teacher Network: Only in When activated, the first Once the training phase is complete and converged, the model parameters are completely frozen and do not participate in gradient updates. They are only used to output the true values ​​of the physics predictions for the playback set, providing a constraint benchmark of the learned knowledge for the student network.

[0046] Step 3: Constructing the dynamic staged total loss function Based on physical constraints and teacher-student consistency constraints, a total loss function is constructed that dynamically adjusts with the training phase to achieve a balance between "learning new knowledge" and "retaining old knowledge," while strictly adhering to the causal evolution law of the physical system. Specifically, the loss function design is divided into two training phases: Step 3.1 First Stage: Dynamic Loss Function for Serialized Temporal Course Learning The kth training phase (current time domain is...) The total loss function consists of four parts: PDE residual loss, initial condition loss, boundary condition loss, and teacher-student consistency constraint loss. Its mathematical expression is: ; The parameters and their physical meanings are clearly defined as follows: Initial condition loss is calculated using Mean Squared Error (MSE) to determine the deviation between the model's predicted values ​​and the true initial values ​​at the initial time step, with weighting coefficients. It is a configurable positive real number; Boundary condition loss is calculated using mean squared error to determine the deviation between the model's predicted boundary values ​​and the actual boundary values, with weighting coefficients. It is a configurable positive real number, ensuring strict satisfaction of boundary conditions; : The PDE residual loss in the current extended time domain is calculated. The mean square value of the residuals of the control equations for all PDE collocation points within the range ensures that the model strictly satisfies the physical control equations in the current time domain; The loss due to the teacher-student consistency constraint only occurs when... Effective immediately, the mean square error is used to calculate the student network and teacher network in the playback set ( The prediction bias on ) is expressed as:

[0047] in This represents the total number of sampling points in the playback set. For student network prediction values, For frozen teacher network predictions; weighting coefficients It is a configurable positive real number used to balance physical constraints and knowledge retention, and to suppress catastrophic forgetting; Step 3.2: Introducing Dual Dynamic Termination Criteria and Training Logic During Training. To further balance training efficiency and physical accuracy in long-term solutions, this method introduces a dual dynamic termination mechanism of "stopping upon reaching the target + budget circuit breaker" at each time-domain expansion stage. This mechanism aims to ensure that the model can advance in a timely manner after mastering the physical laws of the current time domain, while avoiding getting stuck in invalid iterations in complex nonlinear regions. The specific criteria are as follows: (1) Condition A: Adaptive achievement criterion driven by physical residuals This criterion is used to evaluate the model's understanding of the physical evolution laws at the current stage. The model is considered to have converged at the current stage when any of the following sub-conditions are met: Relative convergence slope: The absolute value of the change in the loss function over multiple consecutive iterations is less than... .

[0048] This criterion references the stability principle in causal sweeping strategies, ensuring that the model has established solid physical constraints within the learned time domain before moving to the next time domain extension.

[0049] (2) Condition B: Calculate the budget-driven mandatory circuit breaker criterion This criterion serves as a "fallback" mechanism for the algorithm, used to handle extreme nonlinear regions such as those containing shock waves. The specific settings are as follows: Iteration limit setting: The number of training epochs in a single stage reaches a preset limit (e.g., 20,000 epochs).

[0050] When dealing with highly nonlinear problems such as the Burgers equation, the drastic gradient changes at the shock front can lead to "gradient deadlock" in the first-order optimizer (Adam). By employing a forced circuit breaker mechanism, the algorithm can interrupt the current inefficient iteration and force a move to the next time step. This strategy activates the gradient flow by introducing new time-step sampling points, helping the model escape local optima and ensuring the causal progression of the overall training process.

[0051] Through the synergistic effect of the above dual criteria, this method achieves on-demand allocation of computational resources along the time axis: rapid advancement in smooth evolution regions and steady exploration in complex shock wave regions, fundamentally balancing the real-time performance and numerical accuracy of the solution.

[0052] Step 3.3 Second Stage: RAR Enhancement and Second-Order Fine-Tuning of the Loss Function After completing the full-time domain course, we enter the second-order fine-tuning stage. At this point, the model has grasped the global physical evolution laws. We remove the teacher-student consistency constraint loss and focus on improving accuracy in local high-error regions. The total loss function simplifies to: ; in To enhance the residual loss of PDE, the mean square value of the residuals for all PDE collocation points, including high-error points automatically extracted by the RAR mechanism, is calculated. and Maintain the same values ​​as in the first stage to ensure the rigid constraints of the physical boundary conditions.

[0053] Step 4: Execution of the two-stage hybrid optimization strategy throughout the entire process A two-stage hybrid optimization strategy of "first-stage serialization temporal exploration + second-stage RAR fine-tuning" is adopted, combining the global exploration capability of the first-order optimizer and the local fine-tuning capability of the second-order optimizer to achieve efficient model optimization and breakthrough in accuracy. The specific implementation is as follows: Step 4.1 First Stage: Serialization Time Domain Course Learning and Teacher-Student Constraint Training Based on Adam Optimizer This stage employs a discrete-time segmented progression mechanism to complete N stages of progressive training. The specific execution flow is as follows: Step 4.1.1: Training Initialization. Load the full dataset to the computing device and initialize the index during the initialization phase. The preceding model path is empty; Step 4.1.2: Single-stage training loop. For the k-th training stage, perform the following operations: (1) Data slicing: based on the current stage and Automatically select the training set and playback set; (2) Model initialization: Initialize the student network, if Loading the first The model weights saved in the first stage are used as initialization parameters; at the same time, the model weights in the second stage are frozen. The convergence model of the stage serves as the teacher network, pre-computing the predicted ground truth values ​​for the playback set; (3) Optimizer and learning rate scheduler configuration: Adam first-order optimizer is adopted, the initial learning rate is a configurable positive real number, and the weight decay coefficient is a configurable positive real number; StepLR step learning rate scheduler is adopted to achieve rapid exploration in the early stage of training and steady convergence in the later stage. (4) Iterative training: The maximum number of iterations in a single stage is set to a configurable positive integer. Each iteration executes the standard process of gradient clearing, forward propagation, residual calculation, back propagation, and parameter update, and regularly prints monitoring indicators to track the training status. (5) Model saving: After a single stage of training is completed, the converged model weights of the current stage are saved as the initial weights and teacher network baselines for the next stage. Phase progression: After completing the k-th phase of training, update the phase index. Repeat the single-stage training process in 4.1.2 until... Complete the course learning in the full time domain [0, T] and obtain the basic model that converges in the full time domain; Step 4.2 Second Stage: RAR Adaptive Enhancement Fine-Tuning Based on L-BFGS Optimizer This stage builds upon the converged model from the first stage by enhancing the training set through the RAR mechanism and combining it with a second-order optimizer to achieve a breakthrough in accuracy. The specific execution flow is as follows: Step 4.2.1: Model and Data Initialization. Load the converged model weights from the final output of the first stage as initialization parameters for the fine-tuning stage; load the full initial conditions, boundary conditions, and the original PDE collocation dataset; Step 4.2.2: RAR residual adaptive resampling is performed.

[0054] (1) Construction of high-density reconnaissance grid: Generate a high-density reconnaissance grid covering the entire time and space domain. The grid density is configured according to the computing resources and solution accuracy requirements. The typical value is tens of thousands to hundreds of thousands of sampling points to achieve full coverage without dead angles. (2) Residual calculation: Input the reconnaissance grid into the initialized model, and calculate the absolute value of the PDE residual at each sampling point through automatic differentiation. The larger the absolute value of the residual, the higher the model prediction error and the more drastic the change in the physical field. (3) High error point extraction: Extract the top K points with the largest absolute residual value in the reconnaissance grid (K is a configurable positive integer, typically hundreds to thousands) as “difficult samples” for areas of drastic local change such as shock waves and interface abrupt changes; (4) Training set enhancement: The K high error points extracted are merged into the original PDE collocation set to form an enhanced PDE training set; Step 4.2.3: Second-order optimizer configuration. The L-BFGS second-order optimizer is adopted. Its learning rate, maximum number of iterations, historical gradient storage, and line search strategy can all be configured according to the problem characteristics. High-precision gradient updates are achieved by utilizing the curvature information of the approximate Hessian matrix. Step 4.2.4: Closure-based fine-tuning training. Define the closure function required by the L-BFGS optimizer. The closure function completes the entire process of gradient zeroing, forward propagation, loss calculation, and backpropagation, and returns the total loss value. Call the optimizer's step method to execute the closure function and complete the convergence optimization. Step 4.2.5: Save the final model. After fine-tuning, save the final converged model weights to obtain the trained MoE-KAN-PINN model.

[0055] Step 5: Generation of Physics Field Prediction Solutions The collocation points obtained in step 1 are input into the trained MoE-KAN-PINN model to complete the end-to-end generation of the physics prediction solution. The specific process is as follows: Step 5.1: Input Normalization. Perform interval normalization on the input spatiotemporal coordinates (x, t) (e.g., normalize to [-1, 1] or [0, 1]) to eliminate the impact of dimensional differences on model training; Step 5.2: Parallel Feature Extraction. The normalized spatiotemporal coordinates are input in parallel to a gating network and M heterogeneous expert networks. Each expert network performs independent nonlinear feature mapping on the input coordinates, generating intermediate feature maps with consistent dimensions. The gated network synchronously calculates and outputs the normalized activation weights of each expert network. ,satisfy ; Step 5.3: Employ a weighted feature-level fusion method. The intermediate feature maps output by each expert network are element-wise weighted with the corresponding activation weights output by the gating network. The weighted feature maps are then summed to obtain the fused comprehensive feature vector. The final prediction solution is obtained by mapping the comprehensive feature vector through the output layer. The mathematical expression is:

[0056] in, Let be the activation weights of the i-th expert network output by the gated network. This is the intermediate feature map output by the i-th expert network. The feature-level fusion mechanism preserves the diversity of feature representations while retaining the underlying shared features of the physical field, which helps improve the model's generalization ability and physical consistency.

[0057] By executing the above five steps in a coordinated manner, the entire PDE solution process is completed. The method of this invention achieves high-precision solution for long-term, strongly nonlinear PDEs, while effectively solving technical problems such as causal failure, catastrophic forgetting, and insufficient local feature capture capability of existing methods.

[0058] Experimental Verification: An Example of Solving the One-Dimensional Burgers Equation The physical problem configuration for the one-dimensional Burgers equation is shown in Table 1.

[0059] Table 1 Physical Problem Configuration for One-Dimensional Burgers Equations

[0060] The solution is performed according to the steps in the embodiments of the present invention, as follows: 1. Data sampling settings Data sampling was performed using the Latin hypercube sampling method, with the following number of sampling points: Initial condition point allocation: 500 points, evenly distributed in The initial moment; Boundary condition points: 500, evenly distributed in and The boundary position; Internal spatiotemporal points: 10,000, evenly distributed across the spatiotemporal domain. , Within; all sampling points are sorted in ascending order by time coordinate and then... Normalization preprocessing.

[0061] 2. Network Structure and Training Parameter Settings The network structure and training parameters of MoE-KAN-PINN are strictly set according to the code configuration. The core parameters are as follows: MoE-KAN Expert Architecture: (1) Heterogeneous KAN experts: including four types of experts: Wav-KAN, RBF-KAN, Chebyshev-KAN and Fourier-KAN, with the hidden layer dimension set to 64; (2) FCN expert: A 3-layer hidden layer structure is used, with the number of neurons being [128, 128, 64], and the activation function is Tanh; (3) Gated network: It adopts a 2-layer fully connected hidden layer structure, each layer contains 64 neurons, and the activation function is Tanh.

[0062] The first training phase: The total time domain is divided into 10 training phases, with Δt = 0.1s and each phase having 20,000 epochs; the Adam optimizer has an initial learning rate of 1 × 10⁻⁶. -3 It decays by 50% every 10,000 rounds; RAR mechanism: The high-density grid is 500×200, and the 2000 high error points with the largest residuals are extracted and expanded into the training set; Second stage training: L-BFGS optimizer learning rate lr=1.0, historical gradient storage 100, maximum number of iterations 5000 rounds.

[0063] 3. Implementation Results After two stages of hybrid optimization training, the MoE-KAN-PINN model constructed by the method of this invention achieves high-precision solution of the one-dimensional Burgers equation, with the core performance being: Accurate capture of shock wave features: The model perfectly captures the extremely steep shock wave cross section in the later stage of time evolution (t≈1.0s) without any numerical oscillations. Compared with the traditional PINN method, the fitting accuracy of the shock wave region is significantly improved. The solution is consistent across the entire domain: the teacher-student consistency constraint mechanism effectively suppresses catastrophic forgetting, and the physical field solution in the entire time and space domain is smooth and continuous in the time direction, without any discontinuity or abrupt change, and strictly follows the evolution law of the physical system. Extremely high solution accuracy: Evaluated on a standard test set of 25,600 grid points, the mean square error (MSE) between the model's predictions and the true reference solution is as low as 1.020256 × 10⁻⁶. -8 The mean absolute error (MAE) is only 7.341510×10 -5 The global relative L2 error was controlled at 1.644173 × 10⁻⁶. -4 (0.016%), as shown in Table 2.

[0064] Table 2 Comparison of accuracy between the present invention and existing technologies

[0065] This invention presents a heterogeneous expert partial differential equation (PDE) solution method based on temporal causality and teacher-student constraints. It achieves strict adherence to causality, suppression of catastrophic forgetting, accurate capture of complex features, and dual improvement in solution accuracy and efficiency. It addresses the technical pain points of existing methods in solving long-term, strongly nonlinear PDEs and provides a more efficient and accurate deep learning solution for simulating complex physical fields in the field of scientific computing.

[0066] It should be noted that, in this document, the terms "comprising," "including," or any other variations thereof are intended to cover non-exclusive inclusion, such that a process, method, article, or apparatus that comprises a list of elements includes not only those elements but also other elements not expressly listed, or elements inherent to such a process, method, article, or apparatus. Unless otherwise specified, an element defined by the phrase "comprising one..." does not exclude the presence of other identical elements in the process, method, article, or apparatus that includes said element.

[0067] The terminology used in the embodiments of this invention is for the purpose of describing particular embodiments only and is not intended to limit the invention. The singular forms “a,” “the,” and “the” as used in the embodiments of this invention and the appended claims are also intended to include the plural forms unless the context clearly indicates otherwise.

[0068] It should be understood that the term "and / or" used in this article is merely a description of the relationship between related objects, indicating that three relationships can exist. For example, A and / or B can represent: A existing alone, A and B existing simultaneously, and B existing alone. Additionally, the character " / " in this article generally indicates that the preceding and following related objects have an "or" relationship.

[0069] Depending on the context, the word "if" as used here can be interpreted as "when," "when," "in response to determination," or "in response to detection." Similarly, depending on the context, the phrase "if determination" or "if detection (of the stated condition or event)" can be interpreted as "when determination," "in response to determination," "when detection (of the stated condition or event)," or "in response to detection (of the stated condition or event)."

[0070] Although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art can still modify the technical solutions described in the foregoing embodiments or make equivalent substitutions for some of the technical features. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

Claims

1. A method for solving heterogeneous expert partial differential equations based on temporal causality and teacher-student constraints, characterized in that, Includes the following steps: S1. Spatiotemporal domain construction and full dataset preloading; S2. Construct a hybrid expert physical information neural network MoE-KAN-PINN and a dual network architecture for teachers and students; S3. Construction of dynamic staged total loss function; S4. Execute a two-stage hybrid optimization strategy to train the MoE-KAN-PINN model; S5. Generation of physical field prediction solutions.

2. The method for solving heterogeneous expert partial differential equations based on temporal causality and teacher-student constraints according to claim 1, characterized in that: Step S1 specifically includes the following: S11. Definition of General Spatiotemporal Domain: Let the spatial solution domain of the target PDE be... d is the spatial dimension, and the total time domain is Where T is the maximum solution time as defined by the physical problem; S12. Generation and Preloading of the Full Dataset: Latin Hypercubic Sampling (LHS) is used to generate internal collocations of the PDE. The set of internal PDE collocations is as follows: In all time and space domains Internal LHS sampling is performed, and the samples are immediately sorted globally in ascending order according to the time coordinate t. S13. Dynamic Adaptive Data Slicing Rules: Based on the time-domain expansion strategy of course learning, a predefined generalized time segmentation plan is used to expand the total time domain. Divided into Each progressive training phase, with a single phase time step. Equal step size settings can be used; for the first Each training phase automatically completes two levels of data slicing; specifically, it includes training set selection and replay set selection. Training set selection: extraction Boundary condition points and PDE collocation points within the range are used as training data for the current stage; Indicates the current maximum training time; Playback set filtering: Extraction PDE collocations within the range, used as the replay dataset for teacher-student consistency constraints, are only used within... Effective immediately; This indicates the maximum time that the preceding sequence has been learned.

3. The method for solving heterogeneous expert partial differential equations based on temporal causality and teacher-student constraints according to claim 1, characterized in that: Step S2 specifically includes the following: S21. Construct a heterogeneous expert library consisting of four parallel heterogeneous KAN expert sub-networks and one FCN expert. Each expert network uses normalized spatiotemporal coordinates (x,t) as input and maps them to intermediate feature maps with consistent dimensions through unique mathematical bases or activation functions to achieve adaptive expression of multi-scale and non-uniform features of physical fields. The four parallel heterogeneous KAN expert sub-networks include: Wav-KAN, RBF-KAN, Chebyshev-KAN, and Fourier-KAN. S22. The gated network adopts a multi-layer fully connected structure with two hidden layers. Each layer has 64 neurons. The activation function is Tanh. The input is the normalized spatiotemporal coordinates (x,t). The output layer is a Softmax layer. Finally, a probability vector with a dimension of 5 is obtained. This vector is the activation weights of the five expert networks and satisfies the normalization condition that the sum of the weights is 1. S23. Definition of Teacher-Student Dual Network Architecture: Student Network: The trainable MoE-KAN-PINN model in the current training phase. During the k-th training phase, the model weights converged in the (k-1)-th phase are loaded as initialization parameters and participate in gradient updates and parameter optimization throughout the process. Teacher Network: Enabled only when k≥2, it completely freezes the model parameters that have been trained and converged in the (k-1)th stage, does not participate in gradient updates, and is only used to output the true values ​​of the physics predictions for the playback set, providing a constraint benchmark of the learned knowledge for the student network.

4. The method for solving heterogeneous expert partial differential equations based on temporal causality and teacher-student constraints according to claim 1, characterized in that: In step S3, based on physical constraints and teacher-student consistency constraints, a total loss function that dynamically adjusts with the training phase is constructed. The loss function design is divided into two training phases, specifically including the following: S31. Dynamic loss function for serialized time-domain course learning: In the k-th training phase, the current time domain is... The total loss function consists of four parts: PDE residual loss, initial condition loss, boundary condition loss, and teacher-student consistency constraint loss. Its mathematical expression is: ； in, Loss due to initial conditions For boundary condition loss, For the PDE residual loss in the current extended time domain, For the loss due to teacher-student consistency constraints, the weighting coefficients are... , , All are configurable positive real numbers; S32, RAR enhancement and second-order fine-tuning loss function: After completing the full-time domain course, we enter the second-order fine-tuning stage. We remove the teacher-student consistency constraint loss and focus on improving accuracy in locally high-error regions. The loss function is simplified to: ； in, To enhance the residual loss of PDE, the mean square value of the residuals for all PDE collocation points, including high error points automatically extracted by the RAR mechanism, is calculated.

5. The method for solving heterogeneous expert partial differential equations based on temporal causality and teacher-student constraints according to claim 4, characterized in that: In step S31, a dual dynamic termination mechanism of "stopping upon reaching the target + budget circuit breaker" is introduced during the training process. The dual dynamic termination criteria are as follows: Condition A, Physical Residue-Driven Adaptive Achievement Criterion: When the relative convergence slope, i.e., the absolute value of the change in the loss function over multiple consecutive iterations, is less than... When the current stage converges, it is considered that the current stage has converged. Condition B, computational budget-driven mandatory circuit breaker criterion: setting an iteration upper limit, i.e., the number of training epochs in a single stage reaches a preset upper limit of 20,000 epochs.

6. The method for solving heterogeneous expert partial differential equations based on temporal causality and teacher-student constraints according to claim 1, characterized in that: In step S4, a two-stage hybrid optimization strategy of "first-stage serialization temporal exploration + second-stage RAR fine-tuning" is adopted to achieve efficient model optimization and accuracy breakthrough. The specific execution includes the following: S41. Serialization Time-Domain Exploration Based on Adam Optimizer: A discrete-time segmented advancement mechanism is adopted to complete N stages of progressive training, resulting in a basic model that converges across the entire time domain. S42. RAR Adaptive Enhancement and Fine-Tuning Based on L-BFGS Optimizer: On the basic model obtained in step S41, the training set is enhanced through the RAR mechanism, and the accuracy is improved by combining the second-order optimizer L-BFGS. By defining the closure function required by the L-BFGS optimizer, the entire process of gradient clearing, forward propagation, loss calculation, and backpropagation is completed within the closure function, and the total loss value is returned. The step method of the optimizer is called to execute the closure function and complete the convergence optimization. S43. After fine-tuning, save the final converged model weights to obtain the trained MoE-KAN-PINN model.

7. The method for solving heterogeneous expert partial differential equations based on temporal causality and teacher-student constraints according to claim 6, characterized in that: In step S41, the training process includes the following: S411, Training Initialization: Load the full dataset to the computing device and initialize the index during the initialization phase. The preceding model path is empty; S412, Single-stage training loop: For the k-th training stage, perform the following operations: (1) Data slicing: Based on the current T_curr and T_prev, complete the automatic selection of training set and replay set; (2) Model initialization: Initialize the student network. If k≥2, load the model weights saved in stage k-1 as initialization parameters; at the same time, freeze the converged model in stage k-1 as the teacher network and pre-calculate the true values ​​of the replay set. (3) Optimizer and learning rate scheduler configuration: The Adam first-order optimizer is adopted, the initial learning rate is a configurable positive real number, and the weight decay coefficient is a configurable positive real number; the StepLR step learning rate scheduler is adopted to achieve rapid exploration in the early stage of training and steady convergence in the later stage. (4) Iterative training: The maximum number of iterations in a single stage is set to a configurable positive integer. Each iteration executes the standard process of gradient clearing, forward propagation, residual calculation, back propagation, and parameter update, and regularly prints monitoring indicators to track the training status. (5) Model saving: After a single stage of training is completed, the converged model weights of the current stage are saved as the initial weights and teacher network baselines for the next stage. S413, Phase Advancement: After completing the k-th phase of training, update the phase index. Repeat step S412 of the single-stage training process until... After completing the course learning in the full time domain [0, T], we obtain the basic model that converges in the full time domain.

8. The method for solving heterogeneous expert partial differential equations based on temporal causality and teacher-student constraints according to claim 6, characterized in that: The enhancement of the training set through the RAR mechanism includes RAR residual adaptive resampling, which is specifically executed as follows: (1) High-density reconnaissance grid construction: Generate a high-density reconnaissance grid covering the entire time and space domain. The grid density is configured according to the computing resources and solution accuracy requirements. The typical value is tens of thousands to hundreds of thousands of sampling points to achieve full coverage without dead angles. (2) Residual calculation: Input the reconnaissance grid into the initialized model, and calculate the absolute value of the PDE residual at each sampling point through automatic differentiation. The larger the absolute value of the residual, the higher the model prediction error and the more drastic the change in the physical field. (3) High error point extraction: Extract the top K points with the largest absolute value of residuals in the reconnaissance grid as "difficult samples" for areas of drastic local change such as shock waves and interface abrupt changes, where K is a configurable positive integer; (4) Training set enhancement: The extracted K high error points are merged into the original PDE collocation set to form an enhanced PDE training set.

9. The method for solving heterogeneous expert partial differential equations based on temporal causality and teacher-student constraints according to claim 1, characterized in that: In step S5, the collocation points obtained in step S1 are input into the trained MoE-KAN-PINN model to complete the end-to-end generation of the physics prediction solution. The specific process includes the following: S51. Input Normalization: Perform interval normalization on the input spatiotemporal coordinates (x, t); S52. Parallel Feature Extraction: Normalized spatiotemporal coordinates are input in parallel to a gated network and M heterogeneous KAN expert sub-networks. Each expert sub-network performs independent nonlinear feature mapping on the input coordinates, generating intermediate feature maps with consistent dimensions. The gated network synchronously calculates and outputs the normalized activation weights of each expert network. ,satisfy ; S53. Weighted Feature-Level Fusion: The intermediate feature maps output by each expert network are element-wise weighted with the corresponding activation weights output by the gating network. The weighted feature maps are then summed to obtain the fused comprehensive feature vector. The final prediction solution is obtained by mapping the comprehensive feature vector through the output layer. The mathematical expression is: ； in, The first output of the gating network i The activation weights of an expert network For the first i Intermediate feature maps output by an expert network.

Images