Grain storage planning method and system based on multi-source data perception

Abstract

The invention relates to the technical field of grain storage, and discloses a grain storage planning method and system based on multi-source data perception, and the method comprises the following steps: synchronously collecting and processing grain pile temperature and environmental meteorological data; the method comprises the following steps: mapping a grain pile into a three-dimensional voxel grid, dividing a thermal inert region and an unstable thermal active region, and calculating an environmental thermal intrusion amount; step ventilation excitation is applied, the cold air wavefront arrival time is determined, and the equivalent permeability is inverted; reconstructing a flow velocity field, reversely solving an unsteady flow thermal coupling energy equation, and separating a biological heat production rate; calculating a dynamic maintenance entropy; and constructing a multi-target operation model, generating a storage planning strategy and issuing the storage planning strategy to an execution mechanism. According to the method, by constructing the boundary heat effect decoupling model and the flow heat coupling inversion mechanism, the endogenous biological heat production rate in the grain pile is separated from environment heat invasion, quantitative positioning of hidden insect and mildew activities is achieved, and an abnormal heat source can be locked when the local temperature does not obviously rise.

Description

Technical Field

[0001] This invention relates to the field of grain storage technology, specifically to a grain storage planning method and system based on multi-source data sensing. Background Technology

[0002] Food security is the cornerstone of social stability, and loss reduction technologies in the storage process are a crucial line of defense for ensuring the quantity and quality of grain. During long-term storage, grain piles are highly susceptible to environmental heat and humidity transfer, as well as their own biological respiration, which can induce localized mold growth, pest infestation, or condensation and dampness, posing safety hazards. To delay grain quality deterioration and reduce losses, modern grain depot management relies on precise monitoring methods and scientific control strategies to perceive and intervene in the temperature and humidity environment within the warehouse in real time, thereby achieving green and low-carbon operations while ensuring grain storage safety.

[0003] Standardized flat or shallow circular grain storage warehouses are typically equipped with digital grain condition monitoring systems and mechanical ventilation facilities. Through matrix-style temperature-measuring cables embedded within the grain pile, the system can collect grain temperature data at different spatial locations, providing a clear picture of the overall temperature level of the grain pile. The existing control logic primarily compares the monitored real-time temperature values ​​with preset safe grain storage standards. If the grain temperature exceeds the limit, it automatically or manually activates the axial flow and centrifugal fans. This temperature threshold-triggered ventilation mechanism utilizes cool outside air to replace the hot and humid air inside the warehouse, playing a positive role in suppressing the overall rise in grain temperature, balancing temperature differences within the warehouse, and preventing condensation, thus improving the standardization of grain management.

[0004] However, the aforementioned existing technologies still have logical limitations when dealing with complex dynamic grain conditions. First, the monitoring data mainly relies on a single physical quantity of temperature, making it difficult to effectively identify the source of heat. Heat conduction from solar radiation on the silo walls often masks weak early biological heating signals deep within the grain pile, causing the system to be unable to distinguish between external environmental interference and internal insect and mold activity, easily leading to false alarms or missed alarms. Second, conventional ventilation control assumes that the pore distribution of the grain pile is uniform, ignoring the high flow resistance characteristics often present in moldy and compacted areas. Airflow tends to flow along low-resistance paths, making it easy for cooling airflow to bypass key heating centers. This results in grain in normal areas losing moisture and weight due to excessive ventilation, while the core diseased areas remain at high temperatures due to airflow short-circuiting. In addition, the existing decision-making model lacks a quantitative assessment of the economics of operation and maintenance. When faced with localized stubborn hotspots, the system often mechanically maintains high-power operation without considering whether the invested electricity costs exceed the value preservation benefits of the grain itself, making it difficult to achieve an optimal balance between energy efficiency and losses. Summary of the Invention

[0005] To address the shortcomings of existing technologies, this invention provides a grain storage planning method and system based on multi-source data perception. In order to solve the problems in existing technologies such as difficulty in effectively identifying the attributes of heat sources, neglecting the heterogeneity of grain pile flow resistance in ventilation control leading to low cooling efficiency, and lack of quantitative assessment of operation and maintenance economics, this invention proposes a grain storage planning method and system based on multi-source data perception.

[0006] To achieve the above objectives, the present invention provides a grain storage planning method based on multi-source data sensing, comprising the following steps:

[0007] Temperature field distribution data inside the grain pile and environmental meteorological data outside the warehouse are collected synchronously at preset time intervals, and the collected raw data sequences are spatiotemporally aligned and denoised and smoothed.

[0008] The physical space of the grain pile is mapped into a three-dimensional voxel grid based on the geometric dimensions of the warehouse. Based on the distance from the three-dimensional voxel grid to the physical boundary and the thermal penetration depth, the voxel set is divided into a thermally inert region and an unsteady thermally active region. The amount of environmental heat intrusion acting on the unsteady thermally active region is calculated using the environmental meteorological data.

[0009] Under the condition of applying step ventilation excitation, the response characteristics of the temperature field distribution data in the time dimension are tracked to determine the time when the cold air wavefront arrives at each voxel, and the equivalent permeability of each voxel is calculated based on the time.

[0010] The equivalent permeability is used to reconstruct the pore velocity vector field inside the grain pile, and the environmental heat intrusion amount is combined to construct an unsteady flow-heat coupling energy equation. By solving the unsteady flow-heat coupling energy equation in reverse, the biothermal production rate generated solely by grain respiration and insect and mold activity is separated.

[0011] Based on the equivalent permeability and the biological heat production rate, the thermodynamic impedance of each voxel and the theoretical energy consumption required to suppress the temperature rise are calculated. Combined with the value of grain dry matter loss, a time integral calculation is performed to generate a dynamic maintenance entropy that characterizes the generalized cost required to maintain the state of food security during the remaining storage period.

[0012] A multi-objective operations model is constructed with the goal of minimizing the total operating cost of the entire warehouse. The dynamic maintenance entropy is mapped to a preset decision logic table, and a storage planning strategy including the operating frequency of ventilation units and the opening and closing status of zone dampers is generated and distributed to the execution agency.

[0013] Preferably, the step of synchronously collecting temperature field distribution data inside the grain pile and environmental meteorological data outside the warehouse at preset time intervals includes:

[0014] The system iterates through the temperature sensor nodes on the internal temperature measuring cable of the grain pile via the communication interface, obtains the real-time temperature value and physical address code of each node, and simultaneously reads environmental meteorological monitoring data.

[0015] The missing nodes in the data sequence are identified, interpolation is performed to complete the sequence based on the data of adjacent nodes, and the sliding time window algorithm is used to smooth the completed data sequence.

[0016] Preferably, the step of calculating the amount of environmental heat intrusion acting on the non-steady-state thermally active zone using the environmental meteorological data includes:

[0017] Calculate the Euclidean distance from the center of each voxel to the nearest physical boundary, and mark the voxels whose Euclidean distance is less than the thermal penetration depth as the unsteady thermally active region;

[0018] Based on the aforementioned environmental meteorological data and the thermal radiation characteristics of the warehouse wall surface, the overall air temperature acting on the outer surface of the warehouse wall is calculated;

[0019] Using a quasi-steady-state heat transfer model, the heat flux per unit volume of each voxel entering the unsteady thermally active zone via conduction through the bin wall is calculated, and the heat flux per unit volume is taken as the amount of environmental heat intrusion.

[0020] Preferably, the calculation of the equivalent permeability of each voxel based on the time inversion includes:

[0021] Calculate the second-order difference of the temperature data sequence of each voxel with respect to time, and determine the time corresponding to the maximum point of the second-order difference after the start time of the step ventilation excitation as the wavefront arrival time.

[0022] Based on Darcy's law, the isotropic equivalent permeability of each voxel is calculated using the time difference between the arrival and start times of the wavefront, aerodynamic viscosity, grain pile porosity, and total pressure difference of the ventilation system.

[0023] Specifically, for voxels whose wavefront arrival time is not detected, the equivalent permeability of the voxel is assigned a preset blocking threshold.

[0024] Preferably, the step of separating the biothermal production rate solely from grain respiration and insect and fungal activity by inversely solving the unsteady flow-heat coupling energy equation includes:

[0025] The pore pressure field is solved based on the equivalent permeability, and the Darcy velocity vector of each voxel is calculated accordingly.

[0026] Establish unsteady-state control equations that include convective heat transfer terms, heat conduction terms, and total heat source generation rate;

[0027] The unsteady-state control equations are discretized using a difference scheme. The amount of environmental heat intrusion is taken as a known background term. The heat transferred from the external environment and the heat exchanged by ventilation and convection are deducted from the total energy change to calculate the biological heat production rate. Non-negative constraint corrections are then applied to the biological heat production rate.

[0028] Preferably, the generation of the dynamic maintenance entropy, which characterizes the generalized cost required to maintain food security over the remaining storage period, includes:

[0029] Based on the equivalent permeability and the difference between the current grain temperature and the air temperature at the inlet, the thermodynamic impedance characterizing the difficulty of ventilation and heat exchange is calculated.

[0030] Based on the biological heat production rate, the thermodynamic impedance, and the flow field distribution correction coefficient, calculate the theoretical ventilation power density required to suppress the temperature rise;

[0031] Based on the stoichiometry of respiration, calculate the rate of loss of dry matter value caused by the said biological heat production rate;

[0032] The dynamic maintenance entropy is obtained by summing the power consumption cost and dry matter loss value within the remaining storage period using a time-iteration algorithm.

[0033] Preferably, the construction of the multi-objective operations research model aimed at minimizing the total operating cost of the entire warehouse includes:

[0034] Define a vector of decision variables that includes the operating frequency of the ventilation unit, the opening and closing status of the zone dampers, and the start and stop status of the refrigeration equipment;

[0035] Establish an objective function to minimize the sum of the total system energy consumption cost, the dry matter loss value under controlled conditions, and the penalty for violating security constraints during the remaining storage cycle;

[0036] Constraints including upper limit of grain temperature, lower limit of moisture content, and equipment operating range are set. Feasible control sequences are traversed, and the sequence that minimizes the objective function is selected as the optimal control strategy.

[0037] Preferably, the step of mapping the dynamic maintenance entropy to a preset decision logic table to generate a storage planning strategy that includes the operating frequency of the ventilation unit and the opening and closing status of the zone dampers includes:

[0038] If it is determined that the thermodynamic impedance of each voxel is greater than the preset impedance threshold and the input-output ratio corresponding to the dynamic maintenance entropy exceeds the preset stop-loss coefficient, and the model predicts that the temperature drop does not meet the requirements, then the local blocking control logic is activated to generate an instruction to close the corresponding high impedance area partition damper.

[0039] If it is determined that the average dynamic maintenance entropy of the entire warehouse is at a moderate level and the second derivative of the biological heat production rate is greater than zero, then the energy efficiency priority control logic is activated, the upper limit of the target control temperature is increased, and priority rotation suggestions are generated.

[0040] If it is determined that the thermodynamic impedance of the entire chamber is uniform and the biological heat production rate is in a low-level stable state, then the low-temperature maintenance control logic is activated, and the intermittent low-volume operation mode is locked at night.

[0041] Preferably, after the storage planning strategy is issued to the execution agency, it further includes:

[0042] Within a preset observation window after the storage planning strategy is implemented, the actual cooling rate of the grain pile is monitored in real time and the deviation between the actual cooling rate and the model prediction value is calculated.

[0043] If the deviation exceeds a preset threshold, the equivalent penetration rate is adaptively corrected using the sign of the deviation, and the optimization calculation process of the multi-objective operations research model is re-executed based on the corrected equivalent penetration rate.

[0044] A grain storage planning system based on multi-source data sensing includes:

[0045] The data acquisition module is used to synchronously collect temperature field distribution data inside the grain pile and environmental meteorological data outside the warehouse at preset time intervals, and to perform noise reduction processing on the collected raw data sequence.

[0046] The spatial discretization module is used to map the physical space of the grain pile into a three-dimensional voxel grid, divide the voxel set into a thermally inert region and an unsteady thermally active region according to the thermal penetration depth, and calculate the amount of environmental heat intrusion acting on the unsteady thermally active region.

[0047] The permeability reconstruction module is used to invert and calculate the equivalent permeability of each voxel based on the time when the temperature wavefront arrives at each voxel under the applied step ventilation excitation.

[0048] The heat source inversion module is used to combine the equivalent permeability and the amount of environmental heat intrusion, and separate the biological heat production rate by solving the unsteady flow-heat coupling energy equation in reverse.

[0049] The maintenance entropy calculation module is used to calculate the thermodynamic impedance and theoretical energy consumption based on the equivalent permeability and the biological heat production rate, and generate a dynamic maintenance entropy that characterizes the generalized cost over the remaining storage period.

[0050] The decision planning module is used to solve the multi-objective operations research model based on the dynamic maintenance entropy, generate storage planning strategies, and distribute them to the execution agency.

[0051] This invention provides a method and system for grain storage planning based on multi-source data sensing. It has the following beneficial effects:

[0052] 1. This invention separates the endogenous biological heat production rate inside the grain pile from the environmental heat intrusion by constructing a boundary thermal effect decoupling model and a flow-heat coupling inversion mechanism. This enables the quantitative location of hidden insect and mold activities and can lock abnormal heat sources before the local temperature rises significantly. In contrast, existing monitoring technologies mostly rely on a single temperature threshold for alarms, which are easily affected by the environmental temperature rise caused by sunlight exposure on the warehouse walls, resulting in false alarms. Furthermore, hidden heat points inside the grain pile caused by poor airflow due to pore compaction are often not detected in time due to the lag in heat transfer, leading to delays in handling.

[0053] 2. This invention introduces a dynamic maintenance entropy assessment index, which integrates thermodynamic impedance, mechanical ventilation energy consumption, and grain dry matter loss value into a unified quantitative decision-making dimension. It can accurately identify areas where the cost of ventilation and cooling is far higher than the value preservation benefits of grain, thereby avoiding ineffective energy input in areas with extremely low heat exchange efficiency. Existing technologies usually adopt a crude whole-warehouse ventilation strategy, ignoring the differences in local airflow resistance inside the grain pile. This approach not only makes it difficult to effectively reduce the temperature in high-impedance areas, but also causes surrounding normal grain to suffer severe water loss and weight reduction due to prolonged excessive ventilation, resulting in significant economic losses.

[0054] 3. This invention utilizes a multi-objective operational model to generate differentiated control strategies, including local blocking, energy efficiency priority, and low-temperature maintenance. This enables precise intervention in different grain conditions. For high-resistance biothermal zones, the system can proactively cut off the corresponding airflow input to prevent heat from spreading through forced convection, rather than blindly increasing fan power. Existing automatic control systems mostly use simple target value start-stop logic, lacking spatially targeted adjustment capabilities. This can easily cause unnecessary fluctuations in the overall grain temperature when dealing with local anomalies, and it also cannot flexibly adjust operating strategies according to electricity price fluctuations to reduce operating costs. Attached Figure Description

[0055] Figure 1 This is a schematic diagram of the system structure according to an embodiment of the present invention;

[0056] Figure 2 This is a schematic diagram of the method flow according to an embodiment of the present invention;

[0057] Figure 3 This is a schematic diagram of the multiphysics coupling inversion logic flow according to an embodiment of the present invention;

[0058] Figure 4 This is a schematic diagram of the strategy decision-making logic based on maintenance entropy in an embodiment of the present invention.

[0059] in:

[0060] 100. Data sensing unit; 110. Temperature detection array; 120. Meteorological monitoring station; 200. Flow field excitation unit; 210. Ventilation actuator; 220. Programmable logic controller; 300. Computing control center; 310. Edge computing server; 311. Data acquisition module; 312. Spatial discretization module; 313. Permeability reconstruction module; 314. Heat source inversion module; 315. Maintenance entropy calculation module; 316. Decision planning module; 320. Communication interface. Detailed Implementation

[0061] The technical solutions in 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.

[0062] See attached document Figure 1 , Figure 1 This is a schematic diagram of a grain storage planning system based on multi-source data sensing according to an embodiment of the present invention. The present invention provides a grain storage planning system based on multi-source data sensing, which is applied to industrial grain storage environments such as flat warehouses or shallow round silos equipped with mechanical ventilation facilities and electronic grain condition detection facilities.

[0063] The grain storage planning system based on multi-source data sensing mainly includes a data sensing unit 100, a flow field excitation unit 200, and a computing control center 300. The data sensing unit 100, the flow field excitation unit 200, and the computing control center 300 are connected through an industrial fieldbus and a local area network, forming a closed-loop control system from physical parameter acquisition and physical field reconstruction to logistics decision generation.

[0064] The data sensing unit 100 is configured to acquire data on the internal physical state of the grain pile and meteorological data of the external environment. The data sensing unit 100 includes a temperature detection array 110 and a meteorological monitoring station 120. The temperature detection array 110 consists of multiple vertically distributed digital temperature measuring cables, each encapsulating temperature sensors arranged at predetermined intervals, thereby constructing a three-dimensional monitoring grid inside the grain pile. The meteorological monitoring station 120 is deployed outside the warehouse and is configured to collect data on ambient dry-bulb temperature, relative humidity, and solar radiation intensity.

[0065] The flow field excitation unit 200 is configured to apply airflow pressure boundary conditions to the grain pile and serves as a signal generation device for fluid characteristic detection. The flow field excitation unit 200 includes a ventilation actuator 210 and a programmable logic controller 220. The ventilation actuator 210 includes an axial flow fan or centrifugal fan connected to a ventilation duct at the bottom of the grain pile. The programmable logic controller 220 is electrically connected to the ventilation actuator 210 and is configured to control the start / stop status and operating power of the fan, and to provide feedback on the fan's operating status signal.

[0066] The computing control center 300 is configured to process sensor data, execute physical inversion algorithms, and generate operations planning instructions. The computing control center 300 includes an edge computing server 310 and a communication interface 320. The edge computing server 310 internally deploys multiple functional modules for executing the grain storage planning method provided by this invention. The communication interface 320 is configured to establish communication connections between the edge computing server 310, the programmable logic controller 220, and the external warehouse management system 400.

[0067] like Figure 1 As shown, the edge computing server 310 is logically divided into the following functional modules: data acquisition module 311, spatial discretization module 312, penetration rate reconstruction module 313, heat source inversion module 314, maintenance entropy calculation module 315, and decision planning module 316.

[0068] The data acquisition module 311 is communicatively connected to the data sensing unit 100 and is configured to synchronously read temperature data from the temperature detection array 110 and meteorological data from the meteorological monitoring station 120 at set time intervals, and perform noise reduction preprocessing on the raw data. The spatial discretization module 312 is configured to establish a three-dimensional voxel model of the grain pile. Based on the geometric dimensions of the warehouse and the spatial coordinates of the temperature measuring cable, the spatial discretization module 312 divides the grain pile into a regular voxel grid and calculates the Euclidean distance of each voxel to the nearest physical boundary. Based on a set heat penetration depth threshold, the spatial discretization module 312 marks the voxel set as an unsteady thermally active region significantly affected by the environment and a relatively adiabatic thermally inert region, and calculates the environmental heat flux acting on the unsteady thermally active region based on meteorological data.

[0069] The permeability reconstruction module 313 is configured to calculate the equivalent permeability distribution inside the grain pile. The permeability reconstruction module 313 receives time-series temperature data transmitted by the data acquisition module 311 and identifies the step cold or heat source signal generated by the flow field excitation unit 200. The permeability reconstruction module 313 calculates the second derivative of the temperature change of each voxel over time, determines the arrival time of the temperature wavefront, and solves for the equivalent permeability of each voxel based on Darcy's law and wavefront propagation time, thereby constructing a three-dimensional scalar field reflecting the airflow resistance characteristics inside the grain pile. The heat source inversion module 314 is configured to calculate the biothermal production rate inside the grain pile. Based on the equivalent permeability output by the permeability reconstruction module 313 and the environmental heat flux output by the spatial discretization module 312, the heat source inversion module 314 constructs an inverse problem model of the unsteady-state heat conduction and convection coupling equations. The heat source inversion module 314 solves in reverse to eliminate the temperature change components caused by environmental heat transfer and ventilation convection, and separates the endogenous biological heat production rate generated only by grain respiration and insect and mold activity.

[0070] The maintenance entropy calculation module 315 is configured to quantify the generalized cost of grain storage. Based on the equivalent permeability and the temperature difference between the current grain temperature and the environment, the maintenance entropy calculation module 315 calculates the thermodynamic impedance of each voxel, characterizing the energy cost of ventilation and cooling in that area. Combining the biological heat production rate, thermodynamic impedance, unit energy consumption price, and quality loss discount factor, the maintenance entropy calculation module 315 calculates the dynamic maintenance entropy required to maintain the grain in the warehouse until the predetermined delivery time through integration. The decision planning module 316 is configured to generate specific storage and logistics strategies. The decision planning module 316 receives the cost matrix output by the maintenance entropy calculation module 315 and order data fed back by the warehouse management system 400. Based on preset circuit breaker thresholds and cost optimization objectives, the decision planning module 316 generates mutually exclusive strategies, including dead zone circuit breaker, turnover priority, or strategic reserves, and converts these strategies into fan control commands and inbound / outbound work orders through the communication interface 320.

[0071] See attached document Figure 2 , Figure 2 This is a general flowchart of a grain storage planning method according to an embodiment of the present invention. The present invention provides a grain storage planning method based on multi-source data sensing, comprising the following steps:

[0072] S100, Perform data acquisition. The data acquisition module 311 periodically acquires data on the internal temperature field distribution of the grain pile and external environmental meteorological data, while recording the operating status parameters of the ventilation equipment;

[0073] S200 performs model construction and boundary decoupling. The spatial discretization module 312 maps the physical space of the grain pile into a three-dimensional voxel mesh, distinguishes between voxels in the thermally inert region and voxels in the unsteady thermally active region, and uses a comprehensive air temperature model to calculate the amount of environmental heat intrusion received by voxels in the unsteady thermally active region;

[0074] S300 performs flow field characteristic detection. When the preset temperature difference condition is met, the flow field excitation unit 200 applies step ventilation excitation, and the permeability reconstruction module 313 tracks the wavefront characteristics in the temperature response curve, calculates the equivalent permeability field of the entire chamber, and identifies high flow resistance areas;

[0075] S400 performs biological heat source location. The heat source inversion module 314 substitutes the equivalent permeability field and environmental heat flux into the inverse operation process of the energy conservation equation to quantify the biological heat production rate of each voxel and realize the location of endogenous heat generation;

[0076] S500, quantification of execution costs. The maintenance entropy calculation module 315 calculates thermodynamic impedance based on physical field parameters and combines them with economic parameters to calculate dynamic maintenance entropy, assessing the overall cost of continued storage;

[0077] S600 performs multi-objective operational decision-making. The decision planning module 316 matches and optimizes the dynamically maintained entropy with external logistics orders, determines the optimal handling plan under the current state, and issues the generated instructions to the corresponding execution agencies.

[0078] In step S100, the data acquisition module 311 performs synchronous acquisition and preprocessing of grain condition and meteorological data. This step is used to generate a raw monitoring data sequence with a unified spatiotemporal reference. Specifically, it includes the following sub-steps:

[0079] S110, the data acquisition module 311 establishes a communication link with the multi-dimensional sensing unit through the communication interface 320 in a master-slave polling mode. For the temperature detection array 110, the data acquisition module 311 issues a read command according to a preset polling period (e.g., 10 to 60 minutes), traversing all temperature sensor nodes on the temperature measuring cables inside the grain pile to obtain the real-time temperature value and corresponding physical address code of each node. For the hardware implementation of the temperature sensor, this embodiment uses a digital single-bus temperature sensor (such as DS18B20 and its compatible models) or a resistance temperature detector (RTD) sensor (such as PT100) that supports multi-point parallel connection, transmitting data via RS485 industrial bus and Modbus RTU protocol. For the meteorological monitoring station 120, the data acquisition module 311 synchronously reads the current dry-bulb temperature, relative humidity, and solar radiation intensity values.

[0080] To eliminate time discrepancies caused by differences in response time between different sensors and transmission delays, the data acquisition module 311 employs timestamp alignment processing, uniformly marking all data acquired within the same polling cycle as the same discrete time point. If individual sensor data packets are lost during a certain polling cycle, the system identifies the physical address of the missing data, retrieves the data of the adjacent nodes in the vertical direction of that physical address, and calculates the complete value using linear interpolation to maintain the integrity of the data matrix.

[0081] In step S120, the data acquisition module 311 performs denoising and smoothing processing on the acquired raw time series data. Since the subsequent step S300 involves calculating the second derivative of temperature over time, high-frequency random noise in the original signal can interfere with the differential calculation results. Therefore, this invention employs a sliding time window filtering algorithm to smooth the raw temperature sequence.

[0082] Specifically, the data acquisition module 311 selects a sliding window of a preset length (e.g., 3 to 5 sampling periods), and performs a weighted summation of the raw temperature data within the window. The weighting coefficient decreases as the time distance increases, giving higher weight to the data at the current time and nearby times, thereby obtaining a smoothed temperature sequence. Temperature sequence after noise reduction This will serve as input data for subsequent physical field inversion.

[0083] S130, the data acquisition module 311 maps the one-dimensional linear sensor data into a three-dimensional discrete matrix. Based on the pre-stored warehouse engineering drawings and temperature measurement cable layout plan, the system pre-generates the physical addresses of the sensors and the three-dimensional spatial coordinates of the grain pile during the configuration phase. The mapping table. The data acquisition module 311 will process the smoothed temperature data. Fill to the corresponding spatial coordinates, construct time. The system collects three-dimensional discrete temperature field data. Simultaneously, it associates and stores the acquired fan operating status (such as on / off, operating frequency, and power feedback value) with the temperature field data at that moment, providing a status basis for subsequent identification of flow field excitation signals. Through these steps, the system completes the conversion from underlying physical signals to structured spatiotemporal data.

[0084] In step S200, the spatial discretization module 312 performs numerical modeling of the physical space of the grain pile and decoupling calculation of boundary thermal effects. Specifically, this includes the following sub-steps:

[0085] S210, the spatial discretization module 312 constructs a three-dimensional voxel model of the grain pile based on the warehouse civil engineering drawings and the physical distribution of the temperature measuring cables. The spatial discretization module 312 establishes a Cartesian coordinate system with the center or corner of the warehouse floor as the origin, defining the effective stacking area of ​​the grain pile. Divided into A regular rectangular voxel unit.

[0086] To address the non-uniform distribution of temperature measurement points, the spatial discretization module 312 maps the three-dimensional discrete temperature field data obtained in step S100 to the center of the voxel grid. For voxel nodes without sensors, the system uses Kriging Interpolation or Inverse Distance Weighted (IDW) to calculate their initial temperature estimates, thereby generating a continuous initial temperature field matrix covering the entire chamber. This temperature field matrix serves as a unified data carrier for subsequent permeability field reconstruction, heat source inversion, and impedance calculation.

[0087] S220, the spatial discretization module 312 marks the region attributes of the voxel set based on thermal conductivity characteristics, and constructs a boundary thermal effect decoupling model. The spatial discretization module 312 traverses all voxels. Calculate the Euclidean distance from the voxel center coordinates to the nearest physical boundary (including the silo walls, silo floor, and grain surface). .

[0088] Spatial Discrete Module 312 calculates thermal penetration depth This serves as the physical basis for dividing the thermally inert region into the unsteady thermally active region. Calculate according to the following formula:

[0089] ;

[0090] in, The thermal conductivity of grain, For grain density, The specific heat capacity of grain; The angular frequency of the change in ambient temperature is used in this embodiment to isolate the thermal interference caused by the diurnal temperature range. Values (Corresponds to a 24-hour cycle).

[0091] Spatial discretization module 312 generates a region indication function based on the following logic. :when When the voxel is determined to be in a non-steady-state thermally active region significantly affected by the external environment, it is marked as a non-steady-state thermally active region. , place ;when When the voxel is determined to be in the thermal inert region, it is marked as such. , place .

[0092] S230, Spatial Discrete Module 312 calculates the action on the unsteady thermally active region. Environmental heat source items of various internal elements The spatial discrete module 312, combined with data collected from the meteorological monitoring station, calculates the current time. Combined air temperature acting on the outer surface of the warehouse wall in different orientations :

[0093] ;

[0094] in, The ambient dry-bulb temperature; This represents the solar radiation intensity corresponding to the orientation. The solar radiation absorptivity of the warehouse wall surface; The convective heat transfer coefficient of the outer surface of the silo wall; The emissivity of the bin wall surface; It is the Stefan-Boltzmann constant; The effective temperature of the sky.

[0095] Based on the calculation The spatial discrete module 312 uses a quasi-steady-state one-dimensional heat transfer model to calculate the voxels that enter the unsteady thermally active region via conduction through the bin walls. Environmental thermal intrusion rate per unit volume :

[0096] ;

[0097] in, The total heat transfer coefficient of the warehouse wall (including the thermal resistance of the wall's thermal conductivity and the convective heat transfer resistance of the inner surface). This represents the effective heat transfer area of ​​the voxel in physical contact with the bin wall (for voxels inside the unsteady thermally active zone that do not directly contact the bin wall, this value is 0 or corrected according to the thermal diffusion equation). This is the current measured temperature of the voxel; This is the volume of the voxel. The values ​​are stored as an environmental heat load matrix, which is used to remove them as background noise in the subsequent heat source inversion step.

[0098] See attached document Figure 3 , Figure 3 This is a logical flowchart of the coupling inversion of permeability field and biothermal source field according to an embodiment of the present invention. In step S300, the permeability reconstruction module 313 performs flow field characteristic detection and parameter reconstruction based on the cold source step response. This step inverts the airflow resistance distribution inside the grain pile by analyzing the thermal response characteristics of the medium under unsteady temperature boundary conditions. Specifically, it includes the following sub-steps:

[0099] S310, the flow field excitation unit 200 determines the triggering timing of the step excitation based on the ambient temperature difference and applies a constant pressure boundary condition. The permeability reconstruction module 313 calculates the average temperature of the current grain pile in real time. With ambient dry bulb temperature The difference. When ,in, The system generates flow field detection commands when the preset effective temperature difference threshold is set, for example, 8°C to 10°C.

[0100] After receiving the command, the flow field excitation unit 200 drives the ventilation actuator 210 to open and locks the fan operating frequency to a fixed value (e.g., 50Hz power frequency) to establish a constant static pressure input at the ventilation opening. The system records the moment when the fan reaches its rated speed from rest and operates stably as the start time of the step signal. At this point, the boundary condition at the bottom of the grain pile is considered as a step function that abruptly changes from a static state to a forced convection state.

[0101] S320, the penetration rate reconstruction module 313 performs temperature wavefront tracking and arrival time determination on all voxels in the storage compartment. This applies to any voxel. The system extracts its in Subsequent time series temperature data Due to the thermal inertia of the grain medium, the temperature change exhibits a gradual characteristic. In order to accurately pinpoint the moment when the convective heat transfer effect arrives, this invention employs the second-order difference extremum method for discrete data.

[0102] The permeability reconstruction module 313 calculates an approximation of the second derivative of the temperature series with respect to time. For discrete time points... Calculate the second difference value using the central difference formula. :

[0103] ;

[0104] in, This is the sampling time interval. System search. exist The first obvious maximum point afterward. This maximum point corresponds to the moment when the cold air front passes through the center of this voxel, causing a sudden change in the local cooling rate due to maximum acceleration. Define the wavefront arrival time. for:

[0105] ;

[0106] If within the preset maximum observation window If no second derivative peak exceeding the preset signal-to-noise ratio threshold is detected, the system marks the voxel as a non-responding node.

[0107] S330, Permeability Reconstruction Module 313 retrieves the isotropic equivalent permeability of each voxel based on wavefront propagation time inversion. According to Darcy's law and the principles of porous media fluid dynamics, under constant pressure difference, the wavefront arrival time is inversely proportional to the media permeability.

[0108] The permeability reconstruction module 313 calculates the equivalent permeability of each voxel according to the following formula. :

[0109] ;

[0110] in, The dynamic viscosity of air is calculated, and this parameter is corrected in real time based on the current air temperature by looking up a table or using the Sutherland formula. The average porosity of the grain accumulation ranges from 0.35 to 0.45. Provides an effective total differential pressure for the ventilation system; The airflow path length is defined as the Euclidean distance from the geometric center of the voxel to the nearest inlet surface. This is a dimensionless geometric correction factor used to compensate for the difference between actual streamline curvature (tortuosity) and straight-line distance. The deviation between the two factors is related to the spatial position of the voxels and is obtained through a preset geometric correction matrix.

[0111] For nodes marked as unresponsive in step S320, the penetration rate reconstruction module 313 directly sets their penetration rate. Assigned value: blocking threshold (e.g., 1.0 × 10) -15 m 2 This area is thus numerically defined as a ventilation dead zone. By traversing all voxels, the system generates a three-dimensional scalar field reflecting the degree of airflow within the grain pile. .

[0112] As attached Figure 3 As shown, in step S400, the heat source inversion module 314 performs bio-heat source localization and quantification based on the flow-heat coupling mechanism. This step establishes a coupled physical model of unsteady heat transfer and porous media seepage inside the grain pile. Using the permeability field obtained in step S300 and the temperature field obtained in step S100 as known conditions, the distribution of internal heat sources generated solely by biological respiration is calculated by inversely solving the energy conservation equation. Specifically, it includes the following sub-steps:

[0113] S410, the heat source inversion module 314 reconstructs the pore wind speed vector field inside the grain pile based on the permeability field data. The heat source inversion module 314 calls the three-dimensional permeability scalar field generated in step S300. The equivalent permeability of each voxel in the field Substituting spatial variables into the porous media fluid control equation, the heat source inversion module 314, based on Darcy's law and the fluid continuity equation, constructs a model for pore pressure. Laplace's equation:

[0114] ;

[0115] Heat source inversion module 314 sets the boundary conditions of the solution domain: the air inlet region at the bottom of the grain pile is set as the first type of boundary condition (Dirichlet boundary), and its pressure value... The total pressure value is taken from the feedback of the ventilation actuator 210; the natural exhaust vents of the grain surface and warehouse are set as the outlet boundary, and the pressure value is taken as atmospheric pressure. Set the warehouse walls as second-type boundary conditions (Neumann boundary), with a pressure gradient of zero (i.e., The finite difference method is used to solve the above equations to obtain the pressure distribution field of the entire chamber. Furthermore, the Darcy velocity vector of each voxel is calculated based on Darcy's law. :

[0116] ;

[0117] in, For inclusion The velocity vector has three directional components; This represents the pressure gradient at that voxel; It is the aerodynamic viscosity.

[0118] S420, the heat source inversion module 314 constructs multiphysics coupled energy control equations containing source terms. The system establishes unsteady-state control equations of the following form:

[0119] ;

[0120] in, For grain bulk density; The specific heat capacity of grain; Temperature of the grain pile; For time; air density; The specific heat capacity of air at constant pressure; The thermal conductivity of grain; This represents the total heat source generation rate per unit volume. The second term on the left-hand side of the equation... The first term on the right-hand side of the equation characterizes the convective heat transfer effect that carries heat in the airflow; the second term on the right-hand side characterizes the heat conduction effect based on Fourier's law.

[0121] S430, the heat source inversion module 314 performs the reverse solution. The heat source inversion module 314 converts the environmental heat source items calculated in step S200 into... and the velocity vector calculated in step S410 Substitute into the equation. The system, through rearranging terms, transforms the biological heat production rate... Solve for it as an unknown:

[0122] ;

[0123] In order to ensure the stability of the numerical calculation, the heat source inversion module 314 uses the upwind difference scheme to discretize the convection terms during the calculation process. That is, according to the velocity vector Choose between forward or backward differencing; use the central difference scheme to discretize the diffusion term. The discrete-time derivative term is obtained using the backward difference scheme. This formula allows the system to subtract the heat transferred from the external environment and the heat exchanged through ventilation and convection from the total energy change; the remaining term represents the endogenous biothermal rate generated by grain respiration and insect and fungal activity.

[0124] S440, the heat source inversion module 314 performs physical constraint corrections on the calculation results. Considering that measurement noise and discretization errors may lead to negative calculated heat production rates, the system... Apply nonnegativity constraints: if the calculation result If the value is less than 0, it is corrected to 0. Furthermore, the system employs local median filtering. The three-dimensional matrix is ​​spatially smoothed to eliminate isolated noise points, and the final output is the distribution field of biological heat production rate.

[0125] In step S500, the maintenance entropy calculation module 315 performs a quantitative calculation of the generalized cost of grain storage. This step, based on the physical field parameters obtained in the previous steps, establishes a comprehensive evaluation model that includes active energy consumption and passive quality loss, and calculates the dynamic maintenance entropy of each voxel and the entire storage facility. Specifically, it includes the following sub-steps:

[0126] S510, Maintenance Entropy Calculation Module 315 calculates the thermodynamic impedance of each voxel. This parameter characterizes the process resistance to removing a unit of heat using mechanical ventilation at a specific spatial location. The maintenance entropy calculation module 315 calculates the equivalent permeability based on the output of step S300. and the grain temperature collected in step S100 Inlet air temperature obtained from meteorological monitoring stations Calculate the thermodynamic impedance:

[0127] ;

[0128] in, Aerodynamic viscosity; A preset minimum heat exchange temperature difference threshold (e.g., 1.0℃) is used to prevent singular denominators. When When the system determines that the area is in thermal equilibrium or unable to exchange heat, it will... Marked as high impedance saturation value.

[0129] S511, Maintenance Entropy Calculation Module 315 calculates the theoretical ventilation power density required to suppress temperature rise. The entropy calculation module 315 maintains the biological heat production rate output in step S400. and the velocity vector magnitude output in step S410 The energy consumption of each voxel can be calculated using the following formula:

[0130] ;

[0131] in, The electro-wind energy conversion efficiency of the ventilation system; Voxel volume; This is the flow field distribution correction factor. In this embodiment, With the local Darcy flow rate of this voxel They are inversely proportional, defined as ,in The average flow rate across the entire warehouse. To prevent minute amounts from being divided by zero, this coefficient is used for energy-weighted compensation of heat transfer efficiency degradation in low-flow-rate regions.

[0132] S512, Maintenance Entropy Calculation Module 315 calculates the value loss rate corresponding to dry matter loss. Based on the stoichiometry of respiration, the entropy calculation module 315 calculates the value of quality loss per unit time:

[0133] ;

[0134] in, This refers to the current market price per unit of grain. The heat-to-mass conversion coefficient is 0.06 × 10⁻⁶. -6 kg / J (calculated based on the enthalpy change of glucose oxidation reaction).

[0135] S513, the entropy calculation module 315 performs discrete iterative calculations to generate a dynamic entropy matrix. The entropy calculation module 315 obtains the predetermined outbound time. Determine the remaining storage period .

[0136] To perform the prediction calculations, the system employs the forward Euler method for time-step iteration. At each time step... Within this timeframe, the system updates the biological heat production rate for the next moment based on the Arrhenius equation. ,in This is the temperature coefficient of the reaction rate.

[0137] The maintenance entropy calculation module 315 calculates the dynamic maintenance entropy of each voxel by accumulating and summing. :

[0138] ;

[0139] in, To obtain the unit electricity price that changes over time, the system reads the time-of-use electricity price table to get the electricity price data for the corresponding time. Characterized voxels Total generalized cost over the remaining storage period.

[0140] S514, Maintenance Entropy Calculation Module 315 calculates the total dynamic maintenance entropy of the entire warehouse. And generate a decision matrix. System statistics. Exceeding the preset threshold A set of voxels. When a voxel satisfies (in When the stop-loss coefficient is set to a preset value (e.g., 1.5), it indicates that the cost of active cooling in the area is significantly higher than its value preservation benefit. The system marks the voxel as an "uneconomically controllable area" and outputs it to the decision planning module 316 as the basis for determining partial circuit breaker or emergency outbound.

[0141] See attached document Figure 4 , Figure 4 This is a flowchart of the control strategy decision-making logic based on dynamic maintenance entropy according to an embodiment of the present invention. In step S600, the decision planning module 316 performs strategy generation and execution mechanism instruction issuance based on multi-objective operations research. This step receives the dynamic maintenance entropy matrix output from step S500 and the physical field states output from S300 and S400, and solves for the control sequence that minimizes the total life cycle operating cost while satisfying the constraints of safe food storage. Specifically, it includes the following sub-steps:

[0142] S610, Decision Planning Module 316 constructs a mathematical operations research model with the objective of minimizing total operating costs. Decision Planning Module 316 defines the decision variable vector. ,in The operating frequency of the ventilation unit, For the first The opening and closing status of the zoned ventilation windows. This indicates the start / stop status of the refrigeration equipment.

[0143] The system establishes an objective function of the following form. Calculate the cumulative generalized cost over the remaining storage cycles:

[0144] ;

[0145] in, Time-of-use pricing; The total energy consumption of the system is modeled according to the similarity law of wind turbines. ,in For power coefficient, For the rated frequency, Rated power; Value of dry matter loss under controlled conditions; This is a penalty coefficient for safety constraints; To predict the number of voxels whose states exceed safety constraints.

[0146] Decision Planning Module 316 sets constraints:

[0147] 1. State Constraints: Predicting Temperature with Arbitrary Voxels And the water content of any voxel .

[0148] 2. Control Constraints: Furthermore, the rate of change of the control variable is limited by the acceleration and deceleration time of the equipment.

[0149] 3. Environmental constraints: Only when the outside temperature is... And external humidity When the moisture level is within the equilibrium range, natural ventilation mode is permitted.

[0150] To address this optimization problem, the decision planning module 316 employs the Discrete Grid Search method. The system divides the future prediction time domain (e.g., 24 hours) into... Each control step will control the frequency. The system is discretized into a finite number of speed levels (e.g., [0, 25, 30, ..., 50] Hz). The system iterates through all feasible control sequence combinations and calculates the corresponding objective function values. Select to make The smallest sequence is the optimal control strategy. .

[0151] S620, the decision planning module 316 matches and activates the corresponding maintenance control logic based on the optimization results and the maintenance entropy characteristics of each region. The system will calculate the maintenance entropy of all voxels in the warehouse. Distribution and thermodynamic impedance The distribution is mapped to a preset logical table, and the following judgment is performed:

[0152] 1. Local blocking control logic (corresponding to the original dead zone circuit breaker):

[0153] When a local region voxel satisfies and If the preset stop-loss threshold is exceeded, and the model prediction shows that the temperature drop in the area will be less than 1.0℃ in the next 48 hours under maximum ventilation power, the system will activate local blocking control.

[0154] The system generates instructions to close the ventilation damper corresponding to the high-impedance area and cut off energy input; at the same time, it generates a local high temperature warning and an emergency outbound suggestion signal and sends it to the host computer.

[0155] 2. Energy efficiency priority control logic (corresponding to the original turnover priority):

[0156] When the total position maintains the average entropy It is at a moderate level, but the biological heat production rate is... When the second derivative is greater than 0 (accelerating heat generation), the system activates energy efficiency priority control.

[0157] The system raises the upper limit of the target control temperature (e.g., from 20°C to 25°C) to reduce cooling energy consumption and recommends rotation before grains are downgraded.

[0158] 3. Low-temperature maintenance control logic (corresponding to the original strategic reserve):

[0159] When the thermodynamic impedance of the entire chamber Uniform and low, and biological heat production rate When the system is in a low-level stable state, it activates cryogenic holding control.

[0160] The system locks into an intermittent low-volume mode, driving the fan at a low frequency (such as 25Hz-30Hz) only during off-peak hours when electricity prices are low at night, to maintain the low-temperature cooling state of the grain pile.

[0161] In step S630, the decision planning module 316 converts the generated control strategy into low-level hardware execution instructions and issues them. The decision planning module 316 has a built-in device driver library that processes the control variables obtained in step S610. It is interpreted as an industrial bus signal.

[0162] For ventilation control, the system writes the frequency setpoint to the inverter via the Modbus-RTU protocol; for ventilation window control, the system sends switching commands to the PLC. Simultaneously, the system activates a closed-loop feedback mechanism. This feedback is observed in the monitoring window after the control command is executed. Within this range, the deviation between the actual cooling rate and the model prediction exceeds a threshold. Decision planning module 316 triggers adaptive parameter correction:

[0163] ;

[0164] in, To correct the step size factor (e.g., 0.1). This is a sign function. The system re-executes the optimization calculation in step S610 based on the corrected permeability parameters to achieve online calibration of the control model.

[0165] To further clarify the collaborative working process of the technical solution of this invention, a specific working scenario example will be used below. This example demonstrates how, when hidden localized heating occurs inside a grain pile, the system avoids the energy waste and delays caused by traditional one-size-fits-all ventilation through multiphysics field inversion and multi-objective operational logic.

[0166] Assume the application scenario is a single-story warehouse with an installed capacity of 5,000 tons, storing corn, currently in the spring when temperatures are rising. The external ambient temperature is a daytime high of 25°C and a nighttime low of 10°C. The warehouse management system is set to a safety alarm temperature threshold of 30°C.

[0167] During the monitoring phase, the system executes steps S100 to S200. Data acquisition module 311 collects temperature data for the entire warehouse at 30-minute intervals. At that moment, the system detected a location in the lower layer of the grain pile with coordinates [missing information]. local area The temperature rises to 32°C. To determine the nature of this temperature rise, the spatial discretization module 312 calculates the distance from the region to the nearest boundary. Calculations show that... (Thermal penetration depth, approximately 0.5m in this example), the system determines that this area belongs to the thermal inert region, and its region indication function... Based on this, the system ruled out the possibility of environmental heat intrusion (such as heat from the warehouse walls) causing the temperature rise, and preliminarily identified an abnormal internal heat source in the area.

[0168] Subsequently, the system triggers step S300 to detect flow field characteristics. The flow field excitation unit 200 utilizes the lower nighttime temperatures (e.g., 2:00 AM) to start the fan and maintain a fixed frequency of 50Hz, applying a step-cold source excitation. The permeability reconstruction module 313 tracks the arrival of the cold air wavefront in this anomalous region. The monitoring data shows that, compared to other normal areas on the same floor, wavefront arrival time There is a significant lag, exceeding 2 hours. The system is then used to calculate the equivalent permeability of the area using the permeability inversion formula. Approximately 2.5 × 10 -11 m 2 This is significantly lower than the average level of the entire portfolio (1.0 × 10). -9 m 2 Based on this, the system determined that the area had severe porosity caking or impurity accumulation, and was therefore a high-resistivity region.

[0169] Based on the above physical parameters, the system executes step S400 to quantify the heat source. The heat source inversion module 314 constructs the flow-heat coupling equation. Although the ventilation system is running at this time, due to... Too low, Darcy velocity in this region Minimal, convective heat transfer term The values ​​were weak. However, temperature monitoring showed that the temperature in the area did not decrease with ventilation; instead, it continued to rise at a rate of 0.2℃ / h. Solving the energy equation in reverse yielded the biothermal production rate of the area. Up to 150W / m 3 This temperature is far above the normal respiratory fever range, indicating that there is intense fungal activity or an outbreak of pests at that location.

[0170] Next, the system executes step S500 to perform a maintenance entropy assessment. The maintenance entropy calculation module 315 calculates the thermodynamic impedance of the region. .because Extremely small and with low local temperature difference utilization. The value is approaching the saturation threshold. Substituting this into the energy consumption calculation formula, in order to suppress this 150W / m... 3 The heat production rate and the required theoretical ventilation power density It will be an astronomical figure. System calculations and comparisons revealed that the cost of the electricity invested... The value of the dry matter loss of this portion of grain has been reached. More than 5 times (i.e.) The calculation results show that attempting to blow through the hardened and heated area by simply increasing the fan power is extremely unreasonable from an economic standpoint and can easily cause surrounding normal grains to lose moisture and break due to excessive ventilation.

[0171] Finally, the system executes step S600 to generate a decision command. Based on the above calculation results, the decision planning module 316 identifies that the area meets the triggering conditions for local blocking control logic. Instead of generating a full-power ventilation command, the system generates a composite control strategy:

[0172] 1. Hardware execution level: Instructions are sent to the frequency converter to reduce the fan frequency to 25Hz, maintaining only basic low flow rate and uniform temperature throughout the chamber, avoiding unnecessary work on high-impedance areas; simultaneously, the area directly corresponding to this zone is shut down. Local ventilation dampers (if zone control is available) prevent heat from spreading to the surroundings due to forced convection.

[0173] 2. Management Interaction Level: Send local high-risk alerts to the warehouse management terminal, clearly marking the coordinates of abnormal points. The fault type is high-impedance bio-heating, and it is recommended to implement physical intervention measures such as mechanical turning or fixed-point chemical fumigation, rather than continuing to rely on ventilation and cooling.

[0174] Through this working process, when faced with complex grain condition anomalies, the present invention can see beyond the apparent high temperature to the essential high resistance and high heat generation, and make the optimal decision that meets both economic and safety requirements based on the quantitative assessment of maintenance entropy, thus achieving a technological leap from passive cooling to proactive and refined management.

[0175] Although embodiments of the invention have been shown and described, it will be understood by those skilled in the art that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the appended claims and their equivalents.

Claims

1. A grain storage planning method based on multi-source data sensing, characterized in that, Includes the following steps: Temperature field distribution data inside the grain pile and environmental meteorological data outside the warehouse are collected synchronously at preset time intervals, and the collected raw data sequences are spatiotemporally aligned and denoised and smoothed. The physical space of the grain pile is mapped into a three-dimensional voxel grid based on the geometric dimensions of the warehouse. Based on the distance from the three-dimensional voxel grid to the physical boundary and the thermal penetration depth, the voxel set is divided into a thermally inert region and an unsteady thermally active region. The amount of environmental heat intrusion acting on the unsteady thermally active region is calculated using the environmental meteorological data. Under the condition of applying step ventilation excitation, the response characteristics of the temperature field distribution data in the time dimension are tracked to determine the time when the cold air wavefront arrives at each voxel, and the equivalent permeability of each voxel is calculated based on the time. The equivalent permeability is used to reconstruct the pore velocity vector field inside the grain pile, and the environmental heat intrusion amount is combined to construct an unsteady flow-heat coupling energy equation. By solving the unsteady flow-heat coupling energy equation in reverse, the biothermal production rate generated solely by grain respiration and insect and mold activity is separated. Based on the equivalent permeability and the biological heat production rate, the thermodynamic impedance of each voxel and the theoretical energy consumption required to suppress the temperature rise are calculated. Combined with the value of grain dry matter loss, a time integral calculation is performed to generate a dynamic maintenance entropy that characterizes the generalized cost required to maintain the state of food security during the remaining storage period. A multi-objective operations model is constructed with the goal of minimizing the total operating cost of the entire warehouse. The dynamic maintenance entropy is mapped to a preset decision logic table, and a storage planning strategy including the operating frequency of ventilation units and the opening and closing status of zone dampers is generated and distributed to the execution agency.

2. The grain storage planning method based on multi-source data sensing according to claim 1, characterized in that, The synchronous collection of temperature field distribution data inside the grain pile and environmental meteorological data outside the warehouse at preset time intervals includes: The system iterates through the temperature sensor nodes on the internal temperature measuring cable of the grain pile via the communication interface, obtains the real-time temperature value and physical address code of each node, and simultaneously reads environmental meteorological monitoring data. The missing nodes in the data sequence are identified, interpolation is performed to complete the sequence based on the data of adjacent nodes, and the sliding time window algorithm is used to smooth the completed data sequence.

3. The grain storage planning method based on multi-source data sensing according to claim 1, characterized in that, The calculation of the environmental heat intrusion amount acting on the non-steady-state thermally active zone using the environmental meteorological data includes: Calculate the Euclidean distance from the center of each voxel to the nearest physical boundary, and mark the voxels whose Euclidean distance is less than the thermal penetration depth as the unsteady thermally active region; Based on the aforementioned environmental meteorological data and the thermal radiation characteristics of the warehouse wall surface, the overall air temperature acting on the outer surface of the warehouse wall is calculated; Using a quasi-steady-state heat transfer model, the heat flux per unit volume of each voxel entering the unsteady thermally active zone via conduction through the bin wall is calculated, and the heat flux per unit volume is taken as the amount of environmental heat intrusion.

4. The grain storage planning method based on multi-source data sensing according to claim 1, characterized in that, The calculation of the equivalent permeability of each voxel based on the time inversion includes: Calculate the second-order difference of the temperature data sequence of each voxel with respect to time, and determine the time corresponding to the maximum point of the second-order difference after the start time of the step ventilation excitation as the wavefront arrival time. Based on Darcy's law, the isotropic equivalent permeability of each voxel is calculated using the time difference between the arrival and start times of the wavefront, aerodynamic viscosity, grain pile porosity, and total pressure difference of the ventilation system. Specifically, for voxels whose wavefront arrival time is not detected, the equivalent permeability of the voxel is assigned a preset blocking threshold.

5. A grain storage planning method based on multi-source data sensing according to claim 1, characterized in that, The process involves solving the unsteady flow-heat coupling energy equation in reverse to separate the biothermal production rate generated solely by grain respiration and insect and fungal activity, including: The pore pressure field is solved based on the equivalent permeability, and the Darcy velocity vector of each voxel is calculated accordingly. Establish unsteady-state control equations that include convective heat transfer terms, heat conduction terms, and total heat source generation rate; The unsteady-state control equations are discretized using a difference scheme. The amount of environmental heat intrusion is taken as a known background term. The heat transferred from the external environment and the heat exchanged by ventilation and convection are deducted from the total energy change to calculate the biological heat production rate. Non-negative constraint corrections are then applied to the biological heat production rate.

6. The grain storage planning method based on multi-source data sensing according to claim 1, characterized in that, The generation of the dynamic maintenance entropy, which characterizes the generalized cost required to maintain food security over the remaining storage period, includes: Based on the equivalent permeability and the difference between the current grain temperature and the air temperature at the inlet, the thermodynamic impedance characterizing the difficulty of ventilation and heat exchange is calculated. Based on the biological heat production rate, the thermodynamic impedance, and the flow field distribution correction coefficient, calculate the theoretical ventilation power density required to suppress the temperature rise; Based on the stoichiometry of respiration, calculate the rate of loss of dry matter value caused by the said biological heat production rate; The dynamic maintenance entropy is obtained by summing the power consumption cost and dry matter loss value within the remaining storage period using a time-iteration algorithm.

7. A grain storage planning method based on multi-source data sensing according to claim 1, characterized in that, The construction of the multi-objective operations research model aimed at minimizing the total operating cost of the entire warehouse includes: Define a vector of decision variables that includes the operating frequency of the ventilation unit, the opening and closing status of the zone dampers, and the start and stop status of the refrigeration equipment; Establish an objective function to minimize the sum of the total system energy consumption cost, the dry matter loss value under controlled conditions, and the penalty for violating security constraints during the remaining storage cycle; Constraints including upper limit of grain temperature, lower limit of moisture content, and equipment operating range are set. Feasible control sequences are traversed, and the sequence that minimizes the objective function is selected as the optimal control strategy.

8. A grain storage planning method based on multi-source data sensing according to claim 7, characterized in that, The step of mapping the dynamic maintenance entropy to a preset decision logic table to generate a storage planning strategy that includes the operating frequency of the ventilation unit and the opening and closing status of the zone dampers includes: If it is determined that the thermodynamic impedance of each voxel is greater than the preset impedance threshold and the input-output ratio corresponding to the dynamic maintenance entropy exceeds the preset stop-loss coefficient, and the model predicts that the temperature drop does not meet the requirements, then the local blocking control logic is activated to generate an instruction to close the corresponding high impedance area partition damper. If it is determined that the average dynamic maintenance entropy of the entire warehouse is at a moderate level and the second derivative of the biological heat production rate is greater than zero, then the energy efficiency priority control logic is activated, the upper limit of the target control temperature is increased, and priority rotation suggestions are generated. If it is determined that the thermodynamic impedance of the entire chamber is uniform and the biological heat production rate is in a low-level stable state, then the low-temperature maintenance control logic is activated, and the intermittent low-volume operation mode is locked at night.

9. A grain storage planning method based on multi-source data sensing according to claim 7, characterized in that, After the storage planning strategy is issued to the execution agency, it also includes: Within a preset observation window after the storage planning strategy is implemented, the actual cooling rate of the grain pile is monitored in real time and the deviation between the actual cooling rate and the model prediction value is calculated. If the deviation exceeds a preset threshold, the equivalent penetration rate is adaptively corrected using the sign of the deviation, and the optimization calculation process of the multi-objective operation research model is re-executed based on the corrected equivalent penetration rate.

10. A grain storage planning system based on multi-source data sensing, applied to the grain storage planning method based on multi-source data sensing as described in any one of claims 1-9, characterized in that, include: The data acquisition module is used to synchronously collect temperature field distribution data inside the grain pile and environmental meteorological data outside the warehouse at preset time intervals, and to perform noise reduction processing on the collected raw data sequence. The spatial discretization module is used to map the physical space of the grain pile into a three-dimensional voxel grid, divide the voxel set into a thermally inert region and an unsteady thermally active region according to the thermal penetration depth, and calculate the amount of environmental heat intrusion acting on the unsteady thermally active region. The permeability reconstruction module is used to invert and calculate the equivalent permeability of each voxel based on the time when the temperature wavefront arrives at each voxel under the applied step ventilation excitation. The heat source inversion module is used to combine the equivalent permeability and the amount of environmental heat intrusion, and separate the biological heat production rate by solving the unsteady flow-heat coupling energy equation in reverse. The maintenance entropy calculation module is used to calculate the thermodynamic impedance and theoretical energy consumption based on the equivalent permeability and the biological heat production rate, and generate a dynamic maintenance entropy that characterizes the generalized cost over the remaining storage period. The decision planning module is used to solve the multi-objective operations research model based on the dynamic maintenance entropy, generate storage planning strategies, and distribute them to the execution agency.

Images