CN107330210B - Time domain analysis method of nonlinear thermoelectric power generation system with parameter changing along with temperature - Google Patents

Time domain analysis method of nonlinear thermoelectric power generation system with parameter changing along with temperature Download PDF

Info

Publication number
CN107330210B
CN107330210B CN201710550251.0A CN201710550251A CN107330210B CN 107330210 B CN107330210 B CN 107330210B CN 201710550251 A CN201710550251 A CN 201710550251A CN 107330210 B CN107330210 B CN 107330210B
Authority
CN
China
Prior art keywords
region
temperature
boundary
iii
generation system
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Expired - Fee Related
Application number
CN201710550251.0A
Other languages
Chinese (zh)
Other versions
CN107330210A (en
Inventor
李秋菊
刘楷安
张洋
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
North China University of Water Resources and Electric Power
Original Assignee
North China University of Water Resources and Electric Power
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by North China University of Water Resources and Electric Power filed Critical North China University of Water Resources and Electric Power
Priority to CN201710550251.0A priority Critical patent/CN107330210B/en
Publication of CN107330210A publication Critical patent/CN107330210A/en
Application granted granted Critical
Publication of CN107330210B publication Critical patent/CN107330210B/en
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/08—Thermal analysis or thermal optimisation

Abstract

The invention provides a time domain analysis method of a nonlinear thermoelectric power generation system with parameters changing along with temperature, which comprises the steps of firstly establishing a thermoelectric partial differential equation model and boundary conditions based on basic physical characteristics of the nonlinear thermoelectric power generation system with parameters changing along with temperature, modeling a Seebeck coefficient in the model as a function changing along with temperature, secondly carrying out more specific analysis on a system hot circuit and a circuit, carrying out time and space discretization treatment, establishing an algebraic equation of physical quantities of nodes and boundary nodes in a region from the angle of energy conservation, starting from a time initial value, carrying out iterative solution, and finally obtaining numerical solutions of the temperature and electric field intensity of all parts in the nonlinear thermoelectric power generation system with parameters changing along with temperature. Through the analysis of the numerical solution of the internal temperature and the electric field intensity of the thermoelectric power generation system, the thermoelectric power generation system can be further analyzed and designed in a more specific mode, and the result is more accurate.

Description

Time domain analysis method of nonlinear thermoelectric power generation system with parameter changing along with temperature
Technical Field
The invention relates to the technical field of semiconductor thermoelectric power generation, in particular to a time domain analysis method of a nonlinear thermoelectric power generation system with parameters changing along with temperature.
Background
In 1821, Seebeck, German scientist, found that in a closed loop of two different metals, when there is a temperature difference between the two joints, the loop will generate a current, which is called Seebeck effect. The thermoelectric power generation utilizes the Seebeck effect, and a certain temperature difference is maintained at two ends of a thermoelectric material, so that a certain voltage and electric power output are generated. Research shows that the thermoelectric figure of merit of semiconductor materials is large, and the thermoelectric materials applied in temperature difference systems by people at present are all semiconductor materials, so the thermoelectric materials are also called semiconductor temperature difference power generation.
With the knowledge of people on the energy crisis, the thermoelectric power generation technology can utilize a large amount of temperature difference and industrial waste heat existing in nature, and has good comprehensive social and economic benefits. Meanwhile, with the interest of people in space exploration, the development of medical physics and the application of large-scale wireless sensors, a power supply system which can supply energy by itself and does not need to be attended to needs to be developed, and obviously, thermoelectric power generation is very suitable for the applications. The thermoelectric power generation is used as an all-solid-state energy conversion mode, has the advantages of no medium leakage, no noise, reliable performance, less maintenance and the like, is a green and environment-friendly energy source, and has more and more obvious application value in the aspects of micro energy, low-grade energy and waste energy utilization. Therefore, the realization of the industrialization of the thermoelectric power generation technology and the application thereof as soon as possible has important practical significance. The analysis and design of the thermoelectric power generation system depend on the actual working condition, so that the reasonable thermoelectric power generation system model aiming at different actual conditions is the basis for more specifically analyzing and designing various parameters of the thermoelectric power system.
Chinese patent document CN201410189306.6 discloses a calculation method of electromotive force of thermoelectric power generation system, in which a mathematical model for analyzing and calculating thermoelectric power system is used to approximate the actual thermoelectric power generation system as follows:
considering that the physical characteristic parameter in the thermoelectric system is a constant value, namely the constant value;
the heat flow or temperature at the boundary is considered to be time-invariant, i.e. constant,
therefore, in the CN201410189306.6 calculation method of electromotive force of the thermoelectric system, the boundary condition at the boundary of different materials is simply approximated to fluid balance at the connection, and is only suitable for the analysis of the steady-state linear thermoelectric system when the thermoelectric generation system is in a normal state and the boundary can be simply approximated to a balanced state, but is not suitable for the analysis of the change of material parameters in the thermoelectric generation system, that is, it cannot be used for analyzing the transient process of the dynamic and nonlinear thermoelectric generation system.
In an actual thermoelectric system, there may be three main transmission coefficients (seebeck coefficient, electrical conductivity and thermal conductivity) which have serious temperature variation in low and medium temperature ranges, and particularly, even if the seebeck coefficient is small, the coefficient value can be greatly changed, so that the three transmission coefficients, particularly the seebeck coefficient, should be modeled as a function of the temperature variation.
Disclosure of Invention
The invention aims to provide a time domain analysis method of a nonlinear temperature difference power generation system with parameters changing along with temperature, the numerical solution of the internal temperature and the electric field intensity of dynamic temperature difference power generation is obtained by the method, the more specific analysis and design can be further carried out on the hot-end heat flow time-varying dynamic temperature difference power generation system, and the result is more accurate.
The technical scheme adopted by the invention is as follows:
a time domain analysis method of a nonlinear temperature difference power generation system with parameters changing along with temperature specifically comprises the following steps:
step one, deriving a mathematical model of a nonlinear thermoelectric generation system thermal circuit and a circuit, wherein the mathematical model is formed by d pairs of semiconductor thermocouples, and parameters of the nonlinear thermoelectric generation system thermal circuit and the circuit change along with temperature according to basic physical characteristics of materials of the thermoelectric generation system and basic laws of thermodynamics and electrics:
the mathematical model of each region in the hot circuit is:
in the region I, the first and second electrodes are arranged in a straight line,
Figure BDA0001344238380000021
in the region (II), the first and second regions,
Figure BDA0001344238380000022
in the region of the region III, the first,
Figure BDA0001344238380000023
the mathematical model of each region of the circuit is:
in the region (II), the first and second regions,
Figure BDA0001344238380000024
Figure BDA0001344238380000025
in the region IV, the light-emitting element,
Figure BDA0001344238380000026
wherein, thermoelectric generation system heat route divides into three regions: region I, region II and region III, region I representingA hot end radiator, wherein the area II represents a pair of semiconductor temperature difference thermocouples, and the area III represents a cold end radiator; the circuit in the thermoelectric power generation system is divided into two areas: region II and region IV; region II represents d versus semiconductor thermoelectric thermocouple, and region IV represents load resistance RL；
Where ρ isIAnd ρIIIDensity of substances in zones I and III, respectively, CvIAnd CvIIIConstant specific heat, k, of the materials in the regions I and III respectivelyIAnd kIIIThermal conductivity of the materials of regions I and III, respectively. T is the temperature of the molten metal,
Figure BDA0001344238380000027
for temperature gradient, E is the electric field strength, pIIRespectively the density of the substance in region II, CvIIIs the constant specific heat of the substance in the area II, epsilonIIIs the dielectric constant of the material in region II, and J is the current density in the circuit, current J0Current density through the circuit at time t-0, αII(T)，σII，kIITotal equivalent Seebeck coefficient, total equivalent electrical conductivity, total equivalent thermal conductivity of the region II material αII(T)＝d*(αP(T)-αN(T))，αP(T)And αN(T)Seebeck coefficients for P-type and N-type semiconductor materials in thermoelectric thermocouples, α, respectivelyP(T)And αN(T)Is a function of temperature, therefore αII(T)Also a function of temperature;
Figure BDA0001344238380000028
kII＝m*(kP+kN) (ii) a Wherein σP，kPThe electric conductivity and the thermal conductivity, sigma, of the P-type semiconductor material in the thermoelectric coupleN，kNThe electric conductivity and the thermal conductivity of the N-type semiconductor material in the thermoelectric couple are respectively; sigmaIVIs the conductivity of the load resistance. Determining boundary conditions of a thermal circuit and a circuit of the nonlinear temperature difference power generation system with parameters changing along with the temperature according to actual working conditions;
the boundary conditions in the hot path are:
boundary A, according to the form of energy provided by heat source, the hot end is equivalent to heat flow density of q0；
Boundary B, the density of the heat flow passing through the boundary is continuous, and the temperature is continuous;
a boundary C, wherein the heat flow density passing through the boundary is continuous, and the temperature is continuous;
boundary D: according to the heat dissipation form of the cold source end, the cold end is equivalent to the temperature of Tl0；
The boundary conditions in the circuit are:
a boundary M, wherein the current density passing through the boundary is continuous, and the electric field intensity is continuous;
the current density passing through the boundary is continuous, and the electric field intensity is continuous;
the boundary of the area I and the heat source in the heat path is a boundary A, the boundary of the area I and the area II is a boundary B, the boundary of the area II and the area III is a boundary C, and the boundary of the area III and the cold source is a boundary D; the boundary of a region II at a heat source end and a region IV in the circuit is a boundary M, wherein the boundary of the region II at a cold source end and the region IV is a boundary N;
performing time discretization treatment on the areas I, II, III and IV to determine time nodes in the areas; performing one-dimensional space discretization processing on the region I, the region II and the region III to determine space nodes of the regions;
step four, determining the temperature of each space node and the iteration initial value of the electric field intensity at the moment when t is 0;
discretizing the mathematical model, and establishing an algebraic equation of the temperature and the electric field intensity of the nodes and boundary nodes in the region;
and step six, starting from the iteration initial values of the space nodes of the temperature and the electric field intensity, carrying out iteration solution according to the algebraic equations of the node temperature and the electric field intensity established in the step five to obtain the numerical solution of the temperature and the electric field intensity of each space node of the nonlinear thermoelectric generation system with the parameters changing along with the temperature, wherein the Seebeck coefficient of the temperature difference thermocouple when the temperature and the electric field intensity of the next time node are obtained is the value corresponding to the temperature of the last time node.
Performing time discretization processing on the area I, the area II, the area III and the area IV in the step three to determine time nodes in the areas; performing one-dimensional space discretization processing on the region I, the region II and the region III, and determining the space nodes of the regions: comprises that
(1) Carrying out time discretization on the region I, the region II, the region III and the region IV, and determining time nodes in the regions:
the time t is divided equally into (N-1) units to obtain N time nodes, and the time interval from one time node to the next time node is delta tau and is called a time step, wherein
Figure BDA0001344238380000031
The nth time node represents the time of N × Δ τ, and N is 0, 1, … …, N;
(2) performing one-dimensional space discretization processing on the region I, the region II and the region III, and determining the space nodes of the regions:
divide region I equally into (N)h-1) units, each unit then having a length Δ xI＝LI/(Nh-1), wherein LIThe length of the region I is defined as the spatial node x at the center of the ith celliThe left and right boundaries are x respectivelyi-1/2And xi+1/21, 2, … …, Nh；
Divide region II equally into (N)w-1) units, each unit then having a length Δ xII＝LII/(Nw-1), wherein LIIFor the length of region II, the center of the j-th cell is the spatial node xjThe left and right boundaries are x respectivelyj-1/2And xj+1/2Taking j as Nh+1，Nh+2，……，Nh+Nw；
Divide region III equally into (N)c-1) units, each unit then having a length Δ xIII＝LIII/(Nc-1), wherein LIIIThe length of region III, the center of the kth element is the spatial node xkThe left and right boundaries arexk-1/2And xk+1/2Taking k as Nh+Nw+1，Nh+Nw+2，……，Nh+Nw+Nc；
And step five, carrying out discretization treatment on the mathematical model, and replacing the values of the space node temperature and the electric field intensity by the average values of the space node temperature and the electric field intensity when establishing an algebraic equation of the temperatures and the electric field intensities of the nodes and the boundary nodes in the region.
The invention provides a time domain analysis method suitable for a nonlinear temperature difference power generation system with parameters changing along with temperature, the method provides a more generally applicable mathematical computation model, firstly establishes a thermoelectric partial differential equation model and boundary conditions based on the basic physical characteristics of the nonlinear temperature difference power generation system with parameters changing along with the temperature, in the model, the Seebeck coefficient is modeled as a function which changes along with the temperature, and the hypothesis of all the functions is related to the original physical characteristics of the actual temperature difference power generation system, secondly, the system hot circuit and the circuit are analyzed in more detail, and are subjected to time and space discretization treatment, an algebraic equation of the physical quantities of the internal nodes and the boundary nodes of the region is established from the energy conservation angle, and then starting from the time initial value, carrying out iterative solution, and finally obtaining the numerical solution of the temperature and the electric field intensity of each position in the nonlinear temperature difference power generation system with the parameters changing along with the temperature. Through the analysis of the numerical solutions of the internal temperature and the electric field intensity of the thermoelectric power generation system, the thermoelectric power generation system can be further analyzed and designed in a more specific mode, so that the model is suitable for time domain analysis of the nonlinear thermoelectric power generation system with the parameters changing along with the temperature, and the result is more accurate.
Drawings
FIG. 1 is a schematic diagram of the basic structure of a thermoelectric power generation system according to the present invention;
FIG. 2 is a simplified schematic diagram of the thermal circuit and electrical circuitry for operation of the thermoelectric generation system of the present invention;
FIG. 3 is a schematic diagram of the discretization process of various regions in accordance with the present invention;
FIG. 4 is a flow chart of the present invention.
Detailed Description
For further explanation of the invention, reference is made to the accompanying drawings
As shown in fig. 1, the thermoelectric power generation system is composed of a hot end radiator 2 in contact with a heat source 1, a ceramic plate 3, a flow deflector 4, a semiconductor thermoelectric couple 5, a cold end radiator 6 in contact with a cold source 7, the ceramic plate 3 and the flow deflector 4. In order to obtain higher electromotive force, the thermoelectric generation system is formed by connecting a pair of semiconductor thermocouples (5) in series with a load resistor (8) through a flow deflector (4) with higher conductivity and a metal wire (9), and is clamped between two ceramic plates (3) which are parallel to each other and have better heat-conducting property and are electrically insulated. The heat flow q flows through the radiator 2 at the hot end, the ceramic plate 3 and the flow deflector 4, and is transmitted to the cold end 7 from the heat source 1 to the semiconductor thermoelectric couple 5, the ceramic plate 3 at the cold end, the flow deflector 4 and the radiator 6 at the cold end, and all the thermoelectric couples are in parallel in the heat flow path. The current J flows into (out) the pair of temperature difference thermocouples 5 through the flow deflector 4 and the metal lead 9 and is connected with the load 8, and all the thermocouples on the current path are connected in series; the thermoelectric coupling is only emitted in the thermoelectric thermocouple.
As shown in fig. 2, since the thermal conductivity of the flow deflector and the ceramic plate is very good, the thermal path of the thermoelectric system can be divided into three regions by neglecting the influence of the flow deflector and the ceramic plate on the thermal path: region I, region II and region III, wherein region I represents a hot side heat sink, region II represents d pairs of semiconductor thermoelectric thermocouples, and region III represents a cold side heat sink; the boundary of the area I and the heat source is a boundary A, the boundary of the area I and the area II is a boundary B, the boundary of the area II and the area III is a boundary C, and the boundary of the area III and the cold source is a boundary D. Because the electric conductivity of the flow deflector and the metal wire is good, the influence of the flow deflector and the metal wire on the circuit is ignored, and the circuit in the system can be divided into two areas: region II and region IV; region II represents d for a semiconductor thermoelectric thermocouple, and region IV represents a load resistor; wherein the boundary between the region II at the heat source end and the region IV is a boundary M, and the boundary between the region II at the heat source end and the region IV is a boundary N.
As shown in fig. 3
(1) Time discretization processing
Tau is a time coordinate, and the time discretization processing is carried out on the area I, the area II, the area III and the area IV
The time t is divided into (N-1) units to obtain N time layers, and the time interval from one time layer to the next time layer is delta tau and is called a time step, wherein
Figure BDA0001344238380000051
The nth time layer represents the time of N × Δ τ, and N is 0, 1, … …, N;
(2) discretization of one-dimensional spatial regions
x is a space coordinate, the one-dimensional space discretization processing is carried out on the region I, the region II and the region III, and the space discretization processing is not carried out on the region IV because the value of the physical quantity is irrelevant to the space size in the region IV
Divide region I equally into (N)h-1) units, each unit then having a length Δ xI＝LI/(Nh-1), wherein LIThe length of the region I is defined as the spatial node x at the center of the ith celliThe left and right boundaries are x respectivelyi-1/2And xi+1/21, 2, … …, Nh；
Divide region II equally into (N)w-1) units, each unit then having a length Δ xII＝LII/(Nw-1), wherein LIIFor the length of region II, the center of the j-th cell is the spatial node xjThe left and right boundaries are x respectivelyj-1/2And xj+1/2Taking j as Nh+1，Nh+2，……，Nh+Nw；
Divide region III equally into (N)c-1) units, each unit then having a length Δ xIII＝LIII/(Nc-1), wherein LIIIThe length of region III, the center of the kth element is the spatial node xkThe left and right boundaries are x respectivelyk-1/2And xk+1/2Taking k as Nh+Nw+1，Nh+Nw+2，……，Nh+Nw+Nc；
(3) Marking of nodes after time and one-dimensional space discretization
For example, (n, m) represents the position in the time-space region, n represents the time node position, and m represents the position of the space node. And m is i, j, k.
As shown in FIG. 4, the time domain analysis method of the nonlinear thermoelectric generation system with the parameter varying with the temperature comprises the following steps:
step 01: according to the basic physical characteristics of the thermoelectric generation system material, the basic laws of thermodynamics and electricity and the specific analysis of the thermoelectric generation system circuit and the thermal circuit, a mathematical model of the nonlinear thermoelectric generation system thermal circuit and the circuit, which is formed by d pairs of semiconductor thermocouples and has the parameter changing with the temperature, is derived:
(1) since zones I and III flow only heat, the mathematical model of the partial differential equation for the temperature change of zones I and III is, according to the law of thermodynamics:
Figure BDA0001344238380000052
and
Figure BDA0001344238380000053
where ρ isIAnd ρIIIDensity of substances in zones I and III, respectively, CvIAnd CvIIIConstant specific heat, k, of the materials in the regions I and III respectivelyIAnd kIIIThermal conductivity of the materials of regions I and III, respectively.
(2) Thermoelectric coupling occurs in a thermoelectric thermocouple, and the heat flow density and current density can be a function of the electric field strength and temperature gradient due to the action of the microparticles in the semiconductor material
Figure BDA0001344238380000061
Figure BDA0001344238380000062
Substituting it into charge storage equation, Gauss's law and energy conservation law
Figure BDA0001344238380000063
Figure BDA0001344238380000064
Wherein, pi is α T,
Figure BDA00013442383800000611
wherein T is the temperature of the molten steel,
Figure BDA00013442383800000612
is a temperature gradient; e is the electric field strength, ρ is the density of the material, CvIs the specific heat of constant volume of the substance, epsilon is the dielectric constant of the substance, J0The current density at time t-0, J, α, σ, k, and pi are the seebeck coefficient, electrical conductivity, thermal conductivity, and peltier coefficient of the material.
So finished, there are:
Figure BDA0001344238380000065
Figure BDA0001344238380000066
where T is the temperature, E is the electric field strength, pIIRespectively, region II material density, CvIIIs the constant specific heat of the substance in the area II, epsilonIIIs the dielectric constant of the region II material, J is the current density, J0Current density at time t-0, αII(T)，σII，kIIThe total equivalent Seebeck coefficient, the total equivalent electric conductivity and the total equivalent thermal conductivity of the substance in the region II, and the region II is formed by connecting d pairs of P-type and N-type semiconductor thermoelectric thermocouples in series electrically and in parallel thermally, so αII(T)＝d*(αP(T)-αN(T))，αP(T)And αN(T)Seebeck coefficients for P-type and N-type semiconductor materials in thermoelectric thermocouples, α, respectivelyP(T)And αN(T)Is a function of temperature variation, αII(T)Also a function of temperature;
Figure BDA0001344238380000067
kII＝m*(kP+kN)；σP，kPthe electric conductivity and the thermal conductivity, sigma, of the P-type semiconductor material in the thermoelectric coupleN，kNThe electric conductivity and the thermal conductivity of the N-type semiconductor material in the thermoelectric thermocouple are respectively.
(3) Due to the load resistance R in the region IVLOnly current flows, according to the basic law of the circuit:
in the region IV, the light-emitting element,
Figure BDA0001344238380000068
σIVis the conductivity of the load resistance.
(4) Therefore, the control equation of each area of the hot circuit is:
in the region I, the first and second electrodes are arranged in a straight line,
Figure BDA0001344238380000069
in the region (II), the first and second regions,
Figure BDA00013442383800000610
in the region of the region III, the first,
Figure BDA0001344238380000071
the control equation for each region of the circuit is:
in the region (II), the first and second regions,
Figure BDA0001344238380000072
Figure BDA0001344238380000073
in the region IV, the light-emitting element,
Figure BDA0001344238380000074
step 02: determining boundary conditions of a thermal circuit and a circuit of the nonlinear thermoelectric power generation system with parameters changing along with the temperature according to actual working conditions:
(1) the boundary conditions in the hot path are:
boundary A, according to the form of energy provided by heat source, the hot end is equivalent to heat flow density of q0；
Boundary B, the density of the heat flow passing through the boundary is continuous, and the temperature is continuous;
a boundary C, wherein the heat flow density passing through the boundary is continuous, and the temperature is continuous;
boundary D: according to the form of the cold end heat dissipation mode, the cold end is equivalent to the temperature Tl0；
(2) The boundary conditions in the circuit are:
a boundary M, wherein the current density passing through the boundary is continuous, and the electric field intensity is continuous;
the current density passing through the boundary is continuous, and the electric field intensity is continuous;
step 03: performing time discretization on the region I, the region II, the region III and the region IV to determine time nodes in the regions, and performing one-dimensional space discretization on the region I, the region II and the region III to determine space nodes in the regions;
(1) time discretization processing
Tau is a time coordinate, and the time discretization processing is carried out on the area I, the area II, the area III and the area IV
The time t is divided equally into (N-1) units to obtain N time nodes, and the time interval from one time node to the next time node is delta tau and is called a time step, wherein
Figure BDA0001344238380000075
The nth time node represents the time of N × Δ τ, and N is 0, 1, … …, N;
(2) discretization of one-dimensional spatial regions
x is a space coordinate, the one-dimensional space discretization processing is carried out on the region I, the region II and the region III, and the space discretization processing is not carried out on the region IV because the value of the physical quantity is irrelevant to the space size in the region IV
Divide region I equally into (N)h-1) units, each unit then having a length Δ xI＝LI/(Nh-1), wherein LIFor the length of region I, N is obtainedhA space node is arranged at the center of the ith unitiThe left and right boundaries are x respectivelyi-1/2And xi+1/21, 2, … …, Nh；
Divide region II equally into (N)w-1) units, each unit then having a length Δ xII＝LII/(Nw-1), wherein LIIFor the length of region II, the center of the j-th cell is the space node xjThe left and right boundaries are x respectivelyj-1/2And xj+1/2Taking j as Nh+1，Nh+2，……，Nh+Nw；
Divide region III equally into (N)c-1) units, each unit then having a length Δ xIII＝LIII/(Nc-1), wherein LIIIFor the length of region III, the center of the k-th element is the spatial node xkThe left and right boundaries are x respectivelyk-1/2And xk+1/2Taking k as Nh+Nw+1，Nh+Nw+2，……，Nh+Nw+Nc；
(3) Marking of node physical quantity after time and one-dimensional space dispersion
EnA value representing the physical quantity E at the nth time node; t ismA value representing a physical quantity T at an m-space node; and m is i, j, k. For example, (n, m) represents the position in the time-space region, n represents the time node position, and m represents the position of the space node.
Figure BDA0001344238380000081
A value representing the nth time mth space node of the physical quantity tth,
Figure BDA0001344238380000082
represents the average value of the mth cell temperature at time n × Δ τ, m ═ i, j, k; likewise, the electric field strength of the node is recorded as
Figure BDA0001344238380000083
Figure BDA0001344238380000084
The mean value of the m-th cell electric field intensity at time n × Δ τ is shown, and m is j. Step 04: determining the iteration initial value of each node temperature and electric field intensity at the moment when t is 0;
firstly, determining the initial values of the temperature T and the electric field intensity E of each unit at the moment when T is equal to 0
Figure BDA0001344238380000085
And
Figure BDA0001344238380000086
setting the heat flux density in the heat path at the initial moment to be the same as q0Cold end temperature of Tl0According to the heat path equation, the initial temperature value of each unit at the initial time is as follows:
Figure BDA0001344238380000087
Figure BDA0001344238380000088
Figure BDA0001344238380000089
the initial value at the boundary is
Figure BDA00013442383800000810
Figure BDA00013442383800000811
Figure BDA00013442383800000812
Figure BDA00013442383800000813
The initial value of the current density is obtained according to the circuit equation
Figure BDA00013442383800000814
Further obtaining the initial value of the electric field intensity of the space node at the initial time
Figure BDA00013442383800000815
Figure BDA0001344238380000091
Wherein, the boundary C and D of the area II in the hot circuit are overlapped with the boundary M and N of the area II in the circuit. Therefore it has the advantages of
Figure BDA0001344238380000092
Figure BDA0001344238380000093
Step 05: discretizing the mathematical model, and establishing an algebraic equation of the temperature and the electric field intensity of the internal nodes and the boundary nodes of the region;
(1) algebraic equation for establishing node temperature and electric field intensity in region
The average value of the physical quantity of the space unit is used for replacing the value of the physical quantity of the node, and according to the Taylor series expansion formula, the method comprises the following steps:
Figure BDA0001344238380000094
Figure BDA0001344238380000095
Figure BDA0001344238380000096
Figure BDA0001344238380000097
Figure BDA0001344238380000098
for the region I, the partial differential equation is subjected to one-dimensional discretization treatment to obtain the final product
Figure BDA0001344238380000099
Figure BDA00013442383800000910
Figure BDA00013442383800000911
Then:
Figure BDA00013442383800000912
in the same way as above, the first and second,
for the region III, the partial differential equation is subjected to one-dimensional discretization treatment to obtain
Figure BDA00013442383800000913
For region II, the one-dimensional form of the partial differential equation is
Figure BDA00013442383800000914
Figure BDA00013442383800000915
Discretized processing can obtain
Figure BDA0001344238380000101
Figure BDA0001344238380000102
Figure BDA0001344238380000103
Figure BDA0001344238380000104
Then
Figure BDA0001344238380000105
Figure BDA0001344238380000106
Figure BDA0001344238380000107
For time node have
Figure BDA0001344238380000108
Figure BDA0001344238380000109
Finally, an algebraic equation of the temperature of the node in the region and the electric field intensity can be sorted out as follows:
in region I
Figure BDA00013442383800001010
Region II
Figure BDA0001344238380000111
Figure BDA0001344238380000112
Figure BDA0001344238380000113
Region III
Figure BDA0001344238380000114
In the region IV, the light-emitting element,
Figure BDA0001344238380000115
(2) establishing an algebraic equation of the physical quantity of the region boundary nodes:
at boundary A, the hot side is equivalent to a heat flux density of q, depending on the form of the energy provided by the heat source0，
Figure BDA0001344238380000116
At the boundary D, the cold source end is equivalent to the temperature T according to the form of the cold end heat dissipation model0；
Figure BDA0001344238380000117
For boundary B, the heat flux density is continuous at any time according to the boundary condition of the thermal circuit, and the temperature is continuous at any time
Figure BDA0001344238380000118
Therefore it has the advantages of
Figure BDA0001344238380000119
Discretizing to obtain i ═ Nh，j＝Nh+1
Figure BDA00013442383800001110
Figure BDA00013442383800001111
For the boundary M, the current density through the boundary at any one time is continuous according to the circuit equation and the boundary conditions in the circuit, and the electric field strength is continuous.
Is provided with
Figure BDA0001344238380000121
Therefore it has the advantages of
Figure BDA0001344238380000122
Since the boundary M of the circuit and the boundary B of the hot circuit are in the same position in region II,
therefore, therein
Figure BDA0001344238380000123
Figure BDA0001344238380000124
The boundary B can be obtained by simultaneous equations (1-1) and (1-2)
Figure BDA0001344238380000125
Figure BDA0001344238380000126
Similarly, for boundary C, the heat flux density is continuous at any time and the temperature is continuous at any time according to the boundary conditions of the thermal circuit
Is provided with
Figure BDA0001344238380000127
Therefore it has the advantages of
Figure BDA0001344238380000128
Discretizing, where j is Nh+Nw，k＝Nh+Nw+1
Figure BDA0001344238380000129
Figure BDA00013442383800001210
For the boundary N, the current density through the boundary at any one time is continuous according to the circuit equation and the boundary conditions in the circuit, and the electric field strength is continuous.
Is provided with
Figure BDA00013442383800001211
Therefore it has the advantages of
Figure BDA00013442383800001212
Since the boundary N of the circuit and the boundary C of the hot circuit are in the same position in the region II,
therefore, therein
Figure BDA00013442383800001213
Figure BDA00013442383800001214
The boundary C can be obtained by simultaneous equations (2-1) and (2-2)
Figure BDA0001344238380000131
Figure BDA0001344238380000132
Finally, the discretization algebraic equation of the temperature and the electric field intensity of each node in the region can be sorted out as follows:
in region I
Figure BDA0001344238380000133
Figure BDA0001344238380000134
Figure BDA0001344238380000135
In the region (II), the first and second regions,
Figure BDA0001344238380000136
Figure BDA0001344238380000137
Figure BDA0001344238380000138
Figure BDA0001344238380000139
Figure BDA0001344238380000141
Figure BDA0001344238380000142
Figure BDA0001344238380000143
Figure BDA0001344238380000144
Figure BDA0001344238380000145
wherein the content of the first and second substances,
Figure BDA0001344238380000146
Figure BDA0001344238380000147
(i＝Nh，j＝Nh+1)
Figure BDA0001344238380000148
Figure BDA0001344238380000149
(j＝Nh+Nw，k＝Nh+Nw+1)
and 06, starting from the iteration initial values of the space units with the temperatures and the electric field intensities, carrying out iteration solution according to the algebraic equations of the node temperatures and the electric field intensities established in the step 05 to obtain numerical solutions of the temperature and the electric field intensities of the space nodes of the nonlinear thermoelectric generation system with the parameters changing along with the temperatures, wherein when the values of the node temperatures and the electric field intensities in (n +1) time are solved, the Seebeck coefficient of the thermoelectric couple takes the corresponding value at the time of n x delta tau under the average temperature of the space node units
Figure BDA00013442383800001410

Claims (3)

1. A time domain analysis method of a nonlinear temperature difference power generation system with parameters changing along with temperature is characterized in that: the method specifically comprises the following steps: step one, deriving a mathematical model of a nonlinear thermoelectric generation system thermal circuit and a circuit, wherein the mathematical model is composed of d pairs of semiconductor thermoelectric thermocouples and has parameters changing with temperature, according to basic physical characteristics of materials of the thermoelectric generation system and basic laws of thermodynamics and electrics:
the mathematical model of each region in the hot circuit is:
in the region I, the first and second electrodes are arranged in a straight line,
Figure FDA0002502228370000011
in the region (II), the first and second regions,
Figure FDA0002502228370000012
in the region of the region III, the first,
Figure FDA0002502228370000013
the mathematical model of each region of the circuit is:
in the region (II), the first and second regions,
Figure FDA0002502228370000014
Figure FDA0002502228370000015
in the region IV, the light-emitting element,
Figure FDA0002502228370000016
wherein, thermoelectric generation system heat route divides into three regions: region I, region II and region III, region I representing a hot side heat sink, region II representing d a semiconductor thermoelectric thermocouple, region III representing a cold side heat sink; the circuit in the thermoelectric power generation system is divided into two areas: region II and region IV; region II represents d versus semiconductor thermoelectric thermocouple, and region IV represents load resistance RL；
Where ρ isIAnd ρIIIDensity of substances in zones I and III, respectively, CvIAnd CvIIIConstant specific heat, k, of the materials in the regions I and III respectivelyIAnd kIIIThermal conductivity of materials in regions I and III, respectively; t is the temperature of the molten metal,
Figure FDA0002502228370000017
is a temperature gradient; e is the electric field strength, ρIIIs the density of the substance of region II, CvIIIs the constant specific heat of the substance in the area II, epsilonIIIs the dielectric constant of the material in region II, and J is the current density in the circuit, current J0Current density through the circuit at time t-0, αII(T)，σII，kIIα is the total equivalent Seebeck coefficient, the total equivalent electrical conductivity and the total equivalent thermal conductivity of the substances in the area IIII(T)＝d*(αP(T)-αN(T))，αP(T)And αN(T)Seebeck coefficients for P-type and N-type semiconductor materials in semiconductor thermoelectric thermocouples, α, respectivelyP(T)And αN(T)Is a function of temperature, therefore αII(T)Also a function of temperature;
Figure FDA0002502228370000018
kII＝m*(kP+kN)；σP，kPthe electric conductivity and the thermal conductivity, sigma, of the P-type semiconductor material in the semiconductor thermoelectric coupleN，kNThe electric conductivity and the thermal conductivity of the N-type semiconductor material in the semiconductor thermoelectric thermocouple are respectively; sigmaIVIs the conductivity of the load resistance;
determining the boundary conditions of the thermal circuit and the circuit of the nonlinear thermoelectric power generation system with parameters changing along with the temperature according to the actual working condition:
the boundary conditions in the hot path are:
boundary A, according to the form of energy provided by heat source, the hot end is equivalent to heat flow density of q0；
Boundary B, the density of the heat flow passing through the boundary is continuous, and the temperature is continuous;
a boundary C, wherein the heat flow density passing through the boundary is continuous, and the temperature is continuous;
boundary D: according to the heat dissipation form of the cold source end, the cold end is equivalent to the temperature of Tl0；
The boundary conditions in the circuit are:
a boundary M, wherein the current density passing through the boundary is continuous, and the electric field intensity is continuous;
the current density passing through the boundary is continuous, and the electric field intensity is continuous;
the boundary of the area I and the heat source in the heat path is a boundary A, the boundary of the area I and the area II is a boundary B, the boundary of the area II and the area III is a boundary C, and the boundary of the area III and the cold source is a boundary D; the boundary of the area II at the heat source end and the area IV in the circuit is a boundary M, and the boundary of the area II at the cold source end and the area IV is a boundary N;
performing time discretization treatment on the areas I, II, III and IV to determine time nodes in the areas; performing one-dimensional space discretization processing on the region I, the region II and the region III to determine space nodes of the regions;
step four, determining the temperature of each space node and the iteration initial value of the electric field intensity at the moment when t is 0;
discretizing the mathematical model, and establishing an algebraic equation of the temperature and the electric field intensity of the nodes and the boundary nodes in the region;
and step six, starting from the iteration initial values of the space nodes of the temperature and the electric field intensity, carrying out iteration solution according to the algebraic equations of the node temperature and the electric field intensity established in the step five to obtain the numerical solution of the temperature and the electric field intensity of each space node of the nonlinear temperature difference power generation system with the parameters changing along with the temperature, wherein when the temperature and the electric field intensity of the next time node are solved, the Seebeck coefficient of the semiconductor thermoelectric couple is the value corresponding to the average temperature of the space unit of the last time node.
2. The time domain analysis method of the nonlinear thermoelectric power generation system with the parameter varying with the temperature according to claim 1, characterized in that: in the third step, the time discretization processing is performed on the area I, the area II, the area III, and the area IV, and the time node in the area is determined, including:
dividing time t into (N-1) units to obtain N time nodes, wherein a time interval from one time node to the next time node is Δ τ, which is called a time step, where Δ τ is t/N, the nth time node represents time N × Δ τ, and N is 0, 1, … …, N;
divide region I equally into (N)h-1) units, each unit then having a length Δ xI＝LI/(Nh-1), wherein LIThe length of the region I is defined as the spatial node x at the center of the ith celliThe left and right boundaries are x respectivelyi-1/2And xi+1/21, 2, … …, Nh；
Divide region II equally into (N)w-1) units, each unit then having a length Δ xII＝LII/(Nw-1), wherein LIIFor the length of region II, the center of the j-th cell is the spatial node xjThe left and right boundaries are x respectivelyj-1/2And xj+1/2Taking j as Nh+1，Nh+2，……，Nh+Nw；
Divide region III equally into (N)c-1) units, each unit then having a length Δ xIII＝LIII/(Nc-1), wherein LIIIThe length of region III, the center of the kth element is the spatial node xkThe left and right boundaries are x respectivelyk-1/2And xk+1/2Taking k as Nh+Nw+1，Nh+Nw+2，……，Nh+Nw+Nc。
3. The time domain analysis method of the nonlinear thermoelectric power generation system with the parameter varying with the temperature according to claim 1, characterized in that: and discretizing the mathematical model, establishing an algebraic equation of the temperature and the electric field intensity of the nodes and the boundary nodes in the region, and replacing the average value of the temperature and the electric field intensity of the space unit with the average value of the temperature and the electric field intensity of the space unit.
CN201710550251.0A 2017-07-07 2017-07-07 Time domain analysis method of nonlinear thermoelectric power generation system with parameter changing along with temperature Expired - Fee Related CN107330210B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710550251.0A CN107330210B (en) 2017-07-07 2017-07-07 Time domain analysis method of nonlinear thermoelectric power generation system with parameter changing along with temperature

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710550251.0A CN107330210B (en) 2017-07-07 2017-07-07 Time domain analysis method of nonlinear thermoelectric power generation system with parameter changing along with temperature

Publications (2)

Publication Number Publication Date
CN107330210A CN107330210A (en) 2017-11-07
CN107330210B true CN107330210B (en) 2020-07-07

Family

ID=60196462

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710550251.0A Expired - Fee Related CN107330210B (en) 2017-07-07 2017-07-07 Time domain analysis method of nonlinear thermoelectric power generation system with parameter changing along with temperature

Country Status (1)

Country Link
CN (1) CN107330210B (en)

Family Cites Families (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2012061763A2 (en) * 2010-11-05 2012-05-10 Bsst Llc Energy management systems and methods with thermoelectric generators
JP5862300B2 (en) * 2011-12-30 2016-02-16 Ｊｆｅスチール株式会社 Thermoelectric generation method
CN103973171B (en) * 2014-05-06 2016-01-27 华北水利水电大学 A kind of thermo-electric generation system electromotive force computational methods
CN104993740B (en) * 2015-07-07 2017-08-08 天津大学 A kind of segmented thermoelectric generator construction design method
CN106533263A (en) * 2016-12-29 2017-03-22 哈尔滨工业大学 Multi-stage semiconductor thermoelectric power generation and cooling integrated system for hypersonic flight vehicle

Also Published As

Publication number Publication date
CN107330210A (en) 2017-11-07

Similar Documents

Publication Publication Date Title
Dongxu et al. Geometry optimization of thermoelectric modules: Simulation and experimental study
Montecucco et al. Constant heat characterisation and geometrical optimisation of thermoelectric generators
Chen et al. Effect of pulsed heat power on the thermal and electrical performances of a thermoelectric generator
Arora et al. Multi-objective and multi-parameter optimization of two-stage thermoelectric generator in electrically series and parallel configurations through NSGA-II
Massaguer et al. Development and validation of a new TRNSYS type for the simulation of thermoelectric generators
Högblom et al. A simulation framework for prediction of thermoelectric generator system performance
Niu et al. Elucidating modeling aspects of thermoelectric generator
Jia et al. Estimation of thermoelectric and mechanical performances of segmented thermoelectric generators under optimal operating conditions
Sandoz-Rosado et al. On the Thomson effect in thermoelectric power devices
Favarel et al. Numerical optimization of the occupancy rate of thermoelectric generators to produce the highest electrical power
Wang A finite element computational scheme for transient and nonlinear coupling thermoelectric fields and the associated thermal stresses in thermoelectric materials
Ahiska et al. Design and implementation of a new portable thermoelectric generator for low geothermal temperatures
CN103973171B (en) A kind of thermo-electric generation system electromotive force computational methods
Wang et al. Temperature and performance modeling of thermoelectric generators
Yilbas et al. Thermoelectric generator performance analysis: Influence of pin tapering on the first and second law efficiencies
Qing et al. An analytical model for performance optimization of thermoelectric generator with temperature dependent materials
Brito et al. Analysis of the effect of module thickness reduction on thermoelectric generator output
Reddy et al. Three-dimensional multiphysics coupled field analysis of an integrated thermoelectric device
He et al. Effects of heat transfer characteristics between fluid channels and thermoelectric modules on optimal thermoelectric performance
Colomer et al. Electrically tunable thermal conductivity in thermoelectric materials: Active and passive control
Kane et al. Temperature dependent analysis of thermoelectric module using Matlab/SIMULINK
Orr et al. Validating an alternative method to predict thermoelectric generator performance
CN107194130B (en) Time domain analysis method for hot-end heat flow time-varying dynamic thermoelectric power generation system
Palaniappan et al. Finite element analysis of thermoelectric refrigeration system
CN107122576B (en) Time domain analysis method for thermoelectric couple material segmented dynamic thermoelectric power generation system

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20200707

Termination date: 20210707