Dynamic energy flow calculation method and system for electric heating comprehensive energy system based on physical conservation residual driving

Abstract

The invention provides a dynamic energy flow calculation method and system for an electric heating comprehensive energy system based on physical conservation residual driving. The method comprises the steps that S1, system parameters are input and initialized; s2, initializing time advance and expansion parameters; s3, time embedding full-pure expansion prediction is carried out; s4, calculating a physical conservation residual error; s5, performing judgment based on physical conservation residual errors; s6, performing time scale rollback and sub-time step adaptive reconstruction; s7, system state updating and time advancing; according to the technical scheme, the physical conservation residual error is constructed based on the power grid node power balance relation and the dynamic energy conservation relation of the hot network management section, the satisfaction degree of the system state to the original physical constraint is converted into the numerical index with the clear physical meaning, and a unified and explainable physical consistency evaluation basis is provided for the dynamic energy flow calculation result.

Description

Technical Field

[0001] This invention relates to the field of integrated electrothermal energy technology, and in particular to a method and system for calculating dynamic energy flow in integrated electrothermal energy systems based on physical conservation residuals. Background Technology

[0002] With the development of integrated energy systems, power grids and heating networks form a coordinated operational relationship through energy conversion and dispatch strategies. The system state involves various physical quantities such as voltage, thermal power, and temperature. In planning analysis, operational assessment, and disturbance response studies, the system operating state often changes over time. Steady-state energy flow analysis methods are insufficient to characterize the continuous changes in the system state over time, which will affect the accuracy of operational assessment results.

[0003] Compared to steady-state energy flow, which only reflects the equilibrium state at a certain moment, dynamic energy flow calculation, in addition to satisfying the power balance constraints of power grid nodes, also needs to characterize the temporal evolution of state variables such as temperature and flow rate in the heating network section, and satisfy the corresponding dynamic energy conservation relationship. Therefore, it is necessary to conduct dynamic energy flow calculations on integrated electrothermal energy systems to describe the continuous evolution of system state variables over time.

[0004] For dynamic energy flow calculations, existing research mainly employs time-stepping methods based on numerical integration and holomorphic function methods based on analytical expansion. Numerical integration methods solve the system's differential equations through discrete-time progression, making the results highly sensitive to the chosen time step. In scenarios involving enhanced nonlinearity or long-duration simulations, it's easy to struggle to balance accuracy and stability. Holomorphic expansion methods, by analytically embedding and progressively expanding the system equations, obtain an analytical expression of the system state over time within a single time step, mitigating the sensitivity of numerical iteration to initial conditions to some extent. However, in engineering applications, these methods still face the challenge of online coordination between the expansion order and the time step, leading to error accumulation or local time step calculation failures during enhanced nonlinearity or abrupt changes in operating state, thus affecting the reliability and engineering applicability of the dynamic energy flow analysis results.

[0005] The dynamic energy flow calculation method based on holomorphic expansion can characterize the continuous evolution of the system state in analytical form and alleviate the problem of numerical iteration being sensitive to initial conditions to some extent, but it still has shortcomings.

[0006] Specifically, existing methods in dynamic calculations typically focus on the numerical recursiveness of power series expansions. The validity of the results is primarily evaluated using numerical criteria such as series truncation error, coefficient recursion error, or the difference in order between adjacent expansions. These criteria mainly reflect the numerical stability of the expansion process itself. However, the power balance relationship between power grid nodes and the dynamic energy conservation constraints of heating network segments are mostly used in existing methods for result verification or auxiliary analysis, and have not yet been used as unified physical consistency criteria to directly drive the dynamic energy flow calculation process. Therefore, although the convergence conditions are met numerically, it does not necessarily guarantee that the system state conforms to the original conservation constraints at the physical level, making it difficult to quantify and control the physical deviations of the calculation results.

[0007] Regarding time-mapping strategies, existing methods typically rely on numerical criteria to adjust the expansion order or computational step size separately. These adjustments are mostly centered around a single parameter dimension, and a unified, coordinated adjustment strategy for both the expansion order and the time scale has not yet been established. This approach fails to effectively utilize the consistency of the system's real-time physical state when adjusting expansion accuracy and time scale, leading to physical-level deviations in the calculation results. This, in turn, affects the reliability and applicability of dynamic energy flow analysis under complex operating conditions.

[0008] Furthermore, in long-term dynamic simulation scenarios, existing methods lack effective coordination and adjustment capabilities for the time progression process. When physical deviations accumulate or the effectiveness of the unfolding decreases during the calculation process, existing methods struggle to make timely adaptive adjustments based on the physical consistency of the system, thus limiting their robustness in complex engineering scenarios. Summary of the Invention

[0009] In view of this, the purpose of this invention is to provide a method and system for calculating dynamic energy flow in an integrated electrothermal energy system based on physical conservation residuals. The method constructs physical conservation residuals based on the power balance relationship of power grid nodes and the dynamic energy conservation relationship of heating network pipe sections, and transforms the degree to which the system state satisfies the original physical constraints into numerical indicators with clear physical meanings, providing a unified and interpretable physical consistency evaluation basis for the dynamic energy flow calculation results.

[0010] To achieve the above objectives, the present invention adopts the following technical solution: a method for calculating the dynamic energy flow of an integrated electrothermal energy system based on physical conservation residuals, comprising the following steps:

[0011] S1: Input the system parameters of the power grid, heating network and electrothermal coupling unit, construct the power balance model of the power grid node and the dynamic energy conservation model of the heating network section; set the time advance step size Δt, the initial order of the holomorphic expansion N, the upper limit of the expansion order Nmax and the physical conservation residual threshold ε, and initialize the system state variables;

[0012] S2: At the current time step [t] k , tk +Δt]in which This is the starting time of the current time step. Set the time advance step size; based on the current time step start state, initialize the time advance structure, and use the initial order N of the holomorphic expansion as the starting order of the expansion recursion within this time step, providing initial conditions for the time embedding holomorphic expansion prediction calculation;

[0013] S3: Introduce embedded parameters within the current time step, perform a holomorphic time series expansion on the system state variables, and obtain the prediction result of the system state at the end of the current time step by solving the problem step by step.

[0014] S4: Substitute the predicted state obtained in S3 back into the power balance equation of the power grid node and the dynamic energy conservation equation of the heating network section, calculate the power balance residual of the power grid and the dynamic energy conservation residual of the heating network respectively, and further calculate the comprehensive physical conservation residual index E.

[0015] S5: Compare the comprehensive physical conservation residual E with the physical conservation residual threshold ε:

[0016] When the comprehensive physical conservation residual E ≤ the physical conservation residual threshold ε, proceed to S7;

[0017] When the comprehensive physical conservation residual E > the physical conservation residual threshold ε, the holomorphic expansion order is increased. When the holomorphic expansion order reaches the preset upper limit and the comprehensive physical conservation residual still cannot meet the threshold requirement, proceed to S6.

[0018] S6: When the comprehensive physical conservation residual exceeds the preset threshold and the order of the holomorphic expansion has reached the preset upper limit, it is considered that the current time step Δt has exceeded the effective time scale of the holomorphic expansion in the physical state of the system, and the time advancement structure needs to be adaptively reconstructed; the time scale rollback and sub-time step subdivision calculation are performed on the current time advancement stage.

[0019] S7: Once the calculation results of the current time advance phase or all its sub-time steps have been accepted, update the system state variables and enter the next time advance phase. Repeat the above steps until the dynamic energy flow calculation of the integrated electrothermal energy system is completed.

[0020] In a preferred embodiment, in step S1, the power grid uses a node power balance model to describe the system operating state; for any node in the power grid, its injected active power and reactive power must satisfy the node power balance relationship, expressed as:

[0021] (1)

[0022] In the formula, i, j∈{1,…,N}, and N is the total number of power grid nodes; Inject complex power into node i; and Let i represent the active power and reactive power of node i, respectively. Let i be the voltage at node i. For the complex conjugate of the voltage at node j; Let be the elements of the node admittance matrix between node i and node j. for Complex conjugate.

[0023] In a preferred embodiment, in S1, the heating network in the dynamic energy conservation model of the heating network segment consists of heating pipes, nodes, and heat transfer fluid flowing in the pipes; the flow and heat exchange process of the heat transfer fluid in the pipes follows the law of energy conservation, and its dynamic mechanism is described by a partial differential equation of temperature with respect to time and space; for any heating pipe, its dynamic energy conservation relationship is expressed as:

[0024] (2)

[0025] In the formula, Let x be the temperature at position x inside the pipe at time t, where x is the spatial position variable and t is time. For the density of the heat transfer medium, For specific heat capacity, The cross-sectional area of ​​the pipe. This refers to the mass flow rate of the heat transfer medium within the pipeline. The heat transfer coefficient per unit length of the pipe; Ambient temperature; This represents the rate of temperature change at position x inside the pipe at time t; This represents the rate of change of temperature at position x within the pipe at time t;

[0026] The heating pipeline is discretized along its length, and it is assumed that the mass flow rate remains constant within a single pipe segment. Spatial integration is performed over the length interval of the pipe segment, and an approximation of the average inlet and outlet temperatures is used. The average temperature of the k-th pipe segment is defined as:

[0027] (3)

[0028] In the formula, Let K be the average temperature of pipe segment k at time t. and Let represent the inlet and outlet temperatures of pipe segment k at time t, respectively, with parentheses (t) indicating that the variable changes with time;

[0029] The original partial differential form of the dynamic energy conservation model is transformed into a first-order ordinary differential equation with respect to time:

[0030] (4)

[0031] In the formula, Let K be the mass flow rate of pipe segment k. , These are the cross-sectional area and length of pipe segment k, respectively; Let k be the equivalent heat transfer coefficient of the pipe section; and Let represent the first derivatives of the inlet and outlet temperatures of pipe segment k with respect to time t, respectively.

[0032] In addition, the heat network nodes satisfy the mass continuity constraint relationship, and the inlet temperature of the pipe section is consistent with the starting node temperature. The above relationships together constitute the dynamic energy conservation model of the heat network pipe section.

[0033] In a preferred embodiment, S3 specifically includes: introducing an embedding parameter z. [0,1]; The embedded parameter z is used to characterize the process of the system evolving from the current physical state to the next physical state in a single time step. Its change reflects the continuous evolution characteristics of the system state variables in the time process.

[0034] In a preferred embodiment, S3 specifically includes: unifying system state variables such as grid node voltage, heating network section temperature, and mass flow rate into a state vector X; and, within the time embedding framework, representing the system state vector as a holomorphic function of the embedding parameter z.

[0035] (5)

[0036] in, This represents a holomorphic function of the system state vector with respect to the embedding parameter z. Let n be the holomorphic expansion coefficients of the system state vector with respect to the embedded parameter z, where n is the index of the expansion order. The nth power of the embedded parameter z represents the basis function term that constitutes the power series expansion;

[0037] Based on the embedding relation (5), the derivative of the system state variable with respect to the embedding parameter z is expressed as:

[0038] (6)

[0039] in, This represents the rate of change of the system state vector X with respect to the embedding parameter z; it combines the embedding parameter z with the physical time advance step size. The mapping relationship is used to construct the recursive update relationship of the state variables in the time process, thereby realizing the stepwise solution of dynamic energy flow;

[0040] Substituting equations (5) and (6) into the power balance model of power grid nodes and the dynamic energy conservation model of heating network sections, and by adjusting the embedded parameters... By matching terms of the same power order, holomorphic expansion coefficients can be established. The recursive relationship between them; following the order of expansion from low to high, the holomorphic expansion coefficients of the system state are solved step by step to obtain the analytical expression of the system state in the current time step;

[0041] After completing the holomorphic expansion recursive solution within the current time step, the system state at the end of the current time progression phase is predicted by setting the embedding parameter z=1. This prediction result is used to characterize the state of the system after completing the physical evolution within a time step. Subsequently, the predicted state is used as the initial state of the next time step, and the above time embedding and holomorphic expansion recursive solution process is repeated until the dynamic energy flow calculation of the electrothermal integrated energy system is completed throughout the entire simulation cycle.

[0042] In a preferred embodiment, S4 specifically includes: substituting the predicted voltage of grid node i into equation (1) to construct the injected power residual of grid node i:

[0043]

[0044] in, Let be the residual of the injected power at grid node i; For a given complex power injection at node i; symbol This indicates taking the complex modulus, i.e., the amplitude, which is used to convert the complex power imbalance into a comparable scalar residual, and is used as a criterion for subsequent physical consistency.

[0045] Substituting the predicted state of the heating network into the dynamic energy conservation model of the heating network, the energy conservation residual of the k-th pipe segment is constructed. for:

[0046] (8)

[0047] Where k is the pipe section number; and These are the inlet and outlet temperatures of pipe section k, respectively. and These represent the rates of change of the inlet and outlet temperatures of pipe segment k with respect to time t, respectively.

[0048] Define the physical conservation residual of the system synthesis This represents the maximum value of the residuals of the power grid power balance and the residuals of the heating network energy conservation:

[0049] (9)

[0050] in, Let i be the power balance residual of the power grid. The residual for energy conservation in the heating network of pipe segment k. This is the power grid reference value. It serves as a benchmark value for the thermal power of the heating network, used to normalize the residuals of the power balance of the power grid and the energy conservation of the heating network, so as to achieve a comprehensive residual evaluation with consistent dimensions.

[0051] In a preferred embodiment, S5 specifically includes: setting a physical conservation error threshold ε, and using the comparison result between the comprehensive physical conservation residual E and the physical conservation error threshold ε as the judgment basis and triggering condition for adaptive adjustment;

[0052] When the comprehensive physical conservation residual E ≤ physical conservation error ε, the calculation result in the current time step is determined to meet the system physical consistency requirements, and the calculation result of the time step is accepted.

[0053] When the comprehensive physical conservation residual E > the physical conservation error ε, the adaptive adjustment process is triggered.

[0054] In a preferred embodiment, triggering the adaptive adjustment process includes:

[0055] Step 1: Roll back the system state to the start time t of the current time step. k Cancel any predictions that have not yet been accepted within the current time step and restore the initial state of that time progression phase;

[0056] Step 2: Divide the original time step Δt evenly into M sub-time steps according to the definition of sub-time step length. The time scale of each sub-time step is represented by Δt. (m) And initialize the time advancement structure and full-pure expansion parameters within the sub-time step based on the sub-time scale;

[0057] Step 3: Within each sub-time step, reintroduce the embedded parameters according to the sub-time progression sequence, perform holomorphic expansion prediction on the system state variables, and calculate the corresponding comprehensive physical conservation residuals; when the residuals meet the preset threshold requirements, accept the calculation results of that sub-time step;

[0058] Step 4: After the calculation result of the current sub-time step is accepted, advance the system state to the next sub-time step and repeat the prediction and decision process within the sub-time step until the calculation of all M sub-time steps is completed;

[0059] Step 5: When the comprehensive physical conservation residual cannot meet the preset threshold requirement under the current sub-time scale, it is considered that the sub-time scale is still insufficient to support the physical consistency description of the system state. The number of sub-time step divisions M is incremented, and the time advancement structure is reconstructed based on the updated time scale. The above sub-time step subdivision calculation process is repeated.

[0060] Step 6: After completing all sub-time step calculations and meeting physical consistency requirements, the system state is advanced to the end time t of the original time step. k +Δt, and return to the main time advance process to continue the dynamic energy flow calculation for subsequent time steps.

[0061] In a preferred embodiment, step 2 involves changing the original time step. The time step is divided into M sub-time steps, and the step size of the m-th sub-time step is defined as follows:

[0062] (10)

[0063] This constitutes the sub-time progression sequence in step 3:

[0064] (11)

[0065] Where m is the sub-time step number; M is the number of sub-time step divisions; This represents the original time step for this time progression phase; Let m be the step size of the m-th sub-time step; It indicates the starting moment of the current time progression phase.

[0066] This invention also provides a dynamic energy flow calculation system for an integrated electrothermal energy system based on physical conservation residuals, comprising a processor, a memory, and a bus. The memory stores machine-readable instructions executed by the processor. When the system is running, the processor communicates with the memory via the bus, and when the machine-readable instructions are executed by the processor, the dynamic energy flow calculation method for an integrated electrothermal energy system based on physical conservation residuals is as described above.

[0067] Compared with the prior art, the present invention has the following beneficial effects:

[0068] (1) Based on the power balance relationship of the power grid nodes and the dynamic energy conservation relationship of the heat network pipe section, a physical conservation residual is constructed, which transforms the degree of satisfaction of the system state with the original physical constraints into a numerical index with clear physical meaning, providing a unified and interpretable physical consistency evaluation basis for the dynamic energy flow calculation results.

[0069] (2) By introducing the physical conservation residual into the dynamic energy flow calculation control process, the coordinated adaptive adjustment of the order of the whole pure expansion and the time advancement step size is realized. While ensuring the physical consistency of the system, the calculation efficiency is taken into account, and the calculation instability or redundancy caused by relying on the adjustment of a single parameter is avoided.

[0070] (3) When the expansion order adjustment cannot meet the physical consistency requirements, a time rollback mechanism based on physical conservation residual feedback is introduced, and the time step is incrementally subdivided. In this mechanism, when the physical conservation residual exceeds a preset threshold, the system rolls back to the previous calculation state and increments the number of sub-time step divisions according to the residual feedback, thereby subdividing the time step. This adjustment ensures the continuous evolution characteristics of the system state and improves the stability and reliability of dynamic energy flow calculation when facing strongly nonlinear and long-term simulation scenarios.

[0071] (4) By maintaining the continuous evolution of the system state between adjacent time stages during the dynamic calculation process, it is beneficial to support the long-term dynamic simulation and operation analysis needs of the integrated electric and thermal energy system. Attached Figure Description

[0072] Figure 1 The main flowchart of the dynamic energy flow calculation method for an electrothermal integrated energy system based on physical conservation residuals according to a preferred embodiment of the present invention is shown below.

[0073] Figure 2 The flowchart of time-scale rollback and sub-time step subdivision based on physical conservation residuals is a preferred embodiment of the present invention. Detailed Implementation

[0074] The present invention will be further described below with reference to the accompanying drawings and embodiments.

[0075] It should be noted that the following detailed descriptions are illustrative and intended to provide further explanation of this application. Unless otherwise specified, all technical and scientific terms used herein have the same meaning as commonly understood by one of ordinary skill in the art to which this application pertains.

[0076] It should be noted that the terminology used herein is for the purpose of describing particular implementations only and is not intended to limit the exemplary implementations according to this application; as used herein, the singular form is intended to include the plural form as well, unless the context clearly indicates otherwise; furthermore, it should be understood that when the terms “comprising” and / or “including” are used in this specification, they indicate the presence of features, steps, operations, devices, components and / or combinations thereof.

[0077] A method for calculating dynamic energy flow in an integrated electrothermal energy system based on physical conservation residuals, referencing Figure 1-2 When physical consistency fails during time progression, the timescale rollback and sub-time step subdivision process is as follows: Figure 2 As shown in the figure, E is the comprehensive physical conservation residual, ε is the physical conservation residual threshold, N is the holomorphic expansion order, and M is the time step subdivision number; the following steps are included:

[0078] S1: System Parameter Input and Initialization

[0079] Input the system parameters of the power grid, heating network and electrothermal coupling unit, construct the power balance model of the power grid node and the dynamic energy conservation model of the heating network section; set the time step Δt, the initial order N of the holomorphic expansion, the upper limit of the expansion order Nmax and the physical conservation residual threshold ε, and initialize the system state variables.

[0080] S2: Time Advancement and Deployment Parameter Initialization

[0081] At the current time step [t] k , t k Within +Δt], based on the initial state of the current time step, the time advancement structure is initialized, and the initial order N of the holomorphic expansion is used as the starting order of the expansion recursion within this time step, providing initial conditions for the time embedding holomorphic expansion prediction calculation.

[0082] S3: Temporal Embedded Holomorphic Unfolding Prediction

[0083] Introduce embedded parameters within the current time step, and perform a holomorphic time series expansion of the system state variables (as shown in the equation). As shown in the figure, the prediction result of the system state at the end of the current time step is obtained by solving the problem step by step.

[0084] S4: Physical Conservation Residual Calculation

[0085] Substitute the predicted state obtained from S3 back into the power balance equation of the grid node (as shown in equation S3). (as shown) and the dynamic energy conservation equation of the heating network pipe section (as shown in equation) As shown), calculate the power balance residuals of the power grid respectively (as shown in the formula). (as shown) and the dynamic energy conservation residual of the heating network (as shown in the formula) (as shown), and further calculate the comprehensive physical conservation residual index E (as shown in the formula). (As shown).

[0086] S5: Judgment based on physical conservation residuals

[0087] The comprehensive physical conservation residual E is compared with the preset threshold ε:

[0088] When E≤ε, the calculation results of the current time advancement stage are considered to meet the system physical consistency requirements, and the process proceeds to step S7;

[0089] When E>ε, it is determined that the current calculation result does not meet the physical consistency requirement. Prioritize improving the time expansion accuracy by increasing the order of the holomorphic expansion while keeping the time advancement structure unchanged. When the expansion order reaches the preset upper limit and the comprehensive physical conservation residual still cannot meet the threshold requirement, it is considered that the current time step has exceeded the effective time scale of the holomorphic expansion in the system state, and proceed to step S6.

[0090] S6: Time Scale Rollback and Sub-Time Step Adaptive Reconstruction

[0091] When the comprehensive physical conservation residual exceeds a preset threshold and the holomorphic expansion order reaches a preset upper limit, it is considered that the current time step Δt has exceeded the effective time scale of the holomorphic expansion in the physical state of the system, and an adaptive reconstruction of the time-progression structure is required. At this time, time scale rollback and sub-time step subdivision calculations are performed on the current time-progression stage, the process of which includes:

[0092] Step 1: Roll back the system state to the start time t of the current time step. k Cancel any predictions that have not yet been accepted within the current time step and restore the initial state of that time progression phase;

[0093] Step 2: Define the sub-time step length according to the formula (e.g.) As shown in the figure, the original time step Δt is evenly divided into M sub-time steps, and the time scale of each sub-time step is represented by Δt. (m) And initialize the time advancement structure and full-pure expansion parameters within the sub-time step based on the sub-time scale;

[0094] Step 3: Within each sub-time step, proceed according to the sub-time progression sequence (as shown in the formula). (As shown) The embedded parameters are reintroduced, and holomorphic expansion prediction is performed on the system state variables. The corresponding integrated physical conservation residuals are calculated. When the residuals meet the preset threshold requirements, the calculation results of the sub-time step are accepted.

[0095] Step 4: After the calculation result of the current sub-time step is accepted, advance the system state to the next sub-time step and repeat the prediction and decision process within the sub-time step until the calculation of all M sub-time steps is completed;

[0096] Step 5: When the comprehensive physical conservation residual cannot meet the preset threshold requirement under the current sub-time scale, it is considered that the sub-time scale is still insufficient to support the physical consistency description of the system state. The number of sub-time step divisions M is incremented, and the time advancement structure is reconstructed based on the updated time scale. The above sub-time step subdivision calculation process is repeated.

[0097] Step 6: After completing all sub-time step calculations and meeting physical consistency requirements, the system state is advanced to the end time t of the original time step. k +Δt, and return to the main time advance process to continue the dynamic energy flow calculation for subsequent time steps.

[0098] S7: System Status Updates and Time Progression

[0099] Once the calculation results of the current time advance phase or all its sub-time steps are accepted, the system state variables are updated, and the next time advance phase is entered. The above steps are repeated until the dynamic energy flow calculation of the integrated electrothermal energy system is completed.

[0100] Specifically, in S1, in order to achieve unified modeling and solution of dynamic energy flow in the integrated electric and thermal energy system, a power balance model of the power grid node and a dynamic energy conservation model of the thermal network section are constructed respectively. These are used to form the original physical constraint equations that the system state variables must satisfy during dynamic operation, providing a foundation for subsequent holomorphic expansion prediction and physical conservation residual construction.

[0101] The power grid section uses a nodal power balance model to describe the system's operating state. For any node in the power grid, the injected active and reactive power must satisfy the nodal power balance relationship, expressed as:

[0102]

[0103] In the formula, i, j∈{1,…,N}, and N is the total number of power grid nodes; Inject complex power into node i; and Let i represent the active power and reactive power of node i, respectively. Let i be the voltage at node i. For the complex conjugate of the voltage at node j; Let be the elements of the node admittance matrix between node i and node j. for Complex conjugate.

[0104] A heating network consists of heating pipes, nodes, and the heat transfer medium flowing within the pipes. The flow and heat exchange process of the heat transfer medium in the pipes follows the law of energy conservation, and its dynamic mechanism can be described by partial differential equations of temperature with respect to time and space. For any heating pipe, its dynamic energy conservation relationship can be expressed as:

[0105]

[0106] In the formula, Let x be the temperature at position x inside the pipe at time t, where x is the spatial position variable and t is time. For the density of the heat transfer medium, For specific heat capacity, The cross-sectional area of ​​the pipe. This refers to the mass flow rate of the heat transfer medium within the pipeline. The heat transfer coefficient per unit length of the pipe; Ambient temperature; This represents the rate of temperature change at position x inside the pipe at time t; This represents the rate of change of temperature at position x at time t inside the pipe.

[0107] To accommodate dynamic energy flow calculations, the heating pipeline is discretized along its length, and it is assumed that the mass flow rate remains constant within a single pipe segment. Spatial integration is performed over the length interval of the pipe segment, and an approximation of the average inlet and outlet temperatures is used. The average temperature of the k-th pipe segment can be defined as:

[0108]

[0109] In the formula, Let K be the average temperature of pipe segment k at time t. and Let represent the inlet and outlet temperatures of pipe segment k at time t, respectively, with parentheses (t) indicating that the variable changes over time.

[0110] Further simplification can transform the original partial differential form of the dynamic energy conservation model into a first-order ordinary differential equation with respect to time:

[0111]

[0112] In the formula, Let K be the mass flow rate of pipe segment k. , These are the cross-sectional area and length of pipe segment k, respectively; Let k be the equivalent heat transfer coefficient of the pipe section; and Let represent the first derivatives of the inlet and outlet temperatures of pipe segment k with respect to time t, respectively.

[0113] In addition, the heat network nodes satisfy the mass continuity constraint relationship, and the inlet temperature of the pipe section is consistent with the starting node temperature. The above relationships together constitute the physical constraint model of the heat network.

[0114] The electrothermal coupling unit connects the power grid and the heating network through energy conversion relationships. Its operating characteristics are described by the physical constraint relationship between electrical power and thermal power, and are uniformly incorporated into the system's physical constraint framework to participate in dynamic energy flow calculations.

[0115] Specifically, in S3, let the current system time step be... ,in This is the starting time of the current time step. Let z be the time step size. To describe the continuous process of the system state evolving from the current time step to the next time step, an embedded parameter z is introduced. [0,1]. The embedding parameter z is used to characterize the process of the system evolving from the current physical state to the next physical state within a single time step. Its changes reflect the continuous evolution characteristics of the system's state variables during the time progression. Although the embedding parameter z itself is not directly equivalent to physical time, its changes correspond to the actual physical evolution process of the system within a time progression step, and therefore have a clear physical meaning.

[0116] System state variables such as grid node voltage, heating network section temperature, and mass flow rate are uniformly denoted as state vector X. Within the time embedding framework, the system state vector is represented as a holomorphic function of the embedding parameter z:

[0117]

[0118] in, This represents a holomorphic function of the system state vector with respect to the embedding parameter z. Let n be the holomorphic expansion coefficients of the system state vector with respect to the embedded parameter z, where n is the index of the expansion order. Let z be the nth power of the embedded parameter z, which is a basis function term that constitutes the power series expansion.

[0119] Based on embedded relation The derivative of the system state variable with respect to the embedded parameter z can be expressed as:

[0120]

[0121] in, This represents the rate of change of the system state vector X with respect to the embedding parameter z. It combines the embedding parameter z with the physical time advance step size. The mapping relationship can be used to construct the recursive update relationship of the state variables in the time process, thereby realizing the stepwise solution of dynamic energy flow.

[0122] The formula Japanese style Substituting these parameters into the power balance model of power grid nodes and the dynamic energy conservation model of heating network sections, the embedded parameters are then analyzed. By matching terms of the same power order, holomorphic expansion coefficients can be established. The recursive relationship between them is established. Following the order of expansion from low to high, the holomorphic expansion coefficients of the system state are solved step by step to obtain the analytical expression of the system state in the current time step.

[0123] After completing the holomorphic expansion recursive solution within the current time step, the system state at the end of the current time progression phase is predicted by setting the embedding parameter z=1. This result characterizes the system state after completing its physical evolution within a time step. Subsequently, this predicted state is used as the initial state for the next time step, and the above time embedding and holomorphic expansion recursive solution process is repeated until the dynamic energy flow calculation of the electrothermal integrated energy system is completed throughout the entire simulation cycle.

[0124] Specifically, in S4, to avoid physical distortion caused by relying solely on numerical convergence criteria, physical conservation residuals are constructed based on the system's physical model and used to verify the physical consistency of dynamic energy flow calculation results.

[0125] Substitute the predicted voltage of grid node i into the equation Construct the injected power residual of grid node i:

[0126]

[0127] in, Let be the residual of the injected power at grid node i; For a given complex power injection at node i; symbol This indicates taking the complex modulus, i.e., the amplitude, which is used to convert the complex power imbalance into a comparable scalar residual, used as a criterion for subsequent physical consistency.

[0128] Substitute the predicted state results of the heating network into the dynamic energy conservation model of the heating network. Construct the energy conservation residual for the k-th pipe segment. for:

[0129]

[0130] Where k is the pipe section number; and These are the inlet and outlet temperatures of pipe section k, respectively. and These represent the rates of change of the inlet and outlet temperatures of pipe segment k with respect to time t, respectively.

[0131] Define the physical conservation residual of the system synthesis This represents the maximum value of the residuals of the power grid power balance and the residuals of the heating network energy conservation:

[0132]

[0133] in, Let i be the power balance residual of the power grid. The residual for energy conservation in the heating network of pipe segment k. This is the power grid reference value. It serves as a benchmark value for the thermal power of the heating network, used to normalize the residuals of the power balance of the power grid and the energy conservation of the heating network, so as to achieve a comprehensive residual evaluation with consistent dimensions.

[0134] Specifically, in S5, an error threshold ε is set, preferably ε = 1 × 10⁻ 4 and in the form The comparison between the defined normalized integrated physical conservation residual E and ε serves as the judgment basis and triggering condition for adaptive adjustment.

[0135] 1. When E≤ε, determine that the calculation result in the current time step meets the system physical consistency requirements, and accept the calculation result of this time step;

[0136] 2. When E>ε, the adaptive adjustment process is triggered.

[0137] Specifically, in S6, the adaptive adjustment process includes:

[0138] (1) Adaptive adjustment of expansion order

[0139] Maintain the current time step Under the condition that remains unchanged, the order N of the holomorphic time series expansion is increased, and the holomorphic expansion prediction and physical conservation residual calculation process in the current time step are re-executed to enhance the accuracy of the time expansion.

[0140] (2) Time step rollback and subdivision control

[0141] When the holomorphic expansion order reaches the preset upper limit Nmax, and the physical conservation residual still satisfies E>ε, the current time step is considered to be... The effective timescale for a fully pure expansion within this system state has been exceeded. At this point, a rollback operation is performed on the current time progression phase, restoring the system state to its initial state at the current time step, and adjusting the time step size. Perform further subdivision processing.

[0142] Specifically, the original time step The time step is divided into M sub-time steps, and the step size of the m-th sub-time step is defined as follows:

[0143]

[0144] This constitutes a new time-progression sequence:

[0145]

[0146] Where m is the sub-time step number; M is the number of sub-time step divisions; This represents the original time step for this time progression phase; Let m be the step size of the m-th sub-time step; It indicates the starting moment of the current time progression phase.

[0147] Within each sub-time step, the time embedding parameters are reintroduced, and a fully virtually expanded prediction calculation based on time embedding is performed to obtain the system state prediction value at the end of the sub-time step. Subsequently, this prediction value is substituted back into the power balance equation of the power grid nodes and the dynamic energy conservation equation of the heating network segment to calculate the corresponding physical conservation residuals. Based on the physical conservation residual judgment results, the expansion order adjustment is performed first. When the expansion order adjustment still cannot meet the physical consistency requirements, the adaptive reconstruction of the sub-time step scale is further performed. The number of sub-time step divisions M can start from the initial value M=2, and is incrementally updated in the manner of M=M+1 if the residual still does not meet the threshold requirements, until the accuracy requirements are met.

[0148] Specifically, in S7, when the calculation result of the current sub-time step meets the physical conservation residual threshold requirement, the final state of that sub-time step is accepted as the initial state of the next sub-time step, and subsequent sub-time steps are continued; after all M sub-time steps are completed, the system state is advanced to the end of the original time step. Then, return to the main time step and proceed to the dynamic energy flow calculation for the next time step.

Claims

1. A method for calculating the dynamic energy flow of an integrated electrothermal energy system based on physical conservation residuals, characterized in that, Includes the following steps: S1: Input the system parameters of the power grid, heating network and electrothermal coupling unit, construct the power balance model of the power grid node and the dynamic energy conservation model of the heating network section; set the time advance step size Δt, the initial order of the holomorphic expansion N, the upper limit of the expansion order Nmax and the physical conservation residual threshold ε, and initialize the system state variables; S2: At the current time step [t] k , t k +Δt]in which This is the starting time of the current time step. Set the time advance step size; based on the current time step start state, initialize the time advance structure, and use the initial order N of the holomorphic expansion as the starting order of the expansion recursion within this time step, providing initial conditions for the time embedding holomorphic expansion prediction calculation; S3: Introduce embedded parameters within the current time step, perform a holomorphic time series expansion on the system state variables, and obtain the prediction result of the system state at the end of the current time step by solving the problem step by step. S4: Substitute the predicted state obtained in S3 back into the power balance equation of the power grid node and the dynamic energy conservation equation of the heating network section, calculate the power balance residual of the power grid and the dynamic energy conservation residual of the heating network respectively, and further calculate the comprehensive physical conservation residual index E. S5: Compare the comprehensive physical conservation residual E with the physical conservation residual threshold ε: When the comprehensive physical conservation residual E ≤ the physical conservation residual threshold ε, proceed to S7; When the comprehensive physical conservation residual E > the physical conservation residual threshold ε, the holomorphic expansion order is increased. When the holomorphic expansion order reaches the preset upper limit and the comprehensive physical conservation residual still cannot meet the threshold requirement, proceed to S6. S6: When the comprehensive physical conservation residual exceeds the preset threshold and the order of the holomorphic expansion has reached the preset upper limit, it is considered that the current time step Δt has exceeded the effective time scale of the holomorphic expansion in the physical state of the system, and the time advancement structure needs to be adaptively reconstructed; the time scale rollback and sub-time step subdivision calculation are performed on the current time advancement stage. S7: Once the calculation results of the current time advance phase or all its sub-time steps have been accepted, update the system state variables and enter the next time advance phase. Repeat the above steps until the dynamic energy flow calculation of the integrated electrothermal energy system is completed.

2. The method for calculating dynamic energy flow in an integrated electrothermal energy system based on physical conservation residuals as described in claim 1, characterized in that, In S1, the power grid uses a node power balance model to describe the system operating state; for any node in the power grid, its injected active power and reactive power must satisfy the node power balance relationship, expressed as: （1） In the formula, i, j∈{1,…,N}, and N is the total number of power grid nodes; Inject complex power into node i; and Let i represent the active power and reactive power of node i, respectively. Let i be the voltage at node i. For the complex conjugate of the voltage at node j; Let be the elements of the node admittance matrix between node i and node j. for Complex conjugate.

3. The method for calculating dynamic energy flow in an electrothermal integrated energy system based on physical conservation residuals as described in claim 1, characterized in that, In S1, the heating network in the dynamic energy conservation model of the heating network section consists of heating pipes, nodes, and heat transfer fluid flowing in the pipes; the flow and heat exchange process of the heat transfer fluid in the pipes follows the law of energy conservation, and its dynamic mechanism is described by partial differential equations of temperature with respect to time and space; for any heating pipe, its dynamic energy conservation relationship is expressed as: （2） In the formula, Let x be the temperature at position x inside the pipe at time t, where x is the spatial position variable and t is time. For the density of the heat transfer medium, For specific heat capacity, The cross-sectional area of ​​the pipe. This refers to the mass flow rate of the heat transfer medium within the pipeline. The heat transfer coefficient per unit length of the pipe; Ambient temperature; This represents the rate of temperature change at position x inside the pipe at time t; This represents the rate of change of temperature at position x within the pipe at time t; The heating pipeline is discretized along its length, and it is assumed that the mass flow rate remains constant within a single pipe segment. Spatial integration is performed over the length interval of the pipe segment, and an approximation of the average inlet and outlet temperatures is used. The average temperature of the k-th pipe segment is defined as: （3） In the formula, Let K be the average temperature of pipe segment k at time t. and Let represent the inlet and outlet temperatures of pipe segment k at time t, respectively, with parentheses (t) indicating that the variable changes with time; The original partial differential form of the dynamic energy conservation model is transformed into a first-order ordinary differential equation with respect to time: （4） In the formula, Let K be the mass flow rate of pipe segment k. , These are the cross-sectional area and length of pipe segment k, respectively; Let k be the equivalent heat transfer coefficient of the pipe section; and Let represent the first derivatives of the inlet and outlet temperatures of pipe segment k with respect to time t, respectively. In addition, the heat network nodes satisfy the mass continuity constraint relationship, and the inlet temperature of the pipe section is consistent with the starting node temperature. The above relationships together constitute the dynamic energy conservation model of the heat network pipe section.

4. The method for calculating dynamic energy flow in an electrothermal integrated energy system based on physical conservation residuals as described in claim 1, characterized in that, S3 specifically includes: introducing an embedding parameter z [0,1]; The embedded parameter z is used to characterize the process of the system evolving from the current physical state to the next physical state in a single time step. Its change reflects the continuous evolution characteristics of the system state variables in the time process.

5. The method for calculating dynamic energy flow in an electrothermal integrated energy system based on physical conservation residuals as described in claim 4, characterized in that, Specifically, S3 includes: unifying system state variables such as grid node voltage, heating network section temperature, and mass flow rate into a state vector X; and, within the time embedding framework, representing the system state vector as a holomorphic function of the embedding parameter z. （5） in, This represents a holomorphic function of the system state vector with respect to the embedding parameter z. Let n be the holomorphic expansion coefficients of the system state vector with respect to the embedded parameter z, where n is the index of the expansion order. The nth power of the embedded parameter z represents the basis function term that constitutes the power series expansion; Based on the embedding relation (5), the derivative of the system state variable with respect to the embedding parameter z is expressed as: （6） in, This represents the rate of change of the system state vector X with respect to the embedding parameter z; it combines the embedding parameter z with the physical time advance step size. The mapping relationship is used to construct the recursive update relationship of the state variables in the time process, thereby realizing the stepwise solution of dynamic energy flow; Substituting equations (5) and (6) into the power balance model of power grid nodes and the dynamic energy conservation model of heating network sections, and by adjusting the embedded parameters... By matching terms of the same power order, holomorphic expansion coefficients can be established. The recursive relationship between them; following the order of expansion from low to high, the holomorphic expansion coefficients of the system state are solved step by step to obtain the analytical expression of the system state in the current time step; After completing the holomorphic expansion recursive solution within the current time step, the system state at the end of the current time progression phase is predicted by setting the embedding parameter z=1. This prediction result is used to characterize the state of the system after completing the physical evolution within a time step. Subsequently, the predicted state is used as the initial state of the next time step, and the above time embedding and holomorphic expansion recursive solution process is repeated until the dynamic energy flow calculation of the electrothermal integrated energy system is completed throughout the entire simulation cycle.

6. The method for calculating dynamic energy flow in an electrothermal integrated energy system based on physical conservation residuals as described in claim 2, characterized in that, Specifically, S4 includes: substituting the predicted voltage of grid node i into equation (1) to construct the injected power residual of grid node i: in, Let be the residual of the injected power at grid node i; For a given complex power injection at node i; symbol This indicates taking the complex modulus, i.e., the amplitude, which is used to convert the complex power imbalance into a comparable scalar residual, and is used as a criterion for subsequent physical consistency. Substituting the predicted state of the heating network into the dynamic energy conservation model of the heating network, the energy conservation residual of the k-th pipe segment is constructed. for: （8） Where k is the pipe section number; and These are the inlet and outlet temperatures of pipe section k, respectively. and These represent the rates of change of the inlet and outlet temperatures of pipe segment k with respect to time t, respectively. Define the physical conservation residual of the system synthesis This represents the maximum value of the residuals of the power grid power balance and the residuals of the heating network energy conservation: （9） in, Let i be the power balance residual of the power grid. The residual for energy conservation in the heating network of pipe segment k. This is the power grid reference value. It serves as a benchmark value for the thermal power of the heating network, used to normalize the residuals of the power balance of the power grid and the energy conservation of the heating network, so as to achieve a comprehensive residual evaluation with consistent dimensions.

7. The method for calculating dynamic energy flow in an integrated electrothermal energy system based on physical conservation residuals as described in claim 1, characterized in that, S5 specifically includes: setting a physical conservation error threshold ε, and using the comparison result between the comprehensive physical conservation residual E and the physical conservation error threshold ε as the judgment basis and triggering condition for adaptive adjustment; When the comprehensive physical conservation residual E ≤ physical conservation error ε, the calculation result in the current time step is determined to meet the system physical consistency requirements, and the calculation result of the time step is accepted. When the comprehensive physical conservation residual E > the physical conservation error ε, the adaptive adjustment process is triggered.

8. The method for calculating dynamic energy flow in an electrothermal integrated energy system based on physical conservation residuals as described in claim 7, characterized in that, Triggering the adaptive adjustment process includes: Step 1: Roll back the system state to the start time t of the current time step. k Cancel any predictions that have not yet been accepted within the current time step and restore the initial state of that time progression phase; Step 2: Divide the original time step Δt evenly into M sub-time steps according to the definition of sub-time step length. The time scale of each sub-time step is represented by Δt. (m) And initialize the time advancement structure and full-pure expansion parameters within the sub-time step based on the sub-time scale; Step 3: Within each sub-time step, reintroduce the embedded parameters according to the sub-time progression sequence, perform holomorphic expansion prediction on the system state variables, and calculate the corresponding comprehensive physical conservation residuals; when the residuals meet the preset threshold requirements, accept the calculation results of that sub-time step; Step 4: After the calculation result of the current sub-time step is accepted, advance the system state to the next sub-time step and repeat the prediction and decision process within the sub-time step until the calculation of all M sub-time steps is completed; Step 5: When the comprehensive physical conservation residual cannot meet the preset threshold requirement under the current sub-time scale, it is considered that the sub-time scale is still insufficient to support the physical consistency description of the system state. The number of sub-time step divisions M is incremented, and the time advancement structure is reconstructed based on the updated time scale. The above sub-time step subdivision calculation process is repeated. Step 6: After completing all sub-time step calculations and meeting physical consistency requirements, the system state is advanced to the end time t of the original time step. k +Δt, and return to the main time advance process to continue the dynamic energy flow calculation for subsequent time steps.

9. The method for calculating dynamic energy flow in an electrothermal integrated energy system based on physical conservation residuals as described in claim 8, characterized in that, In step 2, the original time step The time step is divided into M sub-time steps, and the step size of the m-th sub-time step is defined as follows: （10） This constitutes the sub-time progression sequence in step 3: （11） Where m is the sub-time step number; M is the number of sub-time step divisions; This represents the original time step for this time progression phase; Let m be the step size of the m-th sub-time step; It indicates the starting moment of the current time progression phase.

10. A dynamic energy flow calculation system for an integrated electrothermal energy system based on physical conservation residuals, comprising a processor, a memory, and a bus, wherein the memory stores machine-readable instructions executed by the processor; characterized in that, When the system is running, the processor and the memory communicate via a bus, and the machine-readable instructions are executed by the processor as described in any one of claims 1 to 9. This is a dynamic energy flow calculation method for an electrothermal integrated energy system based on physical conservation residual drive.

Images